NAMD
Public Member Functions | Protected Member Functions | Protected Attributes | Friends | List of all members
colvarproxy_namd Class Reference

Communication between colvars and NAMD (implementation of colvarproxy) More...

#include <colvarproxy_namd.h>

Inheritance diagram for colvarproxy_namd:
GlobalMaster

Public Member Functions

 colvarproxy_namd ()
 
 ~colvarproxy_namd ()
 
int setup ()
 
int reset ()
 
int update_atoms_map (AtomIDList::const_iterator begin, AtomIDList::const_iterator end)
 
void calculate ()
 
void log (std::string const &message)
 
void error (std::string const &message)
 
int set_unit_system (std::string const &units_in, bool check_only)
 
void exit (std::string const &message)
 
void add_energy (cvm::real energy)
 
void request_total_force (bool yesno)
 
bool total_forces_enabled () const
 
int run_force_callback ()
 
int run_colvar_callback (std::string const &name, std::vector< const colvarvalue * > const &cvcs, colvarvalue &value)
 
int run_colvar_gradient_callback (std::string const &name, std::vector< const colvarvalue * > const &cvcs, std::vector< cvm::matrix2d< cvm::real > > &gradient)
 
cvm::real backend_angstrom_value ()
 
cvm::real boltzmann ()
 
cvm::real temperature ()
 
cvm::real rand_gaussian ()
 
cvm::real dt ()
 
virtual int replica_enabled ()
 
virtual int replica_index ()
 
virtual int num_replicas ()
 
virtual void replica_comm_barrier ()
 
virtual int replica_comm_recv (char *msg_data, int buf_len, int src_rep)
 
virtual int replica_comm_send (char *msg_data, int msg_len, int dest_rep)
 
int init_atom (int atom_number)
 
int check_atom_id (int atom_number)
 
int init_atom (cvm::residue_id const &residue, std::string const &atom_name, std::string const &segment_id)
 
int check_atom_id (cvm::residue_id const &residue, std::string const &atom_name, std::string const &segment_id)
 
void clear_atom (int index)
 
void update_atom_properties (int index)
 
cvm::rvector position_distance (cvm::atom_pos const &pos1, cvm::atom_pos const &pos2) const
 
int load_atoms (char const *filename, cvm::atom_group &atoms, std::string const &pdb_field, double const pdb_field_value=0.0)
 
int load_coords (char const *filename, std::vector< cvm::atom_pos > &pos, const std::vector< int > &indices, std::string const &pdb_field, double const pdb_field_value=0.0)
 
int scalable_group_coms ()
 
int init_atom_group (std::vector< int > const &atoms_ids)
 
void clear_atom_group (int index)
 
int update_group_properties (int index)
 
int init_volmap (int volmap_id)
 
int init_volmap (const char *volmap_name)
 
void clear_volmap (int index)
 
std::ostream * output_stream (std::string const &output_name, std::ios_base::openmode mode)
 
int flush_output_stream (std::ostream *os)
 
int close_output_stream (std::string const &output_name)
 
int backup_file (char const *filename)
 
char const * script_obj_to_str (unsigned char *obj)
 
std::vector< std::string > script_obj_to_str_vector (unsigned char *obj)
 
- Public Member Functions inherited from GlobalMaster
void processData (AtomIDList::iterator a_i, AtomIDList::iterator a_e, PositionList::iterator p_i, PositionList::iterator g_i, PositionList::iterator g_e, BigRealList::iterator gm_i, BigRealList::iterator gm_e, ForceList::iterator gtf_i, ForceList::iterator gtf_e, IntList::iterator goi_i, IntList::iterator goi_e, BigRealList::iterator gov_i, BigRealList::iterator gov_e, AtomIDList::iterator last_atoms_forced_i, AtomIDList::iterator last_atoms_forced_e, ForceList::iterator last_forces_i, AtomIDList::iterator, AtomIDList::iterator, ForceList::iterator)
 
bool changedAtoms ()
 
const AtomIDListrequestedAtoms ()
 
bool changedForces ()
 
const AtomIDListforcedAtoms ()
 
const ForceListappliedForces ()
 
bool changedGroups ()
 
const ResizeArray< AtomIDList > & requestedGroups ()
 
const ForceListgroupForces ()
 
bool changedGridObjs ()
 
const IntListrequestedGridObjs ()
 
const BigRealListgridObjForces ()
 
bool requestedTotalForces ()
 
void clearChanged ()
 
virtual ~GlobalMaster ()
 
void check () const
 
void setLattice (const Lattice *lat)
 

Protected Member Functions

void init_tcl_pointers ()
 
- Protected Member Functions inherited from GlobalMaster
 GlobalMaster ()
 
AtomIDListmodifyRequestedAtoms ()
 
AtomIDListmodifyForcedAtoms ()
 
ForceListmodifyAppliedForces ()
 
ResizeArray< AtomIDList > & modifyRequestedGroups ()
 
ForceListmodifyGroupForces ()
 
IntListmodifyRequestedGridObjects ()
 
BigRealListmodifyGridObjForces ()
 
AtomIDList::const_iterator getAtomIdBegin ()
 
AtomIDList::const_iterator getAtomIdEnd ()
 
PositionList::const_iterator getAtomPositionBegin ()
 
PositionList::const_iterator getGroupPositionBegin ()
 
PositionList::const_iterator getGroupPositionEnd ()
 
ForceList::const_iterator getGroupTotalForceBegin ()
 
ForceList::const_iterator getGroupTotalForceEnd ()
 
IntList::const_iterator getGridObjIndexBegin ()
 
IntList::const_iterator getGridObjIndexEnd ()
 
BigRealList::const_iterator getGridObjValueBegin ()
 
BigRealList::const_iterator getGridObjValueEnd ()
 
AtomIDList::const_iterator getLastAtomsForcedBegin ()
 
AtomIDList::const_iterator getLastAtomsForcedEnd ()
 
ForceList::const_iterator getLastForcesBegin ()
 
AtomIDList::const_iterator getForceIdBegin ()
 
AtomIDList::const_iterator getForceIdEnd ()
 
ForceList::const_iterator getTotalForce ()
 
void requestTotalForce (bool yesno=true)
 
BigRealList::const_iterator getGroupMassBegin ()
 
BigRealList::const_iterator getGroupMassEnd ()
 

Protected Attributes

std::vector< int > atoms_map
 Array of atom indices (relative to the colvarproxy arrays), usedfor faster copy of atomic data. More...
 
SimParameterssimparams
 Pointer to the NAMD simulation input object. More...
 
BigReal thermostat_temperature
 Self-explained. More...
 
Random random
 NAMD-style PRNG object. More...
 
bool first_timestep
 
size_t previous_NAMD_step
 
SubmitReductionreduction
 Used to submit restraint energy as MISC. More...
 
- Protected Attributes inherited from GlobalMaster
bool totalForceRequested
 
const Latticelattice
 
AtomIDList::iterator atomIdBegin
 
AtomIDList::iterator atomIdEnd
 
PositionList::iterator atomPositionBegin
 
PositionList::iterator groupPositionBegin
 
PositionList::iterator groupPositionEnd
 
BigRealList::iterator groupMassBegin
 
BigRealList::iterator groupMassEnd
 
ForceList::iterator groupTotalForceBegin
 
ForceList::iterator groupTotalForceEnd
 
IntList::iterator gridObjIndexBegin
 
IntList::iterator gridObjIndexEnd
 
BigRealList::iterator gridObjValueBegin
 
BigRealList::iterator gridObjValueEnd
 
AtomIDList::iterator lastAtomsForcedBegin
 
ForceList::iterator lastForcesBegin
 
AtomIDList::iterator lastAtomsForcedEnd
 
AtomIDList::iterator forceIdBegin
 
AtomIDList::iterator forceIdEnd
 
ForceList::iterator totalForceBegin
 
bool reqAtomsChanged
 
AtomIDList reqAtoms
 
bool appForcesChanged
 
AtomIDList fAtoms
 
ForceList appForces
 
bool reqGroupsChanged
 
ResizeArray< AtomIDListreqGroups
 
ForceList grpForces
 
bool reqGridObjsChanged
 
IntList reqGridObjs
 
BigRealList gridobjForces
 

Friends

class cvm::atom
 

Additional Inherited Members

- Public Attributes inherited from GlobalMaster
int step
 
int old_num_groups_requested
 

Detailed Description

Communication between colvars and NAMD (implementation of colvarproxy)

Definition at line 34 of file colvarproxy_namd.h.

Constructor & Destructor Documentation

colvarproxy_namd::colvarproxy_namd ( )

Definition at line 38 of file colvarproxy_namd.C.

References COLVARPROXY_VERSION, Node::colvars, Node::configList, StringList::data, debug, endi(), error(), ConfigList::find(), first_timestep, SimParameters::firstTimestep, init_tcl_pointers(), iout, SimParameters::langevinOn, SimParameters::langevinTemp, SimParameters::loweAndersenOn, SimParameters::loweAndersenTemp, Node::Object(), ReductionMgr::Object(), SimParameters::outputFilename, random, SimParameters::randomSeed, SimParameters::reassignFreq, SimParameters::reassignTemp, reduction, REDUCTIONS_BASIC, GlobalMaster::requestTotalForce(), SimParameters::rescaleFreq, SimParameters::rescaleTemp, SimParameters::restartFilename, SimParameters::restartFrequency, Node::simParameters, simparams, SimParameters::stochRescaleOn, SimParameters::stochRescaleTemp, SimParameters::tCoupleOn, SimParameters::tCoupleTemp, thermostat_temperature, and ReductionMgr::willSubmit().

39 {
40  version_int = get_version_from_string(COLVARPROXY_VERSION);
41 
42  first_timestep = true;
43  requestTotalForce(total_force_requested);
44 
45  angstrom_value = 1.;
46 
47  // initialize pointers to NAMD configuration data
49 
50  if (cvm::debug())
51  iout << "Info: initializing the colvars proxy object.\n" << endi;
52 
53  // find the configuration file, if provided
54  StringList *config = Node::Object()->configList->find("colvarsConfig");
55 
56  // find the input state file
57  StringList *input_restart = Node::Object()->configList->find("colvarsInput");
58  input_prefix_str = std::string(input_restart ? input_restart->data : "");
59  if (input_prefix_str.rfind(".colvars.state") != std::string::npos) {
60  // strip the extension, if present
61  input_prefix_str.erase(input_prefix_str.rfind(".colvars.state"),
62  std::string(".colvars.state").size());
63  }
64 
65  // get the thermostat temperature
66  if (simparams->rescaleFreq > 0)
68  else if (simparams->reassignFreq > 0)
70  else if (simparams->langevinOn)
72  else if (simparams->tCoupleOn)
74  else if (simparams->loweAndersenOn)
76  else if (simparams->stochRescaleOn)
78  else
80 
82 
83  // both fields are taken from data structures already available
84  updated_masses_ = updated_charges_ = true;
85 
86  // take the output prefixes from the namd input
87  output_prefix_str = std::string(simparams->outputFilename);
88  restart_output_prefix_str = std::string(simparams->restartFilename);
89  restart_frequency_engine = simparams->restartFrequency;
90 
91  // check if it is possible to save output configuration
92  if ((!output_prefix_str.size()) && (!restart_output_prefix_str.size())) {
93  error("Error: neither the final output state file or "
94  "the output restart file could be defined, exiting.\n");
95  }
96 
97 
98 #ifdef NAMD_TCL
99  have_scripts = true;
100 
102 
103  // See is user-scripted forces are defined
104  if (Tcl_FindCommand(reinterpret_cast<Tcl_Interp *>(tcl_interp_),
105  "calc_colvar_forces", NULL, 0) == NULL) {
106  force_script_defined = false;
107  } else {
108  force_script_defined = true;
109  }
110 #else
111  force_script_defined = false;
112  have_scripts = false;
113 #endif
114 
115 
116  // initialize module: this object will be the communication proxy
117  colvars = new colvarmodule(this);
118  cvm::log("Using NAMD interface, version "+
119  cvm::to_str(COLVARPROXY_VERSION)+".\n");
120 
121  errno = 0;
122  if (config) {
123  colvars->read_config_file(config->data);
124  }
125 
126  colvars->setup();
127  colvars->setup_input();
128  colvars->setup_output();
129 
130  // save to Node for Tcl script access
131  Node::Object()->colvars = colvars;
132 
133 #ifdef NAMD_TCL
134  // Construct instance of colvars scripting interface
135  script = new colvarscript(this);
136 #endif
137 
138  if (simparams->firstTimestep != 0) {
139  cvm::log("Initializing step number as firstTimestep.\n");
140  colvars->it = colvars->it_restart =
141  static_cast<cvm::step_number>(simparams->firstTimestep);
142  }
143 
145 
146  if (cvm::debug())
147  iout << "Info: done initializing the colvars proxy object.\n" << endi;
148 }
static Node * Object()
Definition: Node.h:86
char restartFilename[128]
SimParameters * simparams
Pointer to the NAMD simulation input object.
Random random
NAMD-style PRNG object.
SimParameters * simParameters
Definition: Node.h:178
BigReal reassignTemp
SubmitReduction * willSubmit(int setID, int size=-1)
Definition: ReductionMgr.C:365
static int debug
Definition: parm.C:36
static ReductionMgr * Object(void)
Definition: ReductionMgr.h:278
#define iout
Definition: InfoStream.h:87
#define COLVARPROXY_VERSION
char outputFilename[128]
Definition: Random.h:37
BigReal rescaleTemp
BigReal langevinTemp
ConfigList * configList
Definition: Node.h:179
void error(std::string const &message)
unsigned int randomSeed
void requestTotalForce(bool yesno=true)
Definition: GlobalMaster.h:133
char * data
Definition: ConfigList.h:48
StringList * find(const char *name) const
Definition: ConfigList.C:341
BigReal thermostat_temperature
Self-explained.
colvarmodule * colvars
Definition: Node.h:184
BigReal tCoupleTemp
infostream & endi(infostream &s)
Definition: InfoStream.C:38
SubmitReduction * reduction
Used to submit restraint energy as MISC.
BigReal stochRescaleTemp
BigReal loweAndersenTemp
colvarproxy_namd::~colvarproxy_namd ( )

Definition at line 151 of file colvarproxy_namd.C.

References reduction.

152 {
153  delete reduction;
154  if (script != NULL) {
155  delete script;
156  script = NULL;
157  }
158  if (colvars != NULL) {
159  delete colvars;
160  colvars = NULL;
161  }
162 }
SubmitReduction * reduction
Used to submit restraint energy as MISC.

Member Function Documentation

void colvarproxy_namd::add_energy ( cvm::real  energy)

Definition at line 596 of file colvarproxy_namd.C.

References SubmitReduction::item(), reduction, and REDUCTION_MISC_ENERGY.

597 {
599 }
BigReal & item(int i)
Definition: ReductionMgr.h:312
SubmitReduction * reduction
Used to submit restraint energy as MISC.
cvm::real colvarproxy_namd::backend_angstrom_value ( )
inline

Definition at line 95 of file colvarproxy_namd.h.

96  {
97  return 1.0;
98  }
int colvarproxy_namd::backup_file ( char const *  filename)

Definition at line 1076 of file colvarproxy_namd.C.

References NAMD_backup_file().

1077 {
1078  if (std::string(filename).rfind(std::string(".colvars.state")) != std::string::npos) {
1079  NAMD_backup_file(filename, ".old");
1080  } else {
1081  NAMD_backup_file(filename, ".BAK");
1082  }
1083  return COLVARS_OK;
1084 }
void NAMD_backup_file(const char *filename, const char *extension)
Definition: common.C:159
cvm::real colvarproxy_namd::boltzmann ( )
inline

Definition at line 100 of file colvarproxy_namd.h.

101  {
102  return 0.001987191;
103  }
void colvarproxy_namd::calculate ( )
virtual

Reimplemented from GlobalMaster.

Definition at line 278 of file colvarproxy_namd.C.

References Lattice::a(), Lattice::a_p(), ResizeArray< T >::add(), atoms_map, Lattice::b(), Lattice::b_p(), ResizeArray< T >::begin(), Lattice::c(), Lattice::c_p(), ResizeArray< T >::clear(), debug, ResizeArray< T >::end(), first_timestep, GlobalMaster::getAtomIdBegin(), GlobalMaster::getAtomIdEnd(), GlobalMaster::getAtomPositionBegin(), GlobalMaster::getForceIdBegin(), GlobalMaster::getForceIdEnd(), GlobalMaster::getGridObjIndexBegin(), GlobalMaster::getGridObjIndexEnd(), GlobalMaster::getGridObjValueBegin(), GlobalMaster::getGridObjValueEnd(), GlobalMaster::getGroupPositionBegin(), GlobalMaster::getGroupPositionEnd(), GlobalMaster::getGroupTotalForceBegin(), GlobalMaster::getGroupTotalForceEnd(), GlobalMaster::getTotalForce(), GlobalMaster::lattice, log(), GlobalMaster::modifyAppliedForces(), GlobalMaster::modifyForcedAtoms(), GlobalMaster::modifyGridObjForces(), GlobalMaster::modifyGroupForces(), Node::molecule, SimParameters::N, Molecule::numAtoms, Node::Object(), Lattice::orthogonal(), previous_NAMD_step, reduction, GlobalMaster::requestedGridObjs(), GlobalMaster::requestedGroups(), ResizeArray< T >::resize(), Vector::set(), ResizeArray< T >::setall(), setup(), simparams, GlobalMaster::step, SubmitReduction::submit(), update_atoms_map(), Vector::x, Vector::y, and Vector::z.

279 {
280  errno = 0;
281 
282  if (first_timestep) {
283 
285  colvars->setup();
286  colvars->setup_input();
287  colvars->setup_output();
288 
289  first_timestep = false;
290 
291  } else {
292  // Use the time step number inherited from GlobalMaster
293  if ( step - previous_NAMD_step == 1 ) {
294  colvars->it++;
295  b_simulation_continuing = false;
296  } else {
297  // Cases covered by this condition:
298  // - run 0
299  // - beginning of a new run statement
300  // The internal counter is not incremented, and the objects are made
301  // aware of this via the following flag
302  b_simulation_continuing = true;
303  }
304  }
305 
307 
308  {
309  Vector const a = lattice->a();
310  Vector const b = lattice->b();
311  Vector const c = lattice->c();
312  unit_cell_x.set(a.x, a.y, a.z);
313  unit_cell_y.set(b.x, b.y, c.z);
314  unit_cell_z.set(c.x, c.y, c.z);
315  }
316 
317  if (!lattice->a_p() && !lattice->b_p() && !lattice->c_p()) {
318  boundaries_type = boundaries_non_periodic;
319  reset_pbc_lattice();
320  } else if (lattice->a_p() && lattice->b_p() && lattice->c_p()) {
321  if (lattice->orthogonal()) {
322  boundaries_type = boundaries_pbc_ortho;
323  } else {
324  boundaries_type = boundaries_pbc_triclinic;
325  }
326  colvarproxy_system::update_pbc_lattice();
327  } else {
328  boundaries_type = boundaries_unsupported;
329  }
330 
331  if (cvm::debug()) {
332  cvm::log(std::string(cvm::line_marker)+
333  "colvarproxy_namd, step no. "+cvm::to_str(colvars->it)+"\n"+
334  "Updating atomic data arrays.\n");
335  }
336 
337  // must delete the forces applied at the previous step: we can do
338  // that because they have already been used and copied to other
339  // memory locations
342 
343  // prepare local arrays
344  for (size_t i = 0; i < atoms_ids.size(); i++) {
345  atoms_positions[i] = cvm::rvector(0.0, 0.0, 0.0);
346  atoms_total_forces[i] = cvm::rvector(0.0, 0.0, 0.0);
347  atoms_new_colvar_forces[i] = cvm::rvector(0.0, 0.0, 0.0);
348  }
349 
350  for (size_t i = 0; i < atom_groups_ids.size(); i++) {
351  atom_groups_total_forces[i] = cvm::rvector(0.0, 0.0, 0.0);
352  atom_groups_new_colvar_forces[i] = cvm::rvector(0.0, 0.0, 0.0);
353  }
354 
355 #if NAMD_VERSION_NUMBER >= 34471681
356  for (int imap = 0; imap < volmaps_ids.size(); imap++) {
357  volmaps_new_colvar_forces[imap] = 0.0;
358  }
359 #endif
360 
361  // create the atom map if needed
362  size_t const n_all_atoms = Node::Object()->molecule->numAtoms;
363  if (atoms_map.size() != n_all_atoms) {
364  atoms_map.resize(n_all_atoms);
365  atoms_map.assign(n_all_atoms, -1);
367  }
368 
369  // if new atomic positions or forces have been communicated by other GlobalMasters, add them to the atom map
370  if ((int(atoms_ids.size()) < (getAtomIdEnd() - getAtomIdBegin())) ||
371  (int(atoms_ids.size()) < (getForceIdEnd() - getForceIdBegin()))) {
374  }
375 
376  {
377  if (cvm::debug()) {
378  cvm::log("Updating positions arrays.\n");
379  }
380  size_t n_positions = 0;
384 
385  for ( ; a_i != a_e; ++a_i, ++p_i ) {
386  atoms_positions[atoms_map[*a_i]] = cvm::rvector((*p_i).x, (*p_i).y, (*p_i).z);
387  n_positions++;
388  }
389 
390  // The following had to be relaxed because some atoms may be forced without their position being requested
391  // if (n_positions < atoms_ids.size()) {
392  // cvm::error("Error: did not receive the positions of all atoms.\n", BUG_ERROR);
393  // }
394  }
395 
396  if (total_force_requested && cvm::step_relative() > 0) {
397 
398  // sort the force arrays the previous step
399  // (but only do so if there *is* a previous step!)
400 
401  {
402  if (cvm::debug()) {
403  cvm::log("Updating total forces arrays.\n");
404  }
405  size_t n_total_forces = 0;
409 
410  for ( ; a_i != a_e; ++a_i, ++f_i ) {
411  atoms_total_forces[atoms_map[*a_i]] = cvm::rvector((*f_i).x, (*f_i).y, (*f_i).z);
412  n_total_forces++;
413  }
414 
415  if (n_total_forces < atoms_ids.size()) {
416  cvm::error("Error: total forces were requested, but total forces "
417  "were not received for all atoms.\n"
418  "The most probable cause is combination of energy "
419  "minimization with a biasing method that requires MD (e.g. ABF).\n"
420  "Always run minimization and ABF separately.", INPUT_ERROR);
421  }
422  }
423 
424  {
425  if (cvm::debug()) {
426  cvm::log("Updating group total forces arrays.\n");
427  }
430  size_t i = 0;
431  if ((f_e - f_i) != ((int) atom_groups_ids.size())) {
432  cvm::error("Error: total forces were requested for scalable groups, "
433  "but they are not in the same number from the number of groups.\n"
434  "The most probable cause is combination of energy "
435  "minimization with a biasing method that requires MD (e.g. ABF).\n"
436  "Always run minimization and ABF separately.", INPUT_ERROR);
437  }
438  for ( ; f_i != f_e; f_i++, i++) {
439  atom_groups_total_forces[i] = cvm::rvector((*f_i).x, (*f_i).y, (*f_i).z);
440  }
441  }
442  }
443 
444  {
445  if (cvm::debug()) {
446  cvm::log("Updating group positions arrays.\n");
447  }
448  // update group data (only coms available so far)
449  size_t ig;
450  // note: getGroupMassBegin() could be used here, but masses and charges
451  // have already been calculated from the last call to setup()
453  for (ig = 0; gp_i != getGroupPositionEnd(); gp_i++, ig++) {
454  atom_groups_coms[ig] = cvm::rvector(gp_i->x, gp_i->y, gp_i->z);
455  }
456  }
457 
458 #if NAMD_VERSION_NUMBER >= 34471681
459  {
460  if (cvm::debug()) {
461  log("Updating grid objects.\n");
462  }
463  // Using a simple nested loop: there probably won't be so many maps that
464  // this becomes performance-limiting
467  for ( ; gov_i != getGridObjValueEnd(); goi_i++, gov_i++) {
468  for (size_t imap = 0; imap < volmaps_ids.size(); imap++) {
469  if (volmaps_ids[imap] == *goi_i) {
470  volmaps_values[imap] = *gov_i;
471  break;
472  }
473  }
474  }
475  }
476 #endif
477 
478  if (cvm::debug()) {
479  cvm::log("atoms_ids = "+cvm::to_str(atoms_ids)+"\n");
480  cvm::log("atoms_ncopies = "+cvm::to_str(atoms_ncopies)+"\n");
481  cvm::log("atoms_masses = "+cvm::to_str(atoms_masses)+"\n");
482  cvm::log("atoms_charges = "+cvm::to_str(atoms_charges)+"\n");
483  cvm::log("atoms_positions = "+cvm::to_str(atoms_positions)+"\n");
484  cvm::log("atoms_total_forces = "+cvm::to_str(atoms_total_forces)+"\n");
485  cvm::log(cvm::line_marker);
486 
487  cvm::log("atom_groups_ids = "+cvm::to_str(atom_groups_ids)+"\n");
488  cvm::log("atom_groups_ncopies = "+cvm::to_str(atom_groups_ncopies)+"\n");
489  cvm::log("atom_groups_masses = "+cvm::to_str(atom_groups_masses)+"\n");
490  cvm::log("atom_groups_charges = "+cvm::to_str(atom_groups_charges)+"\n");
491  cvm::log("atom_groups_coms = "+cvm::to_str(atom_groups_coms)+"\n");
492  cvm::log("atom_groups_total_forces = "+cvm::to_str(atom_groups_total_forces)+"\n");
493  cvm::log(cvm::line_marker);
494 
495 #if NAMD_VERSION_NUMBER >= 34471681
496  cvm::log("volmaps_ids = "+cvm::to_str(volmaps_ids)+"\n");
497  cvm::log("volmaps_values = "+cvm::to_str(volmaps_values)+"\n");
498  cvm::log(cvm::line_marker);
499 #endif
500  }
501 
502  // call the collective variable module
503  if (colvars->calc() != COLVARS_OK) {
504  cvm::error("Error in the collective variables module.\n", COLVARS_ERROR);
505  }
506 
507  if (cvm::debug()) {
508  cvm::log(cvm::line_marker);
509  cvm::log("atoms_new_colvar_forces = "+cvm::to_str(atoms_new_colvar_forces)+"\n");
510  cvm::log(cvm::line_marker);
511  cvm::log("atom_groups_new_colvar_forces = "+cvm::to_str(atom_groups_new_colvar_forces)+"\n");
512  cvm::log(cvm::line_marker);
513 #if NAMD_VERSION_NUMBER >= 34471681
514  cvm::log("volmaps_new_colvar_forces = "+cvm::to_str(volmaps_new_colvar_forces)+"\n");
515  cvm::log(cvm::line_marker);
516 #endif
517  }
518 
519  // communicate all forces to the MD integrator
520  for (size_t i = 0; i < atoms_ids.size(); i++) {
521  cvm::rvector const &f = atoms_new_colvar_forces[i];
522  modifyForcedAtoms().add(atoms_ids[i]);
523  modifyAppliedForces().add(Vector(f.x, f.y, f.z));
524  }
525 
526  if (atom_groups_new_colvar_forces.size() > 0) {
529  for (int ig = 0; gf_i != modifyGroupForces().end(); gf_i++, ig++) {
530  cvm::rvector const &f = atom_groups_new_colvar_forces[ig];
531  *gf_i = Vector(f.x, f.y, f.z);
532  }
533  }
534 
535 #if NAMD_VERSION_NUMBER >= 34471681
536  if (volmaps_new_colvar_forces.size() > 0) {
541  for ( ; goi_i != getGridObjIndexEnd(); goi_i++, gof_i++) {
542  for (size_t imap = 0; imap < volmaps_ids.size(); imap++) {
543  if (volmaps_ids[imap] == *goi_i) {
544  *gof_i = volmaps_new_colvar_forces[imap];
545  break;
546  }
547  }
548  }
549  }
550 #endif
551 
552  // send MISC energy
553  reduction->submit();
554 
555  // NAMD does not destruct GlobalMaster objects, so we must remember
556  // to write all output files at the end of a run
557  if (step == simparams->N) {
558  post_run();
559  }
560 }
static Node * Object()
Definition: Node.h:86
ForceList & modifyAppliedForces()
Definition: GlobalMaster.C:162
BigRealList & modifyGridObjForces()
Definition: GlobalMaster.C:179
SimParameters * simparams
Pointer to the NAMD simulation input object.
const Elem * const_iterator
Definition: ResizeArray.h:38
AtomIDList::const_iterator getForceIdEnd()
Definition: GlobalMaster.C:260
AtomIDList::const_iterator getAtomIdBegin()
Definition: GlobalMaster.C:190
PositionList::const_iterator getGroupPositionEnd()
Definition: GlobalMaster.C:206
Definition: Vector.h:64
BigReal z
Definition: Vector.h:66
ForceList::const_iterator getGroupTotalForceBegin()
Definition: GlobalMaster.C:210
int orthogonal() const
Definition: Lattice.h:257
static int debug
Definition: parm.C:36
PositionList::const_iterator getAtomPositionBegin()
Definition: GlobalMaster.C:198
const ResizeArray< AtomIDList > & requestedGroups()
Definition: GlobalMaster.C:149
AtomIDList::const_iterator getForceIdBegin()
Definition: GlobalMaster.C:255
std::vector< int > atoms_map
Array of atom indices (relative to the colvarproxy arrays), usedfor faster copy of atomic data...
IntList::const_iterator getGridObjIndexBegin()
Definition: GlobalMaster.C:218
iterator end(void)
Definition: ResizeArray.h:37
void clear()
Definition: ResizeArray.h:87
void setall(const Elem &elem)
Definition: ResizeArray.h:90
BigReal x
Definition: Vector.h:66
void log(std::string const &message)
AtomIDList & modifyForcedAtoms()
Definition: GlobalMaster.C:157
int numAtoms
Definition: Molecule.h:556
ForceList & modifyGroupForces()
Definition: GlobalMaster.C:167
Bool set(const char *s)
Definition: Vector.h:228
BigRealList::const_iterator getGridObjValueBegin()
Definition: GlobalMaster.C:226
int add(const Elem &elem)
Definition: ResizeArray.h:97
IntList::const_iterator getGridObjIndexEnd()
Definition: GlobalMaster.C:222
const Lattice * lattice
Definition: GlobalMaster.h:141
void resize(int i)
Definition: ResizeArray.h:84
AtomIDList::const_iterator getAtomIdEnd()
Definition: GlobalMaster.C:194
BigReal y
Definition: Vector.h:66
Vector b() const
Definition: Lattice.h:253
void submit(void)
Definition: ReductionMgr.h:323
const IntList & requestedGridObjs()
Definition: GlobalMaster.C:153
SubmitReduction * reduction
Used to submit restraint energy as MISC.
int b_p() const
Definition: Lattice.h:274
int update_atoms_map(AtomIDList::const_iterator begin, AtomIDList::const_iterator end)
PositionList::const_iterator getGroupPositionBegin()
Definition: GlobalMaster.C:202
int a_p() const
Definition: Lattice.h:273
Molecule * molecule
Definition: Node.h:176
Vector a() const
Definition: Lattice.h:252
ForceList::const_iterator getTotalForce()
Definition: GlobalMaster.C:265
Vector c() const
Definition: Lattice.h:254
BigRealList::const_iterator getGridObjValueEnd()
Definition: GlobalMaster.C:230
int c_p() const
Definition: Lattice.h:275
ForceList::const_iterator getGroupTotalForceEnd()
Definition: GlobalMaster.C:214
iterator begin(void)
Definition: ResizeArray.h:36
int colvarproxy_namd::check_atom_id ( int  atom_number)

Definition at line 652 of file colvarproxy_namd.C.

References debug, Node::molecule, Molecule::numAtoms, and Node::Object().

Referenced by init_atom().

653 {
654  // NAMD's internal numbering starts from zero
655  int const aid = (atom_number-1);
656 
657  if (cvm::debug())
658  cvm::log("Adding atom "+cvm::to_str(atom_number)+
659  " for collective variables calculation.\n");
660 
661  if ( (aid < 0) || (aid >= Node::Object()->molecule->numAtoms) ) {
662  cvm::error("Error: invalid atom number specified, "+
663  cvm::to_str(atom_number)+"\n", INPUT_ERROR);
664  return INPUT_ERROR;
665  }
666 
667  return aid;
668 }
static Node * Object()
Definition: Node.h:86
static int debug
Definition: parm.C:36
int numAtoms
Definition: Molecule.h:556
Molecule * molecule
Definition: Node.h:176
int colvarproxy_namd::check_atom_id ( cvm::residue_id const &  residue,
std::string const &  atom_name,
std::string const &  segment_id 
)

Definition at line 698 of file colvarproxy_namd.C.

References Molecule::get_atom_from_name(), Node::molecule, and Node::Object().

701 {
702  int const aid =
703  (segment_id.size() ?
704  Node::Object()->molecule->get_atom_from_name(segment_id.c_str(),
705  residue,
706  atom_name.c_str()) :
708  residue,
709  atom_name.c_str()));
710 
711  if (aid < 0) {
712  // get_atom_from_name() has returned an error value
713  cvm::error("Error: could not find atom \""+
714  atom_name+"\" in residue "+
715  cvm::to_str(residue)+
716  ( (segment_id != "MAIN") ?
717  (", segment \""+segment_id+"\"") :
718  ("") )+
719  "\n", INPUT_ERROR);
720  return INPUT_ERROR;
721  }
722 
723  return aid;
724 }
static Node * Object()
Definition: Node.h:86
int get_atom_from_name(const char *segid, int resid, const char *aname) const
Definition: Molecule.C:121
Molecule * molecule
Definition: Node.h:176
void colvarproxy_namd::clear_atom ( int  index)

Definition at line 759 of file colvarproxy_namd.C.

760 {
761  colvarproxy::clear_atom(index);
762  // TODO remove it from GlobalMaster arrays?
763 }
void colvarproxy_namd::clear_atom_group ( int  index)

Definition at line 1194 of file colvarproxy_namd.C.

1195 {
1196  // do nothing, keep the NAMD arrays in sync with the colvarproxy ones
1197  colvarproxy::clear_atom_group(index);
1198 }
void colvarproxy_namd::clear_volmap ( int  index)

Definition at line 1290 of file colvarproxy_namd.C.

1291 {
1292  colvarproxy::clear_volmap(index);
1293 }
int colvarproxy_namd::close_output_stream ( std::string const &  output_name)

Definition at line 1057 of file colvarproxy_namd.C.

1058 {
1059  std::list<std::ostream *>::iterator osi = output_files.begin();
1060  std::list<std::string>::iterator osni = output_stream_names.begin();
1061  for ( ; osi != output_files.end(); osi++, osni++) {
1062  if (*osni == output_name) {
1063  if (((ofstream_namd *) *osi)->is_open()) {
1064  ((ofstream_namd *) *osi)->close();
1065  }
1066  delete *osi;
1067  output_files.erase(osi);
1068  output_stream_names.erase(osni);
1069  return COLVARS_OK;
1070  }
1071  }
1072  return COLVARS_ERROR;
1073 }
cvm::real colvarproxy_namd::dt ( )
inline
void colvarproxy_namd::error ( std::string const &  message)

Definition at line 624 of file colvarproxy_namd.C.

References log(), NAMD_die(), and NAMD_err().

Referenced by colvarproxy_namd().

625 {
626  log(message);
627  switch (cvm::get_error()) {
628  case FILE_ERROR:
629  errno = EIO; break;
630  case COLVARS_NOT_IMPLEMENTED:
631  errno = ENOSYS; break;
632  case MEMORY_ERROR:
633  errno = ENOMEM; break;
634  }
635  char const *msg = "Error in the collective variables module "
636  "(see above for details)";
637  if (errno) {
638  NAMD_err(msg);
639  } else {
640  NAMD_die(msg);
641  }
642 }
void NAMD_err(const char *err_msg)
Definition: common.C:102
void log(std::string const &message)
void NAMD_die(const char *err_msg)
Definition: common.C:83
void colvarproxy_namd::exit ( std::string const &  message)

Definition at line 645 of file colvarproxy_namd.C.

References BackEnd::exit(), and log().

646 {
647  log(message);
648  BackEnd::exit();
649 }
static void exit(int status=0)
Definition: BackEnd.C:276
void log(std::string const &message)
int colvarproxy_namd::flush_output_stream ( std::ostream *  os)

Definition at line 1043 of file colvarproxy_namd.C.

1044 {
1045  std::list<std::ostream *>::iterator osi = output_files.begin();
1046  std::list<std::string>::iterator osni = output_stream_names.begin();
1047  for ( ; osi != output_files.end(); osi++, osni++) {
1048  if (*osi == os) {
1049  ((ofstream_namd *) *osi)->flush();
1050  return COLVARS_OK;
1051  }
1052  }
1053  return COLVARS_ERROR;
1054 }
int colvarproxy_namd::init_atom ( int  atom_number)

Definition at line 671 of file colvarproxy_namd.C.

References ResizeArray< T >::add(), check_atom_id(), GlobalMaster::modifyRequestedAtoms(), and update_atom_properties().

672 {
673  // save time by checking first whether this atom has been requested before
674  // (this is more common than a non-valid atom number)
675  int aid = (atom_number-1);
676 
677  for (size_t i = 0; i < atoms_ids.size(); i++) {
678  if (atoms_ids[i] == aid) {
679  // this atom id was already recorded
680  atoms_ncopies[i] += 1;
681  return i;
682  }
683  }
684 
685  aid = check_atom_id(atom_number);
686 
687  if (aid < 0) {
688  return INPUT_ERROR;
689  }
690 
691  int const index = add_atom_slot(aid);
692  modifyRequestedAtoms().add(aid);
693  update_atom_properties(index);
694  return index;
695 }
AtomIDList & modifyRequestedAtoms()
Definition: GlobalMaster.C:127
void update_atom_properties(int index)
int add(const Elem &elem)
Definition: ResizeArray.h:97
int check_atom_id(int atom_number)
int colvarproxy_namd::init_atom ( cvm::residue_id const &  residue,
std::string const &  atom_name,
std::string const &  segment_id 
)

For AMBER topologies, the segment id is automatically set to "MAIN" (the segment id assigned by NAMD's AMBER topology parser), and is therefore optional when an AMBER topology is used

Definition at line 731 of file colvarproxy_namd.C.

References ResizeArray< T >::add(), check_atom_id(), debug, GlobalMaster::modifyRequestedAtoms(), and update_atom_properties().

734 {
735  int const aid = check_atom_id(residue, atom_name, segment_id);
736 
737  for (size_t i = 0; i < atoms_ids.size(); i++) {
738  if (atoms_ids[i] == aid) {
739  // this atom id was already recorded
740  atoms_ncopies[i] += 1;
741  return i;
742  }
743  }
744 
745  if (cvm::debug())
746  cvm::log("Adding atom \""+
747  atom_name+"\" in residue "+
748  cvm::to_str(residue)+
749  " (index "+cvm::to_str(aid)+
750  ") for collective variables calculation.\n");
751 
752  int const index = add_atom_slot(aid);
753  modifyRequestedAtoms().add(aid);
754  update_atom_properties(index);
755  return index;
756 }
AtomIDList & modifyRequestedAtoms()
Definition: GlobalMaster.C:127
static int debug
Definition: parm.C:36
void update_atom_properties(int index)
int add(const Elem &elem)
Definition: ResizeArray.h:97
int check_atom_id(int atom_number)
int colvarproxy_namd::init_atom_group ( std::vector< int > const &  atoms_ids)

Definition at line 1124 of file colvarproxy_namd.C.

References debug, GlobalMaster::modifyGroupForces(), GlobalMaster::modifyRequestedGroups(), Node::molecule, Molecule::numAtoms, Node::Object(), ResizeArray< T >::resize(), ResizeArray< T >::size(), and update_group_properties().

1125 {
1126  if (cvm::debug())
1127  cvm::log("Reguesting from NAMD a group of size "+cvm::to_str(atoms_ids.size())+
1128  " for collective variables calculation.\n");
1129 
1130  // Note: modifyRequestedGroups is supposed to be in sync with the colvarproxy arrays,
1131  // and to stay that way during a simulation
1132 
1133  // compare this new group to those already allocated inside GlobalMaster
1134  int ig;
1135  for (ig = 0; ig < modifyRequestedGroups().size(); ig++) {
1136  AtomIDList const &namd_group = modifyRequestedGroups()[ig];
1137  bool b_match = true;
1138 
1139  if (namd_group.size() != ((int) atoms_ids.size())) {
1140  b_match = false;
1141  } else {
1142  int ia;
1143  for (ia = 0; ia < namd_group.size(); ia++) {
1144  int const aid = atoms_ids[ia];
1145  if (namd_group[ia] != aid) {
1146  b_match = false;
1147  break;
1148  }
1149  }
1150  }
1151 
1152  if (b_match) {
1153  if (cvm::debug())
1154  cvm::log("Group was already added.\n");
1155  // this group already exists
1156  atom_groups_ncopies[ig] += 1;
1157  return ig;
1158  }
1159  }
1160 
1161  // add this group (note: the argument of add_atom_group_slot() is redundant for NAMD, and provided only for consistency)
1162  size_t const index = add_atom_group_slot(atom_groups_ids.size());
1163  modifyRequestedGroups().resize(atom_groups_ids.size());
1164  // the following is done in calculate()
1165  // modifyGroupForces().resize(atom_groups_ids.size());
1166  AtomIDList &namd_group = modifyRequestedGroups()[index];
1167  namd_group.resize(atoms_ids.size());
1168  int const n_all_atoms = Node::Object()->molecule->numAtoms;
1169  for (size_t ia = 0; ia < atoms_ids.size(); ia++) {
1170  int const aid = atoms_ids[ia];
1171  if (cvm::debug())
1172  cvm::log("Adding atom "+cvm::to_str(aid+1)+
1173  " for collective variables calculation.\n");
1174  if ( (aid < 0) || (aid >= n_all_atoms) ) {
1175  cvm::error("Error: invalid atom number specified, "+
1176  cvm::to_str(aid+1)+"\n", INPUT_ERROR);
1177  return -1;
1178  }
1179  namd_group[ia] = aid;
1180  }
1181 
1182  update_group_properties(index);
1183 
1184  if (cvm::debug()) {
1185  cvm::log("Group has index "+cvm::to_str(index)+"\n");
1186  cvm::log("modifyRequestedGroups length = "+cvm::to_str(modifyRequestedGroups().size())+
1187  ", modifyGroupForces length = "+cvm::to_str(modifyGroupForces().size())+"\n");
1188  }
1189 
1190  return index;
1191 }
static Node * Object()
Definition: Node.h:86
static int debug
Definition: parm.C:36
int numAtoms
Definition: Molecule.h:556
ResizeArray< AtomIDList > & modifyRequestedGroups()
Definition: GlobalMaster.C:184
int update_group_properties(int index)
ForceList & modifyGroupForces()
Definition: GlobalMaster.C:167
void resize(int i)
Definition: ResizeArray.h:84
int size(void) const
Definition: ResizeArray.h:127
Molecule * molecule
Definition: Node.h:176
void colvarproxy_namd::init_tcl_pointers ( )
protected

Definition at line 565 of file colvarproxy_namd.C.

References Node::getScript(), and Node::Object().

Referenced by colvarproxy_namd().

566 {
567 #ifdef NAMD_TCL
568  // Store pointer to NAMD's Tcl interpreter
569  tcl_interp_ = reinterpret_cast<void *>(Node::Object()->getScript()->interp);
570 #endif
571 }
static Node * Object()
Definition: Node.h:86
ScriptTcl * getScript(void)
Definition: Node.h:192
int colvarproxy_namd::init_volmap ( int  volmap_id)

Definition at line 1240 of file colvarproxy_namd.C.

References ResizeArray< T >::add(), GlobalMaster::modifyRequestedGridObjects(), Node::molecule, Molecule::numGridforceGrids, and Node::Object().

Referenced by init_volmap().

1241 {
1242  for (size_t i = 0; i < volmaps_ids.size(); i++) {
1243  if (volmaps_ids[i] == volmap_id) {
1244  // this map has already been requested
1245  volmaps_ncopies[i] += 1;
1246  return i;
1247  }
1248  }
1249 
1250  Molecule *mol = Node::Object()->molecule;
1251 
1252  if ((volmap_id < 0) || (volmap_id >= mol->numGridforceGrids)) {
1253  return cvm::error("Error: invalid numeric ID ("+cvm::to_str(volmap_id)+
1254  ") for map.\n", INPUT_ERROR);
1255  }
1256 
1257  int const index = add_volmap_slot(volmap_id);
1258  modifyRequestedGridObjects().add(volmap_id);
1259 
1260  return index;
1261 }
static Node * Object()
Definition: Node.h:86
int numGridforceGrids
Definition: Molecule.h:590
IntList & modifyRequestedGridObjects()
Definition: GlobalMaster.C:173
int add(const Elem &elem)
Definition: ResizeArray.h:97
Molecule * molecule
Definition: Node.h:176
int colvarproxy_namd::init_volmap ( const char *  volmap_name)

Definition at line 1264 of file colvarproxy_namd.C.

References Molecule::get_gridfrc_grid(), GridforceGrid::get_scale(), MGridforceParamsList::index_for_key(), init_volmap(), SimParameters::mgridforcelist, Node::molecule, Node::Object(), simparams, Vector::x, Vector::y, and Vector::z.

1265 {
1266  if (volmap_name == NULL) {
1267  return cvm::error("Error: no grid object name provided.", INPUT_ERROR);
1268  }
1269  int volmap_id = simparams->mgridforcelist.index_for_key(volmap_name);
1270 
1271  int error_code = init_volmap(volmap_id);
1272 
1273  if (error_code == COLVARS_OK) {
1274  // Check that the scale factor is correctly set to zero
1275  Molecule *mol = Node::Object()->molecule;
1276  GridforceGrid const *grid = mol->get_gridfrc_grid(volmap_id);
1277  Vector const gfScale = grid->get_scale();
1278  if ((gfScale.x != 0.0) || (gfScale.y != 0.0) || (gfScale.z != 0.0)) {
1279  error_code |= cvm::error("Error: GridForce map \""+
1280  std::string(volmap_name)+
1281  "\" has non-zero scale factors.\n",
1282  INPUT_ERROR);
1283  }
1284  }
1285 
1286  return error_code;
1287 }
static Node * Object()
Definition: Node.h:86
SimParameters * simparams
Pointer to the NAMD simulation input object.
Definition: Vector.h:64
BigReal z
Definition: Vector.h:66
int index_for_key(const char *key)
BigReal x
Definition: Vector.h:66
GridforceGrid * get_gridfrc_grid(int gridnum) const
Definition: Molecule.h:1276
MGridforceParamsList mgridforcelist
virtual Vector get_scale(void) const =0
BigReal y
Definition: Vector.h:66
int init_volmap(int volmap_id)
Molecule * molecule
Definition: Node.h:176
int colvarproxy_namd::load_atoms ( char const *  filename,
cvm::atom_group &  atoms,
std::string const &  pdb_field,
double const  pdb_field_value = 0.0 
)

Definition at line 961 of file colvarproxy_namd.C.

References PDB::atom(), cvm::atom, e_pdb_beta, e_pdb_occ, e_pdb_x, e_pdb_y, e_pdb_z, PDB::num_atoms(), and pdb_field_str2enum().

965 {
966  if (pdb_field_str.size() == 0)
967  cvm::error("Error: must define which PDB field to use "
968  "in order to define atoms from a PDB file.\n", INPUT_ERROR);
969 
970  PDB *pdb = new PDB(pdb_filename);
971  size_t const pdb_natoms = pdb->num_atoms();
972 
973  e_pdb_field pdb_field_index = pdb_field_str2enum(pdb_field_str);
974 
975  for (size_t ipdb = 0; ipdb < pdb_natoms; ipdb++) {
976 
977  double atom_pdb_field_value = 0.0;
978 
979  switch (pdb_field_index) {
980  case e_pdb_occ:
981  atom_pdb_field_value = (pdb->atom(ipdb))->occupancy();
982  break;
983  case e_pdb_beta:
984  atom_pdb_field_value = (pdb->atom(ipdb))->temperaturefactor();
985  break;
986  case e_pdb_x:
987  atom_pdb_field_value = (pdb->atom(ipdb))->xcoor();
988  break;
989  case e_pdb_y:
990  atom_pdb_field_value = (pdb->atom(ipdb))->ycoor();
991  break;
992  case e_pdb_z:
993  atom_pdb_field_value = (pdb->atom(ipdb))->zcoor();
994  break;
995  default:
996  break;
997  }
998 
999  if ( (pdb_field_value) &&
1000  (atom_pdb_field_value != pdb_field_value) ) {
1001  continue;
1002  } else if (atom_pdb_field_value == 0.0) {
1003  continue;
1004  }
1005 
1006  if (atoms.is_enabled(colvardeps::f_ag_scalable)) {
1007  atoms.add_atom_id(ipdb);
1008  } else {
1009  atoms.add_atom(cvm::atom(ipdb+1));
1010  }
1011  }
1012 
1013  delete pdb;
1014  return (cvm::get_error() ? COLVARS_ERROR : COLVARS_OK);
1015 }
Definition: PDB.h:35
static __thread atom * atoms
e_pdb_field pdb_field_str2enum(std::string const &pdb_field_str)
int num_atoms(void)
Definition: PDB.C:323
e_pdb_field
PDBAtom * atom(int place)
Definition: PDB.C:394
int colvarproxy_namd::load_coords ( char const *  filename,
std::vector< cvm::atom_pos > &  pos,
const std::vector< int > &  indices,
std::string const &  pdb_field,
double const  pdb_field_value = 0.0 
)

Definition at line 849 of file colvarproxy_namd.C.

References PDB::atom(), e_pdb_beta, e_pdb_occ, e_pdb_x, e_pdb_y, e_pdb_z, PDB::num_atoms(), and pdb_field_str2enum().

854 {
855  if (pdb_field_str.size() == 0 && indices.size() == 0) {
856  cvm::error("Bug alert: either PDB field should be defined or list of "
857  "atom IDs should be available when loading atom coordinates!\n", BUG_ERROR);
858  }
859 
860  e_pdb_field pdb_field_index;
861  bool const use_pdb_field = (pdb_field_str.size() > 0);
862  if (use_pdb_field) {
863  pdb_field_index = pdb_field_str2enum(pdb_field_str);
864  }
865 
866  // next index to be looked up in PDB file (if list is supplied)
867  std::vector<int>::const_iterator current_index = indices.begin();
868 
869  PDB *pdb = new PDB(pdb_filename);
870  size_t const pdb_natoms = pdb->num_atoms();
871 
872  if (pos.size() != pdb_natoms) {
873 
874  bool const pos_allocated = (pos.size() > 0);
875 
876  size_t ipos = 0, ipdb = 0;
877  for ( ; ipdb < pdb_natoms; ipdb++) {
878 
879  if (use_pdb_field) {
880  // PDB field mode: skip atoms with wrong value in PDB field
881  double atom_pdb_field_value = 0.0;
882 
883  switch (pdb_field_index) {
884  case e_pdb_occ:
885  atom_pdb_field_value = (pdb->atom(ipdb))->occupancy();
886  break;
887  case e_pdb_beta:
888  atom_pdb_field_value = (pdb->atom(ipdb))->temperaturefactor();
889  break;
890  case e_pdb_x:
891  atom_pdb_field_value = (pdb->atom(ipdb))->xcoor();
892  break;
893  case e_pdb_y:
894  atom_pdb_field_value = (pdb->atom(ipdb))->ycoor();
895  break;
896  case e_pdb_z:
897  atom_pdb_field_value = (pdb->atom(ipdb))->zcoor();
898  break;
899  default:
900  break;
901  }
902 
903  if ( (pdb_field_value) &&
904  (atom_pdb_field_value != pdb_field_value) ) {
905  continue;
906  } else if (atom_pdb_field_value == 0.0) {
907  continue;
908  }
909 
910  } else {
911  // Atom ID mode: use predefined atom IDs from the atom group
912  if (((int) ipdb) != *current_index) {
913  // Skip atoms not in the list
914  continue;
915  } else {
916  current_index++;
917  }
918  }
919 
920  if (!pos_allocated) {
921  pos.push_back(cvm::atom_pos(0.0, 0.0, 0.0));
922  } else if (ipos >= pos.size()) {
923  cvm::error("Error: the PDB file \""+
924  std::string(pdb_filename)+
925  "\" contains coordinates for "
926  "more atoms than needed.\n", BUG_ERROR);
927  }
928 
929  pos[ipos] = cvm::atom_pos((pdb->atom(ipdb))->xcoor(),
930  (pdb->atom(ipdb))->ycoor(),
931  (pdb->atom(ipdb))->zcoor());
932  ipos++;
933  if (!use_pdb_field && current_index == indices.end())
934  break;
935  }
936 
937  if (ipos < pos.size() || (!use_pdb_field && current_index != indices.end())) {
938  size_t n_requested = use_pdb_field ? pos.size() : indices.size();
939  cvm::error("Error: number of matching records in the PDB file \""+
940  std::string(pdb_filename)+"\" ("+cvm::to_str(ipos)+
941  ") does not match the number of requested coordinates ("+
942  cvm::to_str(n_requested)+").\n", INPUT_ERROR);
943  return COLVARS_ERROR;
944  }
945  } else {
946 
947  // when the PDB contains exactly the number of atoms of the array,
948  // ignore the fields and just read coordinates
949  for (size_t ia = 0; ia < pos.size(); ia++) {
950  pos[ia] = cvm::atom_pos((pdb->atom(ia))->xcoor(),
951  (pdb->atom(ia))->ycoor(),
952  (pdb->atom(ia))->zcoor());
953  }
954  }
955 
956  delete pdb;
957  return COLVARS_OK;
958 }
Definition: PDB.h:35
e_pdb_field pdb_field_str2enum(std::string const &pdb_field_str)
int num_atoms(void)
Definition: PDB.C:323
e_pdb_field
PDBAtom * atom(int place)
Definition: PDB.C:394
void colvarproxy_namd::log ( std::string const &  message)

Definition at line 614 of file colvarproxy_namd.C.

References endi(), and iout.

Referenced by calculate(), error(), exit(), setup(), and update_atom_properties().

615 {
616  std::istringstream is(message);
617  std::string line;
618  while (std::getline(is, line))
619  iout << "colvars: " << line << "\n";
620  iout << endi;
621 }
#define iout
Definition: InfoStream.h:87
infostream & endi(infostream &s)
Definition: InfoStream.C:38
int colvarproxy_namd::num_replicas ( )
virtual

Definition at line 1399 of file colvarproxy_namd.C.

1399  {
1400  return CmiNumPartitions();
1401 }
std::ostream * colvarproxy_namd::output_stream ( std::string const &  output_name,
std::ios_base::openmode  mode 
)

Definition at line 1018 of file colvarproxy_namd.C.

References debug, and ofstream_namd::is_open().

1020 {
1021  if (cvm::debug()) {
1022  cvm::log("Using colvarproxy_namd::output_stream()\n");
1023  }
1024 
1025  std::ostream *os = get_output_stream(output_name);
1026  if (os != NULL) return os;
1027 
1028  if (!(mode & (std::ios_base::app | std::ios_base::ate))) {
1029  colvarproxy::backup_file(output_name);
1030  }
1031  ofstream_namd *osf = new ofstream_namd(output_name.c_str(), mode);
1032  if (!osf->is_open()) {
1033  cvm::error("Error: cannot write to file \""+output_name+"\".\n",
1034  FILE_ERROR);
1035  return NULL;
1036  }
1037  output_stream_names.push_back(output_name);
1038  output_files.push_back(osf);
1039  return osf;
1040 }
static int debug
Definition: parm.C:36
bool is_open() const
Definition: fstream_namd.h:30
cvm::rvector colvarproxy_namd::position_distance ( cvm::atom_pos const &  pos1,
cvm::atom_pos const &  pos2 
) const

Definition at line 782 of file colvarproxy_namd.C.

References Lattice::delta(), GlobalMaster::lattice, Vector::x, Vector::y, and Vector::z.

785 {
786  Position const p1(pos1.x, pos1.y, pos1.z);
787  Position const p2(pos2.x, pos2.y, pos2.z);
788  // return p2 - p1
789  if (this->lattice != NULL) {
790  Vector const d = this->lattice->delta(p2, p1);
791  return cvm::rvector(d.x, d.y, d.z);
792  }
793  else {
794  return colvarproxy_system::position_distance(pos1, pos2);
795  }
796 }
Definition: Vector.h:64
BigReal z
Definition: Vector.h:66
Vector delta(const Position &pos1, const Position &pos2) const
Definition: Lattice.h:144
BigReal x
Definition: Vector.h:66
const Lattice * lattice
Definition: GlobalMaster.h:141
BigReal y
Definition: Vector.h:66
cvm::real colvarproxy_namd::rand_gaussian ( )
inline

Definition at line 110 of file colvarproxy_namd.h.

References Random::gaussian(), and random.

111  {
112  return random.gaussian();
113  }
Random random
NAMD-style PRNG object.
BigReal gaussian(void)
Definition: Random.h:116
void colvarproxy_namd::replica_comm_barrier ( )
virtual

Definition at line 1404 of file colvarproxy_namd.C.

References replica_barrier().

1404  {
1405  replica_barrier();
1406 }
void replica_barrier()
int colvarproxy_namd::replica_comm_recv ( char *  msg_data,
int  buf_len,
int  src_rep 
)
virtual

Definition at line 1409 of file colvarproxy_namd.C.

References DataMessage::data, replica_recv(), and DataMessage::size.

1410  {
1411  DataMessage *recvMsg = NULL;
1412  replica_recv(&recvMsg, src_rep, CkMyPe());
1413  CmiAssert(recvMsg != NULL);
1414  int retval = recvMsg->size;
1415  if (buf_len >= retval) {
1416  memcpy(msg_data,recvMsg->data,retval);
1417  } else {
1418  retval = 0;
1419  }
1420  CmiFree(recvMsg);
1421  return retval;
1422 }
void replica_recv(DataMessage **precvMsg, int srcPart, int srcPE)
char data[1]
Definition: DataExchanger.h:23
int colvarproxy_namd::replica_comm_send ( char *  msg_data,
int  msg_len,
int  dest_rep 
)
virtual

Definition at line 1425 of file colvarproxy_namd.C.

References replica_send().

1426  {
1427  replica_send(msg_data, msg_len, dest_rep, CkMyPe());
1428  return msg_len;
1429 }
void replica_send(const char *sndbuf, int sendcount, int destPart, int destPE)
int colvarproxy_namd::replica_enabled ( )
virtual

Definition at line 1385 of file colvarproxy_namd.C.

1385  {
1386 #if CMK_HAS_PARTITION
1387  return COLVARS_OK;
1388 #else
1389  return COLVARS_NOT_IMPLEMENTED;
1390 #endif
1391 }
int colvarproxy_namd::replica_index ( )
virtual

Definition at line 1394 of file colvarproxy_namd.C.

1394  {
1395  return CmiMyPartition();
1396 }
void colvarproxy_namd::request_total_force ( bool  yesno)

Definition at line 601 of file colvarproxy_namd.C.

References debug, and GlobalMaster::requestTotalForce().

602 {
603  if (cvm::debug()) {
604  cvm::log("colvarproxy_namd::request_total_force()\n");
605  }
606  total_force_requested = yesno;
607  requestTotalForce(total_force_requested);
608  if (cvm::debug()) {
609  cvm::log("colvarproxy_namd::request_total_force() end\n");
610  }
611 }
static int debug
Definition: parm.C:36
void requestTotalForce(bool yesno=true)
Definition: GlobalMaster.h:133
int colvarproxy_namd::reset ( )

Definition at line 258 of file colvarproxy_namd.C.

References atoms_map, ResizeArray< T >::clear(), GlobalMaster::modifyRequestedAtoms(), GlobalMaster::modifyRequestedGridObjects(), and GlobalMaster::modifyRequestedGroups().

259 {
260  int error_code = COLVARS_OK;
261 
262  // Unrequest all atoms and group from NAMD
264  modifyRequestedGroups().clear();
265 #if NAMD_VERSION_NUMBER >= 34471681
267 #endif
268 
269  atoms_map.clear();
270 
271  // Clear internal Proxy records
272  error_code |= colvarproxy::reset();
273 
274  return error_code;
275 }
AtomIDList & modifyRequestedAtoms()
Definition: GlobalMaster.C:127
std::vector< int > atoms_map
Array of atom indices (relative to the colvarproxy arrays), usedfor faster copy of atomic data...
IntList & modifyRequestedGridObjects()
Definition: GlobalMaster.C:173
void clear()
Definition: ResizeArray.h:87
ResizeArray< AtomIDList > & modifyRequestedGroups()
Definition: GlobalMaster.C:184
int colvarproxy_namd::run_colvar_callback ( std::string const &  name,
std::vector< const colvarvalue * > const &  cvcs,
colvarvalue &  value 
)

Definition at line 578 of file colvarproxy_namd.C.

582 {
583  return colvarproxy::tcl_run_colvar_callback(name, cvc_values, value);
584 }
int colvarproxy_namd::run_colvar_gradient_callback ( std::string const &  name,
std::vector< const colvarvalue * > const &  cvcs,
std::vector< cvm::matrix2d< cvm::real > > &  gradient 
)

Definition at line 586 of file colvarproxy_namd.C.

590 {
591  return colvarproxy::tcl_run_colvar_gradient_callback(name, cvc_values,
592  gradient);
593 }
int colvarproxy_namd::run_force_callback ( )

Definition at line 573 of file colvarproxy_namd.C.

574 {
575  return colvarproxy::tcl_run_force_callback();
576 }
int colvarproxy_namd::scalable_group_coms ( )
inline

Definition at line 208 of file colvarproxy_namd.h.

209  {
210  return COLVARS_OK;
211  }
char const * colvarproxy_namd::script_obj_to_str ( unsigned char *  obj)

Definition at line 1087 of file colvarproxy_namd.C.

1088 {
1089 #ifdef NAMD_TCL
1090  return colvarproxy_tcl::tcl_get_str(obj);
1091 #else
1092  // This is most likely not going to be executed
1093  return colvarproxy::script_obj_to_str(obj);
1094 #endif
1095 }
std::vector< std::string > colvarproxy_namd::script_obj_to_str_vector ( unsigned char *  obj)

Definition at line 1098 of file colvarproxy_namd.C.

References debug.

1099 {
1100  if (cvm::debug()) {
1101  cvm::log("Called colvarproxy_namd::script_obj_to_str_vector().\n");
1102  }
1103  std::vector<std::string> result;
1104 #ifdef NAMD_TCL
1105  Tcl_Interp *interp = reinterpret_cast<Tcl_Interp *>(tcl_interp_);
1106  Tcl_Obj *tcl_obj = reinterpret_cast<Tcl_Obj *>(obj);
1107  Tcl_Obj **tcl_list_elems = NULL;
1108  int count = 0;
1109  if (Tcl_ListObjGetElements(interp, tcl_obj, &count, &tcl_list_elems) ==
1110  TCL_OK) {
1111  result.reserve(count);
1112  for (int i = 0; i < count; i++) {
1113  result.push_back(Tcl_GetString(tcl_list_elems[i]));
1114  }
1115  } else {
1116  Tcl_SetResult(interp,
1117  const_cast<char *>("Cannot parse Tcl list."), TCL_STATIC);
1118  }
1119 #endif
1120  return result;
1121 }
static int debug
Definition: parm.C:36
int colvarproxy_namd::set_unit_system ( std::string const &  units_in,
bool  check_only 
)

Definition at line 1228 of file colvarproxy_namd.C.

1229 {
1230  if (units_in != "real") {
1231  cvm::error("Error: Specified unit system \"" + units_in + "\" is unsupported in NAMD. Supported units are \"real\" (A, kcal/mol).\n");
1232  return COLVARS_ERROR;
1233  }
1234  return COLVARS_OK;
1235 }
int colvarproxy_namd::setup ( )

Definition at line 200 of file colvarproxy_namd.C.

References log(), GlobalMaster::modifyGridObjForces(), GlobalMaster::modifyRequestedGroups(), simparams, ResizeArray< T >::size(), update_atom_properties(), update_group_properties(), and SimParameters::wrapAll.

Referenced by calculate().

201 {
202  if (colvars->size() == 0) return COLVARS_OK;
203 
204  cvm::log("Updating NAMD interface:\n");
205 
206  errno = 0;
207 
208  if (simparams->wrapAll) {
209  cvm::log("Warning: enabling wrapAll can lead to inconsistent results "
210  "for Colvars calculations: please disable wrapAll, "
211  "as is the default option in NAMD.\n");
212  }
213 
214  cvm::log("updating atomic data ("+cvm::to_str(atoms_ids.size())+" atoms).\n");
215 
216  size_t i;
217  for (i = 0; i < atoms_ids.size(); i++) {
219 
220  // zero out mutable arrays
221  atoms_positions[i] = cvm::rvector(0.0, 0.0, 0.0);
222  atoms_total_forces[i] = cvm::rvector(0.0, 0.0, 0.0);
223  atoms_new_colvar_forces[i] = cvm::rvector(0.0, 0.0, 0.0);
224  }
225 
226  size_t n_group_atoms = 0;
227  for (int ig = 0; ig < modifyRequestedGroups().size(); ig++) {
228  n_group_atoms += modifyRequestedGroups()[ig].size();
229  }
230 
231  cvm::log("updating group data ("+cvm::to_str(atom_groups_ids.size())+
232  " scalable groups, "+
233  cvm::to_str(n_group_atoms)+" atoms in total).\n");
234 
235  // Note: groupMassBegin, groupMassEnd may be used here, but they won't work for charges
236  for (int ig = 0; ig < modifyRequestedGroups().size(); ig++) {
237 
238  // update mass and charge
240 
241  atom_groups_coms[ig] = cvm::rvector(0.0, 0.0, 0.0);
242  atom_groups_total_forces[ig] = cvm::rvector(0.0, 0.0, 0.0);
243  atom_groups_new_colvar_forces[ig] = cvm::rvector(0.0, 0.0, 0.0);
244  }
245 
246 #if NAMD_VERSION_NUMBER >= 34471681
247  log("updating grid object data ("+cvm::to_str(volmaps_ids.size())+
248  " grid objects in total).\n");
249  for (int imap = 0; imap < modifyGridObjForces().size(); imap++) {
250  volmaps_new_colvar_forces[imap] = 0.0;
251  }
252 #endif
253 
254  return COLVARS_OK;
255 }
BigRealList & modifyGridObjForces()
Definition: GlobalMaster.C:179
SimParameters * simparams
Pointer to the NAMD simulation input object.
void log(std::string const &message)
ResizeArray< AtomIDList > & modifyRequestedGroups()
Definition: GlobalMaster.C:184
int update_group_properties(int index)
void update_atom_properties(int index)
int size(void) const
Definition: ResizeArray.h:127
cvm::real colvarproxy_namd::temperature ( )
inline

Definition at line 105 of file colvarproxy_namd.h.

References thermostat_temperature.

106  {
107  return thermostat_temperature;
108  }
BigReal thermostat_temperature
Self-explained.
bool colvarproxy_namd::total_forces_enabled ( ) const
inline

Definition at line 82 of file colvarproxy_namd.h.

83  {
84  return total_force_requested;
85  }
void colvarproxy_namd::update_atom_properties ( int  index)

Definition at line 766 of file colvarproxy_namd.C.

References Molecule::atomcharge(), Molecule::atommass(), log(), Node::molecule, and Node::Object().

Referenced by init_atom(), setup(), and update_atoms_map().

767 {
768  Molecule *mol = Node::Object()->molecule;
769  // update mass
770  double const mass = mol->atommass(atoms_ids[index]);
771  if (mass <= 0.001) {
772  this->log("Warning: near-zero mass for atom "+
773  cvm::to_str(atoms_ids[index]+1)+
774  "; expect unstable dynamics if you apply forces to it.\n");
775  }
776  atoms_masses[index] = mass;
777  // update charge
778  atoms_charges[index] = mol->atomcharge(atoms_ids[index]);
779 }
static Node * Object()
Definition: Node.h:86
void log(std::string const &message)
Real atomcharge(int anum) const
Definition: Molecule.h:1048
Real atommass(int anum) const
Definition: Molecule.h:1038
Molecule * molecule
Definition: Node.h:176
int colvarproxy_namd::update_atoms_map ( AtomIDList::const_iterator  begin,
AtomIDList::const_iterator  end 
)

Definition at line 165 of file colvarproxy_namd.C.

References atoms_map, debug, and update_atom_properties().

Referenced by calculate().

167 {
168  for (AtomIDList::const_iterator a_i = begin; a_i != end; a_i++) {
169 
170  if (cvm::debug()) {
171  cvm::log("Updating atoms_map for atom ID "+cvm::to_str(*a_i)+"\n");
172  }
173 
174  if (atoms_map[*a_i] >= 0) continue;
175 
176  for (size_t i = 0; i < atoms_ids.size(); i++) {
177  if (atoms_ids[i] == *a_i) {
178  atoms_map[*a_i] = i;
179  break;
180  }
181  }
182 
183  if (atoms_map[*a_i] < 0) {
184  // this atom is probably managed by another GlobalMaster:
185  // add it here anyway to avoid having to test for array boundaries at each step
186  int const index = add_atom_slot(*a_i);
187  atoms_map[*a_i] = index;
188  update_atom_properties(index);
189  }
190  }
191 
192  if (cvm::debug()) {
193  cvm::log("atoms_map = "+cvm::to_str(atoms_map)+".\n");
194  }
195 
196  return COLVARS_OK;
197 }
const Elem * const_iterator
Definition: ResizeArray.h:38
static int debug
Definition: parm.C:36
std::vector< int > atoms_map
Array of atom indices (relative to the colvarproxy arrays), usedfor faster copy of atomic data...
void update_atom_properties(int index)
int colvarproxy_namd::update_group_properties ( int  index)

Definition at line 1201 of file colvarproxy_namd.C.

References Molecule::atomcharge(), Molecule::atommass(), debug, GlobalMaster::modifyRequestedGroups(), Node::molecule, Node::Object(), and ResizeArray< T >::size().

Referenced by init_atom_group(), and setup().

1202 {
1203  AtomIDList const &namd_group = modifyRequestedGroups()[index];
1204  if (cvm::debug()) {
1205  cvm::log("Re-calculating total mass and charge for scalable group no. "+cvm::to_str(index+1)+" ("+
1206  cvm::to_str(namd_group.size())+" atoms).\n");
1207  }
1208 
1209  cvm::real total_mass = 0.0;
1210  cvm::real total_charge = 0.0;
1211  for (int i = 0; i < namd_group.size(); i++) {
1212  total_mass += Node::Object()->molecule->atommass(namd_group[i]);
1213  total_charge += Node::Object()->molecule->atomcharge(namd_group[i]);
1214  }
1215  atom_groups_masses[index] = total_mass;
1216  atom_groups_charges[index] = total_charge;
1217 
1218  if (cvm::debug()) {
1219  cvm::log("total mass = "+cvm::to_str(total_mass)+
1220  ", total charge = "+cvm::to_str(total_charge)+"\n");
1221  }
1222 
1223  return COLVARS_OK;
1224 }
static Node * Object()
Definition: Node.h:86
static int debug
Definition: parm.C:36
ResizeArray< AtomIDList > & modifyRequestedGroups()
Definition: GlobalMaster.C:184
Real atomcharge(int anum) const
Definition: Molecule.h:1048
Real atommass(int anum) const
Definition: Molecule.h:1038
int size(void) const
Definition: ResizeArray.h:127
Molecule * molecule
Definition: Node.h:176

Friends And Related Function Documentation

friend class cvm::atom
friend

Definition at line 61 of file colvarproxy_namd.h.

Referenced by load_atoms().

Member Data Documentation

std::vector<int> colvarproxy_namd::atoms_map
protected

Array of atom indices (relative to the colvarproxy arrays), usedfor faster copy of atomic data.

Definition at line 40 of file colvarproxy_namd.h.

Referenced by calculate(), reset(), and update_atoms_map().

bool colvarproxy_namd::first_timestep
protected

Definition at line 51 of file colvarproxy_namd.h.

Referenced by calculate(), and colvarproxy_namd().

size_t colvarproxy_namd::previous_NAMD_step
protected

Definition at line 52 of file colvarproxy_namd.h.

Referenced by calculate().

Random colvarproxy_namd::random
protected

NAMD-style PRNG object.

Definition at line 49 of file colvarproxy_namd.h.

Referenced by colvarproxy_namd(), and rand_gaussian().

SubmitReduction* colvarproxy_namd::reduction
protected

Used to submit restraint energy as MISC.

Definition at line 55 of file colvarproxy_namd.h.

Referenced by add_energy(), calculate(), colvarproxy_namd(), and ~colvarproxy_namd().

SimParameters* colvarproxy_namd::simparams
protected

Pointer to the NAMD simulation input object.

Definition at line 43 of file colvarproxy_namd.h.

Referenced by calculate(), colvarproxy_namd(), dt(), init_volmap(), and setup().

BigReal colvarproxy_namd::thermostat_temperature
protected

Self-explained.

Definition at line 46 of file colvarproxy_namd.h.

Referenced by colvarproxy_namd(), and temperature().


The documentation for this class was generated from the following files: