NAMD
colvarproxy_namd.C
Go to the documentation of this file.
1 // -*- c++ -*-
2 
3 // This file is part of the Collective Variables module (Colvars).
4 // The original version of Colvars and its updates are located at:
5 // https://github.com/Colvars/colvars
6 // Please update all Colvars source files before making any changes.
7 // If you wish to distribute your changes, please submit them to the
8 // Colvars repository at GitHub.
9 
10 #include <errno.h>
11 
12 #include "common.h"
13 #include "fstream_namd.h"
14 #include "BackEnd.h"
15 #include "InfoStream.h"
16 #include "Node.h"
17 #include "Molecule.h"
18 #include "PDB.h"
19 #include "PDBData.h"
20 #include "ReductionMgr.h"
21 #include "ScriptTcl.h"
22 
23 #ifdef NAMD_TCL
24 #include <tcl.h>
25 #endif
26 
27 // For replica exchange
28 #include "converse.h"
29 #include "DataExchanger.h"
30 
31 #include "colvarmodule.h"
32 #include "colvaratoms.h"
33 #include "colvarproxy.h"
34 #include "colvarproxy_namd.h"
35 #include "colvarscript.h"
36 
37 
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 }
149 
150 
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 }
163 
164 
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 }
198 
199 
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 }
256 
257 
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 }
276 
277 
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 }
561 
562 
563 // Callback functions
564 
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 }
572 
574 {
575  return colvarproxy::tcl_run_force_callback();
576 }
577 
579  std::string const &name,
580  std::vector<const colvarvalue *> const &cvc_values,
581  colvarvalue &value)
582 {
583  return colvarproxy::tcl_run_colvar_callback(name, cvc_values, value);
584 }
585 
587  std::string const &name,
588  std::vector<const colvarvalue *> const &cvc_values,
589  std::vector<cvm::matrix2d<cvm::real> > &gradient)
590 {
591  return colvarproxy::tcl_run_colvar_gradient_callback(name, cvc_values,
592  gradient);
593 }
594 
595 
596 void colvarproxy_namd::add_energy(cvm::real energy)
597 {
599 }
600 
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 }
612 
613 
614 void colvarproxy_namd::log(std::string const &message)
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 }
622 
623 
624 void colvarproxy_namd::error(std::string const &message)
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 }
643 
644 
645 void colvarproxy_namd::exit(std::string const &message)
646 {
647  log(message);
648  BackEnd::exit();
649 }
650 
651 
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 }
669 
670 
671 int colvarproxy_namd::init_atom(int atom_number)
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 }
696 
697 
698 int colvarproxy_namd::check_atom_id(cvm::residue_id const &residue,
699  std::string const &atom_name,
700  std::string const &segment_id)
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 }
725 
726 
727 
731 int colvarproxy_namd::init_atom(cvm::residue_id const &residue,
732  std::string const &atom_name,
733  std::string const &segment_id)
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 }
757 
758 
760 {
761  colvarproxy::clear_atom(index);
762  // TODO remove it from GlobalMaster arrays?
763 }
764 
765 
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 }
780 
781 
782 cvm::rvector colvarproxy_namd::position_distance(cvm::atom_pos const &pos1,
783  cvm::atom_pos const &pos2)
784  const
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 }
797 
798 
799 
808 };
809 
810 
811 e_pdb_field pdb_field_str2enum(std::string const &pdb_field_str)
812 {
813  e_pdb_field pdb_field = e_pdb_none;
814 
815  if (colvarparse::to_lower_cppstr(pdb_field_str) ==
816  colvarparse::to_lower_cppstr("O")) {
817  pdb_field = e_pdb_occ;
818  }
819 
820  if (colvarparse::to_lower_cppstr(pdb_field_str) ==
821  colvarparse::to_lower_cppstr("B")) {
822  pdb_field = e_pdb_beta;
823  }
824 
825  if (colvarparse::to_lower_cppstr(pdb_field_str) ==
826  colvarparse::to_lower_cppstr("X")) {
827  pdb_field = e_pdb_x;
828  }
829 
830  if (colvarparse::to_lower_cppstr(pdb_field_str) ==
831  colvarparse::to_lower_cppstr("Y")) {
832  pdb_field = e_pdb_y;
833  }
834 
835  if (colvarparse::to_lower_cppstr(pdb_field_str) ==
836  colvarparse::to_lower_cppstr("Z")) {
837  pdb_field = e_pdb_z;
838  }
839 
840  if (pdb_field == e_pdb_none) {
841  cvm::error("Error: unsupported PDB field, \""+
842  pdb_field_str+"\".\n", INPUT_ERROR);
843  }
844 
845  return pdb_field;
846 }
847 
848 
849 int colvarproxy_namd::load_coords(char const *pdb_filename,
850  std::vector<cvm::atom_pos> &pos,
851  const std::vector<int> &indices,
852  std::string const &pdb_field_str,
853  double const pdb_field_value)
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 }
959 
960 
961 int colvarproxy_namd::load_atoms(char const *pdb_filename,
962  cvm::atom_group &atoms,
963  std::string const &pdb_field_str,
964  double const pdb_field_value)
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 }
1016 
1017 
1018 std::ostream * colvarproxy_namd::output_stream(std::string const &output_name,
1019  std::ios_base::openmode mode)
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 }
1041 
1042 
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 }
1055 
1056 
1057 int colvarproxy_namd::close_output_stream(std::string const &output_name)
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 }
1074 
1075 
1076 int colvarproxy_namd::backup_file(char const *filename)
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 }
1085 
1086 
1087 char const *colvarproxy_namd::script_obj_to_str(unsigned char *obj)
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 }
1096 
1097 
1098 std::vector<std::string> colvarproxy_namd::script_obj_to_str_vector(unsigned char *obj)
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 }
1122 
1123 
1124 int colvarproxy_namd::init_atom_group(std::vector<int> const &atoms_ids)
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 }
1192 
1193 
1195 {
1196  // do nothing, keep the NAMD arrays in sync with the colvarproxy ones
1197  colvarproxy::clear_atom_group(index);
1198 }
1199 
1200 
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 }
1225 
1226 
1227 
1228 int colvarproxy_namd::set_unit_system(std::string const &units_in, bool /*check_only*/)
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 }
1236 
1237 
1238 #if NAMD_VERSION_NUMBER >= 34471681
1239 
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 }
1262 
1263 
1264 int colvarproxy_namd::init_volmap(const char *volmap_name)
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 }
1288 
1289 
1291 {
1292  colvarproxy::clear_volmap(index);
1293 }
1294 
1295 #endif
1296 
1297 
1298 #if CMK_SMP && USE_CKLOOP // SMP only
1299 
1300 void calc_colvars_items_smp(int first, int last, void *result, int paramNum, void *param)
1301 {
1302  colvarproxy_namd *proxy = (colvarproxy_namd *) param;
1303  colvarmodule *cv = proxy->colvars;
1304 
1305  cvm::increase_depth();
1306  for (int i = first; i <= last; i++) {
1307  colvar *x = (*(cv->variables_active_smp()))[i];
1308  int x_item = (*(cv->variables_active_smp_items()))[i];
1309  if (cvm::debug()) {
1310  cvm::log("["+cvm::to_str(proxy->smp_thread_id())+"/"+cvm::to_str(proxy->smp_num_threads())+
1311  "]: calc_colvars_items_smp(), first = "+cvm::to_str(first)+
1312  ", last = "+cvm::to_str(last)+", cv = "+
1313  x->name+", cvc = "+cvm::to_str(x_item)+"\n");
1314  }
1315  x->calc_cvcs(x_item, 1);
1316  }
1317  cvm::decrease_depth();
1318 }
1319 
1320 
1321 int colvarproxy_namd::smp_colvars_loop()
1322 {
1323  colvarmodule *cv = this->colvars;
1324  CkLoop_Parallelize(calc_colvars_items_smp, 1, this,
1325  cv->variables_active_smp()->size(),
1326  0, cv->variables_active_smp()->size()-1);
1327  return cvm::get_error();
1328 }
1329 
1330 
1331 void calc_cv_biases_smp(int first, int last, void *result, int paramNum, void *param)
1332 {
1333  colvarproxy_namd *proxy = (colvarproxy_namd *) param;
1334  colvarmodule *cv = proxy->colvars;
1335 
1336  cvm::increase_depth();
1337  for (int i = first; i <= last; i++) {
1338  colvarbias *b = (*(cv->biases_active()))[i];
1339  if (cvm::debug()) {
1340  cvm::log("["+cvm::to_str(proxy->smp_thread_id())+"/"+cvm::to_str(proxy->smp_num_threads())+
1341  "]: calc_cv_biases_smp(), first = "+cvm::to_str(first)+
1342  ", last = "+cvm::to_str(last)+", bias = "+
1343  b->name+"\n");
1344  }
1345  b->update();
1346  }
1347  cvm::decrease_depth();
1348 }
1349 
1350 
1351 int colvarproxy_namd::smp_biases_loop()
1352 {
1353  colvarmodule *cv = this->colvars;
1354  CkLoop_Parallelize(calc_cv_biases_smp, 1, this,
1355  cv->biases_active()->size(), 0, cv->biases_active()->size()-1);
1356  return cvm::get_error();
1357 }
1358 
1359 
1360 void calc_cv_scripted_forces(int paramNum, void *param)
1361 {
1362  colvarproxy_namd *proxy = (colvarproxy_namd *) param;
1363  colvarmodule *cv = proxy->colvars;
1364  if (cvm::debug()) {
1365  cvm::log("["+cvm::to_str(proxy->smp_thread_id())+"/"+cvm::to_str(proxy->smp_num_threads())+
1366  "]: calc_cv_scripted_forces()\n");
1367  }
1368  cv->calc_scripted_forces();
1369 }
1370 
1371 
1372 int colvarproxy_namd::smp_biases_script_loop()
1373 {
1374  colvarmodule *cv = this->colvars;
1375  CkLoop_Parallelize(calc_cv_biases_smp, 1, this,
1376  cv->biases_active()->size(), 0, cv->biases_active()->size()-1,
1377  1, NULL, CKLOOP_NONE,
1378  calc_cv_scripted_forces, 1, this);
1379  return cvm::get_error();
1380 }
1381 
1382 #endif // #if CMK_SMP && USE_CKLOOP
1383 
1384 
1386 #if CMK_HAS_PARTITION
1387  return COLVARS_OK;
1388 #else
1389  return COLVARS_NOT_IMPLEMENTED;
1390 #endif
1391 }
1392 
1393 
1395  return CmiMyPartition();
1396 }
1397 
1398 
1400  return CmiNumPartitions();
1401 }
1402 
1403 
1405  replica_barrier();
1406 }
1407 
1408 
1409 int colvarproxy_namd::replica_comm_recv(char* msg_data, int buf_len,
1410  int src_rep) {
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 }
1423 
1424 
1425 int colvarproxy_namd::replica_comm_send(char* msg_data, int msg_len,
1426  int dest_rep) {
1427  replica_send(msg_data, msg_len, dest_rep, CkMyPe());
1428  return msg_len;
1429 }
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
ForceList & modifyAppliedForces()
Definition: GlobalMaster.C:162
char restartFilename[128]
AtomIDList & modifyRequestedAtoms()
Definition: GlobalMaster.C:127
cvm::rvector position_distance(cvm::atom_pos const &pos1, cvm::atom_pos const &pos2) const
BigRealList & modifyGridObjForces()
Definition: GlobalMaster.C:179
Definition: PDB.h:35
SimParameters * simparams
Pointer to the NAMD simulation input object.
Random random
NAMD-style PRNG object.
const Elem * const_iterator
Definition: ResizeArray.h:38
void NAMD_err(const char *err_msg)
Definition: common.C:102
ScriptTcl * getScript(void)
Definition: Node.h:192
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)
AtomIDList::const_iterator getForceIdEnd()
Definition: GlobalMaster.C:260
AtomIDList::const_iterator getAtomIdBegin()
Definition: GlobalMaster.C:190
PositionList::const_iterator getGroupPositionEnd()
Definition: GlobalMaster.C:206
Communication between colvars and NAMD (implementation of colvarproxy)
static void exit(int status=0)
Definition: BackEnd.C:276
Definition: Vector.h:64
SimParameters * simParameters
Definition: Node.h:178
static __thread atom * atoms
e_pdb_field pdb_field_str2enum(std::string const &pdb_field_str)
BigReal reassignTemp
BigReal & item(int i)
Definition: ReductionMgr.h:312
void replica_send(const char *sndbuf, int sendcount, int destPart, int destPE)
BigReal z
Definition: Vector.h:66
int init_atom(int atom_number)
ForceList::const_iterator getGroupTotalForceBegin()
Definition: GlobalMaster.C:210
int numGridforceGrids
Definition: Molecule.h:590
SubmitReduction * willSubmit(int setID, int size=-1)
Definition: ReductionMgr.C:365
int orthogonal() const
Definition: Lattice.h:257
static int debug
Definition: parm.C:36
int close_output_stream(std::string const &output_name)
static ReductionMgr * Object(void)
Definition: ReductionMgr.h:278
#define iout
Definition: InfoStream.h:87
int num_atoms(void)
Definition: PDB.C:323
PositionList::const_iterator getAtomPositionBegin()
Definition: GlobalMaster.C:198
const ResizeArray< AtomIDList > & requestedGroups()
Definition: GlobalMaster.C:149
AtomIDList::const_iterator getForceIdBegin()
Definition: GlobalMaster.C:255
int init_atom_group(std::vector< int > const &atoms_ids)
int run_colvar_callback(std::string const &name, std::vector< const colvarvalue * > const &cvcs, colvarvalue &value)
#define COLVARPROXY_VERSION
char const * script_obj_to_str(unsigned char *obj)
std::vector< int > atoms_map
Array of atom indices (relative to the colvarproxy arrays), usedfor faster copy of atomic data...
void replica_recv(DataMessage **precvMsg, int srcPart, int srcPE)
char outputFilename[128]
Definition: Random.h:37
int index_for_key(const char *key)
e_pdb_field
void clear_atom(int index)
IntList::const_iterator getGridObjIndexBegin()
Definition: GlobalMaster.C:218
PDBAtom * atom(int place)
Definition: PDB.C:394
BigReal rescaleTemp
int run_colvar_gradient_callback(std::string const &name, std::vector< const colvarvalue * > const &cvcs, std::vector< cvm::matrix2d< cvm::real > > &gradient)
iterator end(void)
Definition: ResizeArray.h:37
void request_total_force(bool yesno)
friend class cvm::atom
bool is_open() const
Definition: fstream_namd.h:30
IntList & modifyRequestedGridObjects()
Definition: GlobalMaster.C:173
int set_unit_system(std::string const &units_in, bool check_only)
void clear_volmap(int index)
void clear()
Definition: ResizeArray.h:87
BigReal langevinTemp
Vector delta(const Position &pos1, const Position &pos2) const
Definition: Lattice.h:144
void setall(const Elem &elem)
Definition: ResizeArray.h:90
BigReal x
Definition: Vector.h:66
std::ostream * output_stream(std::string const &output_name, std::ios_base::openmode mode)
void log(std::string const &message)
void replica_barrier()
AtomIDList & modifyForcedAtoms()
Definition: GlobalMaster.C:157
int numAtoms
Definition: Molecule.h:556
void exit(std::string const &message)
ResizeArray< AtomIDList > & modifyRequestedGroups()
Definition: GlobalMaster.C:184
int update_group_properties(int index)
GridforceGrid * get_gridfrc_grid(int gridnum) const
Definition: Molecule.h:1276
void NAMD_die(const char *err_msg)
Definition: common.C:83
ForceList & modifyGroupForces()
Definition: GlobalMaster.C:167
ConfigList * configList
Definition: Node.h:179
void error(std::string const &message)
virtual int replica_comm_recv(char *msg_data, int buf_len, int src_rep)
virtual int replica_index()
MGridforceParamsList mgridforcelist
void update_atom_properties(int index)
void NAMD_backup_file(const char *filename, const char *extension)
Definition: common.C:159
virtual Vector get_scale(void) const =0
Bool set(const char *s)
Definition: Vector.h:228
BigRealList::const_iterator getGridObjValueBegin()
Definition: GlobalMaster.C:226
std::vector< std::string > script_obj_to_str_vector(unsigned char *obj)
int add(const Elem &elem)
Definition: ResizeArray.h:97
void clear_atom_group(int index)
IntList::const_iterator getGridObjIndexEnd()
Definition: GlobalMaster.C:222
const Lattice * lattice
Definition: GlobalMaster.h:141
unsigned int randomSeed
void requestTotalForce(bool yesno=true)
Definition: GlobalMaster.h:133
char data[1]
Definition: DataExchanger.h:23
void resize(int i)
Definition: ResizeArray.h:84
AtomIDList::const_iterator getAtomIdEnd()
Definition: GlobalMaster.C:194
char * data
Definition: ConfigList.h:48
int check_atom_id(int atom_number)
virtual int replica_enabled()
virtual void replica_comm_barrier()
BigReal y
Definition: Vector.h:66
Vector b() const
Definition: Lattice.h:253
Real atomcharge(int anum) const
Definition: Molecule.h:1048
StringList * find(const char *name) const
Definition: ConfigList.C:341
BigReal thermostat_temperature
Self-explained.
colvarmodule * colvars
Definition: Node.h:184
BigReal tCoupleTemp
int load_atoms(char const *filename, cvm::atom_group &atoms, std::string const &pdb_field, double const pdb_field_value=0.0)
void add_energy(cvm::real energy)
Real atommass(int anum) const
Definition: Molecule.h:1038
int flush_output_stream(std::ostream *os)
void submit(void)
Definition: ReductionMgr.h:323
int size(void) const
Definition: ResizeArray.h:127
infostream & endi(infostream &s)
Definition: InfoStream.C:38
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 backup_file(char const *filename)
gridSize x
int init_volmap(int volmap_id)
virtual int replica_comm_send(char *msg_data, int msg_len, int dest_rep)
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
virtual int num_replicas()
Vector a() const
Definition: Lattice.h:252
ForceList::const_iterator getTotalForce()
Definition: GlobalMaster.C:265
Vector c() const
Definition: Lattice.h:254
BigReal stochRescaleTemp
BigRealList::const_iterator getGridObjValueEnd()
Definition: GlobalMaster.C:230
int c_p() const
Definition: Lattice.h:275
BigReal loweAndersenTemp
ForceList::const_iterator getGroupTotalForceEnd()
Definition: GlobalMaster.C:214
iterator begin(void)
Definition: ResizeArray.h:36