colvarproxy_namd Class Reference

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

#include <colvarproxy_namd.h>

Inheritance diagram for colvarproxy_namd:
GlobalMaster

List of all members.

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)

Protected Member Functions

void init_tcl_pointers ()

Protected Attributes

std::vector< int > atoms_map
 Array of atom indices (relative to the colvarproxy arrays), usedfor faster copy of atomic data.
SimParameterssimparams
 Pointer to the NAMD simulation input object.
BigReal thermostat_temperature
 Self-explained.
Random random
 NAMD-style PRNG object.
bool first_timestep
size_t previous_NAMD_step
SubmitReductionreduction
 Used to submit restraint energy as MISC.

Friends

class cvm::atom

Detailed Description

Communication between colvars and NAMD (implementation of colvarproxy).

Definition at line 29 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, log(), SimParameters::loweAndersenOn, SimParameters::loweAndersenTemp, ReductionMgr::Object(), Node::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().

00039 {
00040   version_int = get_version_from_string(COLVARPROXY_VERSION);
00041 
00042   first_timestep = true;
00043   requestTotalForce(total_force_requested);
00044 
00045   angstrom_value = 1.;
00046 
00047   // initialize pointers to NAMD configuration data
00048   simparams = Node::Object()->simParameters;
00049 
00050   if (cvm::debug())
00051     iout << "Info: initializing the colvars proxy object.\n" << endi;
00052 
00053   // find the configuration file, if provided
00054   StringList *config = Node::Object()->configList->find("colvarsConfig");
00055 
00056   // find the input state file
00057   StringList *input_restart = Node::Object()->configList->find("colvarsInput");
00058   input_prefix_str = std::string(input_restart ? input_restart->data : "");
00059   if (input_prefix_str.rfind(".colvars.state") != std::string::npos) {
00060     // strip the extension, if present
00061     input_prefix_str.erase(input_prefix_str.rfind(".colvars.state"),
00062                            std::string(".colvars.state").size());
00063   }
00064 
00065   // get the thermostat temperature
00066   if (simparams->rescaleFreq > 0)
00067     thermostat_temperature = simparams->rescaleTemp;
00068   else if (simparams->reassignFreq > 0)
00069     thermostat_temperature = simparams->reassignTemp;
00070   else if (simparams->langevinOn)
00071     thermostat_temperature = simparams->langevinTemp;
00072   else if (simparams->tCoupleOn)
00073     thermostat_temperature = simparams->tCoupleTemp;
00074   else if (simparams->loweAndersenOn)
00075     thermostat_temperature = simparams->loweAndersenTemp;
00076   else if (simparams->stochRescaleOn)
00077     thermostat_temperature = simparams->stochRescaleTemp;
00078   else
00079     thermostat_temperature = 0.0;
00080 
00081   random = Random(simparams->randomSeed);
00082 
00083   // both fields are taken from data structures already available
00084   updated_masses_ = updated_charges_ = true;
00085 
00086   // take the output prefixes from the namd input
00087   output_prefix_str = std::string(simparams->outputFilename);
00088   restart_output_prefix_str = std::string(simparams->restartFilename);
00089   restart_frequency_engine = simparams->restartFrequency;
00090 
00091   // check if it is possible to save output configuration
00092   if ((!output_prefix_str.size()) && (!restart_output_prefix_str.size())) {
00093     error("Error: neither the final output state file or "
00094           "the output restart file could be defined, exiting.\n");
00095   }
00096 
00097 
00098 #ifdef NAMD_TCL
00099   have_scripts = true;
00100 
00101   init_tcl_pointers();
00102 
00103   // See is user-scripted forces are defined
00104   if (Tcl_FindCommand(reinterpret_cast<Tcl_Interp *>(tcl_interp_),
00105                       "calc_colvar_forces", NULL, 0) == NULL) {
00106     force_script_defined = false;
00107   } else {
00108     force_script_defined = true;
00109   }
00110 #else
00111   force_script_defined = false;
00112   have_scripts = false;
00113 #endif
00114 
00115 
00116   // initialize module: this object will be the communication proxy
00117   colvars = new colvarmodule(this);
00118   cvm::log("Using NAMD interface, version "+
00119            cvm::to_str(COLVARPROXY_VERSION)+".\n");
00120 
00121   errno = 0;
00122   if (config) {
00123     colvars->read_config_file(config->data);
00124   }
00125 
00126   colvars->setup();
00127   colvars->setup_input();
00128   colvars->setup_output();
00129 
00130   // save to Node for Tcl script access
00131   Node::Object()->colvars = colvars;
00132 
00133 #ifdef NAMD_TCL
00134   // Construct instance of colvars scripting interface
00135   script = new colvarscript(this);
00136 #endif
00137 
00138   if (simparams->firstTimestep != 0) {
00139     cvm::log("Initializing step number as firstTimestep.\n");
00140     colvars->it = colvars->it_restart =
00141       static_cast<cvm::step_number>(simparams->firstTimestep);
00142   }
00143 
00144   reduction = ReductionMgr::Object()->willSubmit(REDUCTIONS_BASIC);
00145 
00146   if (cvm::debug())
00147     iout << "Info: done initializing the colvars proxy object.\n" << endi;
00148 }

colvarproxy_namd::~colvarproxy_namd (  ) 

Definition at line 151 of file colvarproxy_namd.C.

References reduction.

00152 {
00153   delete reduction;
00154   if (script != NULL) {
00155     delete script;
00156     script = NULL;
00157   }
00158   if (colvars != NULL) {
00159     delete colvars;
00160     colvars = NULL;
00161   }
00162 }


Member Function Documentation

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

Definition at line 578 of file colvarproxy_namd.C.

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

00579 {
00580   reduction->item(REDUCTION_MISC_ENERGY) += energy;
00581 }

cvm::real colvarproxy_namd::backend_angstrom_value (  )  [inline]

Definition at line 90 of file colvarproxy_namd.h.

00091   {
00092     return 1.0;
00093   }

int colvarproxy_namd::backup_file ( char const *  filename  ) 

Definition at line 1058 of file colvarproxy_namd.C.

References NAMD_backup_file().

Referenced by output_stream().

01059 {
01060   if (std::string(filename).rfind(std::string(".colvars.state")) != std::string::npos) {
01061     NAMD_backup_file(filename, ".old");
01062   } else {
01063     NAMD_backup_file(filename, ".BAK");
01064   }
01065   return COLVARS_OK;
01066 }

cvm::real colvarproxy_namd::boltzmann (  )  [inline]

Definition at line 95 of file colvarproxy_namd.h.

00096   {
00097     return 0.001987191;
00098   }

void colvarproxy_namd::calculate (  )  [virtual]

Reimplemented from GlobalMaster.

Definition at line 274 of file colvarproxy_namd.C.

References Lattice::a(), Lattice::a_p(), ResizeArray< Elem >::add(), atoms_map, Lattice::b(), Lattice::b_p(), ResizeArray< Elem >::begin(), Lattice::c(), Lattice::c_p(), ResizeArray< Elem >::clear(), debug, ResizeArray< Elem >::end(), error(), 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< Elem >::resize(), Vector::set(), ResizeArray< Elem >::setall(), setup(), simparams, GlobalMaster::step, SubmitReduction::submit(), update_atoms_map(), Vector::x, Vector::y, and Vector::z.

00275 {
00276   errno = 0;
00277 
00278   if (first_timestep) {
00279 
00280     colvarproxy_namd::setup();
00281     colvars->setup();
00282     colvars->setup_input();
00283     colvars->setup_output();
00284 
00285     first_timestep = false;
00286 
00287   } else {
00288     // Use the time step number inherited from GlobalMaster
00289     if ( step - previous_NAMD_step == 1 ) {
00290       colvars->it++;
00291     }
00292     // Other cases could mean:
00293     // - run 0
00294     // - beginning of a new run statement
00295     // then the internal counter should not be incremented
00296   }
00297 
00298   previous_NAMD_step = step;
00299 
00300   {
00301     Vector const a = lattice->a();
00302     Vector const b = lattice->b();
00303     Vector const c = lattice->c();
00304     unit_cell_x.set(a.x, a.y, a.z);
00305     unit_cell_y.set(b.x, b.y, c.z);
00306     unit_cell_z.set(c.x, c.y, c.z);
00307   }
00308 
00309   if (!lattice->a_p() && !lattice->b_p() && !lattice->c_p()) {
00310     boundaries_type = boundaries_non_periodic;
00311     reset_pbc_lattice();
00312   } else if (lattice->a_p() && lattice->b_p() && lattice->c_p()) {
00313     if (lattice->orthogonal()) {
00314       boundaries_type = boundaries_pbc_ortho;
00315     } else {
00316       boundaries_type = boundaries_pbc_triclinic;
00317     }
00318     colvarproxy_system::update_pbc_lattice();
00319   } else {
00320     boundaries_type = boundaries_unsupported;
00321   }
00322 
00323   if (cvm::debug()) {
00324     cvm::log(std::string(cvm::line_marker)+
00325              "colvarproxy_namd, step no. "+cvm::to_str(colvars->it)+"\n"+
00326              "Updating atomic data arrays.\n");
00327   }
00328 
00329   // must delete the forces applied at the previous step: we can do
00330   // that because they have already been used and copied to other
00331   // memory locations
00332   modifyForcedAtoms().clear();
00333   modifyAppliedForces().clear();
00334 
00335   // prepare local arrays
00336   for (size_t i = 0; i < atoms_ids.size(); i++) {
00337     atoms_positions[i] = cvm::rvector(0.0, 0.0, 0.0);
00338     atoms_total_forces[i] = cvm::rvector(0.0, 0.0, 0.0);
00339     atoms_new_colvar_forces[i] = cvm::rvector(0.0, 0.0, 0.0);
00340   }
00341 
00342   for (size_t i = 0; i < atom_groups_ids.size(); i++) {
00343     atom_groups_total_forces[i] = cvm::rvector(0.0, 0.0, 0.0);
00344     atom_groups_new_colvar_forces[i] = cvm::rvector(0.0, 0.0, 0.0);
00345   }
00346 
00347   for (int imap = 0; imap < volmaps_ids.size(); imap++) {
00348     volmaps_new_colvar_forces[imap] = 0.0;
00349   }
00350 
00351   // create the atom map if needed
00352   size_t const n_all_atoms = Node::Object()->molecule->numAtoms;
00353   if (atoms_map.size() != n_all_atoms) {
00354     atoms_map.resize(n_all_atoms);
00355     atoms_map.assign(n_all_atoms, -1);
00356     update_atoms_map(getAtomIdBegin(), getAtomIdEnd());
00357   }
00358 
00359   // if new atomic positions or forces have been communicated by other GlobalMasters, add them to the atom map
00360   if ((int(atoms_ids.size()) < (getAtomIdEnd() - getAtomIdBegin())) ||
00361       (int(atoms_ids.size()) < (getForceIdEnd() - getForceIdBegin()))) {
00362     update_atoms_map(getAtomIdBegin(), getAtomIdEnd());
00363     update_atoms_map(getForceIdBegin(), getForceIdEnd());
00364   }
00365 
00366   {
00367     if (cvm::debug()) {
00368       cvm::log("Updating positions arrays.\n");
00369     }
00370     size_t n_positions = 0;
00371     AtomIDList::const_iterator a_i = getAtomIdBegin();
00372     AtomIDList::const_iterator a_e = getAtomIdEnd();
00373     PositionList::const_iterator p_i = getAtomPositionBegin();
00374 
00375     for ( ; a_i != a_e; ++a_i, ++p_i ) {
00376       atoms_positions[atoms_map[*a_i]] = cvm::rvector((*p_i).x, (*p_i).y, (*p_i).z);
00377       n_positions++;
00378     }
00379 
00380     // The following had to be relaxed because some atoms may be forced without their position being requested
00381     // if (n_positions < atoms_ids.size()) {
00382     //   cvm::error("Error: did not receive the positions of all atoms.\n", BUG_ERROR);
00383     // }
00384   }
00385 
00386   if (total_force_requested && cvm::step_relative() > 0) {
00387 
00388     // sort the force arrays the previous step
00389     // (but only do so if there *is* a previous step!)
00390 
00391     {
00392       if (cvm::debug()) {
00393         cvm::log("Updating total forces arrays.\n");
00394       }
00395       size_t n_total_forces = 0;
00396       AtomIDList::const_iterator a_i = getForceIdBegin();
00397       AtomIDList::const_iterator a_e = getForceIdEnd();
00398       ForceList::const_iterator f_i = getTotalForce();
00399 
00400       for ( ; a_i != a_e; ++a_i, ++f_i ) {
00401         atoms_total_forces[atoms_map[*a_i]] = cvm::rvector((*f_i).x, (*f_i).y, (*f_i).z);
00402         n_total_forces++;
00403       }
00404 
00405       if (n_total_forces < atoms_ids.size()) {
00406         cvm::error("Error: total forces were requested, but total forces "
00407                    "were not received for all atoms.\n"
00408                    "The most probable cause is combination of energy "
00409                    "minimization with a biasing method that requires MD (e.g. ABF).\n"
00410                    "Always run minimization and ABF separately.", INPUT_ERROR);
00411       }
00412     }
00413 
00414     {
00415       if (cvm::debug()) {
00416         cvm::log("Updating group total forces arrays.\n");
00417       }
00418       ForceList::const_iterator f_i = getGroupTotalForceBegin();
00419       ForceList::const_iterator f_e = getGroupTotalForceEnd();
00420       size_t i = 0;
00421       if ((f_e - f_i) != ((int) atom_groups_ids.size())) {
00422         cvm::error("Error: total forces were requested for scalable groups, "
00423                    "but they are not in the same number from the number of groups.\n"
00424                    "The most probable cause is combination of energy "
00425                    "minimization with a biasing method that requires MD (e.g. ABF).\n"
00426                    "Always run minimization and ABF separately.", INPUT_ERROR);
00427       }
00428       for ( ; f_i != f_e; f_i++, i++) {
00429         atom_groups_total_forces[i] = cvm::rvector((*f_i).x, (*f_i).y, (*f_i).z);
00430       }
00431     }
00432   }
00433 
00434   {
00435     if (cvm::debug()) {
00436       cvm::log("Updating group positions arrays.\n");
00437     }
00438     // update group data (only coms available so far)
00439     size_t ig;
00440     // note: getGroupMassBegin() could be used here, but masses and charges
00441     // have already been calculated from the last call to setup()
00442     PositionList::const_iterator gp_i = getGroupPositionBegin();
00443     for (ig = 0; gp_i != getGroupPositionEnd(); gp_i++, ig++) {
00444       atom_groups_coms[ig] = cvm::rvector(gp_i->x, gp_i->y, gp_i->z);
00445     }
00446   }
00447 
00448   {
00449     if (cvm::debug()) {
00450       log("Updating grid objects.\n");
00451     }
00452     // Using a simple nested loop: there probably won't be so many maps that
00453     // this becomes performance-limiting
00454     IntList::const_iterator goi_i = getGridObjIndexBegin();
00455     BigRealList::const_iterator gov_i = getGridObjValueBegin();
00456     for ( ; gov_i != getGridObjValueEnd(); goi_i++, gov_i++) {
00457       for (size_t imap = 0; imap < volmaps_ids.size(); imap++) {
00458         if (volmaps_ids[imap] == *goi_i) {
00459           volmaps_values[imap] = *gov_i;
00460           break;
00461         }
00462       }
00463     }
00464   }
00465 
00466   if (cvm::debug()) {
00467     cvm::log("atoms_ids = "+cvm::to_str(atoms_ids)+"\n");
00468     cvm::log("atoms_ncopies = "+cvm::to_str(atoms_ncopies)+"\n");
00469     cvm::log("atoms_masses = "+cvm::to_str(atoms_masses)+"\n");
00470     cvm::log("atoms_charges = "+cvm::to_str(atoms_charges)+"\n");
00471     cvm::log("atoms_positions = "+cvm::to_str(atoms_positions)+"\n");
00472     cvm::log("atoms_total_forces = "+cvm::to_str(atoms_total_forces)+"\n");
00473     cvm::log(cvm::line_marker);
00474 
00475     cvm::log("atom_groups_ids = "+cvm::to_str(atom_groups_ids)+"\n");
00476     cvm::log("atom_groups_ncopies = "+cvm::to_str(atom_groups_ncopies)+"\n");
00477     cvm::log("atom_groups_masses = "+cvm::to_str(atom_groups_masses)+"\n");
00478     cvm::log("atom_groups_charges = "+cvm::to_str(atom_groups_charges)+"\n");
00479     cvm::log("atom_groups_coms = "+cvm::to_str(atom_groups_coms)+"\n");
00480     cvm::log("atom_groups_total_forces = "+cvm::to_str(atom_groups_total_forces)+"\n");
00481     cvm::log(cvm::line_marker);
00482 
00483     cvm::log("volmaps_ids = "+cvm::to_str(volmaps_ids)+"\n");
00484     cvm::log("volmaps_values = "+cvm::to_str(volmaps_values)+"\n");
00485     cvm::log(cvm::line_marker);
00486   }
00487 
00488   // call the collective variable module
00489   if (colvars->calc() != COLVARS_OK) {
00490     cvm::error("Error in the collective variables module.\n", COLVARS_ERROR);
00491   }
00492 
00493   if (cvm::debug()) {
00494     cvm::log(cvm::line_marker);
00495     cvm::log("atoms_new_colvar_forces = "+cvm::to_str(atoms_new_colvar_forces)+"\n");
00496     cvm::log(cvm::line_marker);
00497     cvm::log("atom_groups_new_colvar_forces = "+cvm::to_str(atom_groups_new_colvar_forces)+"\n");
00498     cvm::log(cvm::line_marker);
00499     cvm::log("volmaps_new_colvar_forces = "+cvm::to_str(volmaps_new_colvar_forces)+"\n");
00500     cvm::log(cvm::line_marker);
00501   }
00502 
00503   // communicate all forces to the MD integrator
00504   for (size_t i = 0; i < atoms_ids.size(); i++) {
00505     cvm::rvector const &f = atoms_new_colvar_forces[i];
00506     modifyForcedAtoms().add(atoms_ids[i]);
00507     modifyAppliedForces().add(Vector(f.x, f.y, f.z));
00508   }
00509 
00510   if (atom_groups_new_colvar_forces.size() > 0) {
00511     modifyGroupForces().resize(requestedGroups().size());
00512     ForceList::iterator gf_i = modifyGroupForces().begin();
00513     for (int ig = 0; gf_i != modifyGroupForces().end(); gf_i++, ig++) {
00514       cvm::rvector const &f = atom_groups_new_colvar_forces[ig];
00515       *gf_i = Vector(f.x, f.y, f.z);
00516     }
00517   }
00518 
00519   if (volmaps_new_colvar_forces.size() > 0) {
00520     modifyGridObjForces().resize(requestedGridObjs().size());
00521     modifyGridObjForces().setall(0.0);
00522     IntList::const_iterator goi_i = getGridObjIndexBegin();
00523     BigRealList::iterator gof_i = modifyGridObjForces().begin();
00524     for ( ; goi_i != getGridObjIndexEnd(); goi_i++, gof_i++) {
00525       for (size_t imap = 0; imap < volmaps_ids.size(); imap++) {
00526         if (volmaps_ids[imap] == *goi_i) {
00527           *gof_i = volmaps_new_colvar_forces[imap];
00528           break;
00529         }
00530       }
00531     }
00532   }
00533 
00534   // send MISC energy
00535   reduction->submit();
00536 
00537   // NAMD does not destruct GlobalMaster objects, so we must remember
00538   // to write all output files at the end of a run
00539   if (step == simparams->N) {
00540     post_run();
00541   }
00542 }

int colvarproxy_namd::check_atom_id ( cvm::residue_id const &  residue,
std::string const &  atom_name,
std::string const &  segment_id 
)

Definition at line 680 of file colvarproxy_namd.C.

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

00683 {
00684   int const aid =
00685     (segment_id.size() ?
00686      Node::Object()->molecule->get_atom_from_name(segment_id.c_str(),
00687                                                   residue,
00688                                                   atom_name.c_str()) :
00689      Node::Object()->molecule->get_atom_from_name("MAIN",
00690                                                   residue,
00691                                                   atom_name.c_str()));
00692 
00693   if (aid < 0) {
00694     // get_atom_from_name() has returned an error value
00695     cvm::error("Error: could not find atom \""+
00696                atom_name+"\" in residue "+
00697                cvm::to_str(residue)+
00698                ( (segment_id != "MAIN") ?
00699                  (", segment \""+segment_id+"\"") :
00700                  ("") )+
00701                "\n", INPUT_ERROR);
00702     return INPUT_ERROR;
00703   }
00704 
00705   return aid;
00706 }

int colvarproxy_namd::check_atom_id ( int  atom_number  ) 

Definition at line 634 of file colvarproxy_namd.C.

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

Referenced by init_atom().

00635 {
00636   // NAMD's internal numbering starts from zero
00637   int const aid = (atom_number-1);
00638 
00639   if (cvm::debug())
00640     cvm::log("Adding atom "+cvm::to_str(atom_number)+
00641         " for collective variables calculation.\n");
00642 
00643   if ( (aid < 0) || (aid >= Node::Object()->molecule->numAtoms) ) {
00644     cvm::error("Error: invalid atom number specified, "+
00645                cvm::to_str(atom_number)+"\n", INPUT_ERROR);
00646     return INPUT_ERROR;
00647   }
00648 
00649   return aid;
00650 }

void colvarproxy_namd::clear_atom ( int  index  ) 

Definition at line 741 of file colvarproxy_namd.C.

00742 {
00743   colvarproxy::clear_atom(index);
00744   // TODO remove it from GlobalMaster arrays?
00745 }

void colvarproxy_namd::clear_atom_group ( int  index  ) 

Definition at line 1176 of file colvarproxy_namd.C.

01177 {
01178   // do nothing, keep the NAMD arrays in sync with the colvarproxy ones
01179   colvarproxy::clear_atom_group(index);
01180 }

void colvarproxy_namd::clear_volmap ( int  index  ) 

Definition at line 1270 of file colvarproxy_namd.C.

01271 {
01272   colvarproxy::clear_volmap(index);
01273 }

int colvarproxy_namd::close_output_stream ( std::string const &  output_name  ) 

Definition at line 1039 of file colvarproxy_namd.C.

01040 {
01041   std::list<std::ostream *>::iterator osi  = output_files.begin();
01042   std::list<std::string>::iterator    osni = output_stream_names.begin();
01043   for ( ; osi != output_files.end(); osi++, osni++) {
01044     if (*osni == output_name) {
01045       if (((ofstream_namd *) *osi)->is_open()) {
01046         ((ofstream_namd *) *osi)->close();
01047       }
01048       delete *osi;
01049       output_files.erase(osi);
01050       output_stream_names.erase(osni);
01051       return COLVARS_OK;
01052     }
01053   }
01054   return COLVARS_ERROR;
01055 }

cvm::real colvarproxy_namd::dt (  )  [inline]
void colvarproxy_namd::error ( std::string const &  message  ) 

Definition at line 606 of file colvarproxy_namd.C.

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

Referenced by calculate(), check_atom_id(), colvarproxy_namd(), init_atom_group(), init_volmap(), load_coords(), output_stream(), and set_unit_system().

00607 {
00608   log(message);
00609   switch (cvm::get_error()) {
00610   case FILE_ERROR:
00611     errno = EIO; break;
00612   case COLVARS_NOT_IMPLEMENTED:
00613     errno = ENOSYS; break;
00614   case MEMORY_ERROR:
00615     errno = ENOMEM; break;
00616   }
00617   char const *msg = "Error in the collective variables module "
00618     "(see above for details)";
00619   if (errno) {
00620     NAMD_err(msg);
00621   } else {
00622     NAMD_die(msg);
00623   }
00624 }

void colvarproxy_namd::exit ( std::string const &  message  ) 

Definition at line 627 of file colvarproxy_namd.C.

References log().

00628 {
00629   log(message);
00630   BackEnd::exit();
00631 }

int colvarproxy_namd::flush_output_stream ( std::ostream *  os  ) 

Definition at line 1025 of file colvarproxy_namd.C.

01026 {
01027   std::list<std::ostream *>::iterator osi  = output_files.begin();
01028   std::list<std::string>::iterator    osni = output_stream_names.begin();
01029   for ( ; osi != output_files.end(); osi++, osni++) {
01030     if (*osi == os) {
01031       ((ofstream_namd *) *osi)->flush();
01032       return COLVARS_OK;
01033     }
01034   }
01035   return COLVARS_ERROR;
01036 }

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 713 of file colvarproxy_namd.C.

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

00716 {
00717   int const aid = check_atom_id(residue, atom_name, segment_id);
00718 
00719   for (size_t i = 0; i < atoms_ids.size(); i++) {
00720     if (atoms_ids[i] == aid) {
00721       // this atom id was already recorded
00722       atoms_ncopies[i] += 1;
00723       return i;
00724     }
00725   }
00726 
00727   if (cvm::debug())
00728     cvm::log("Adding atom \""+
00729         atom_name+"\" in residue "+
00730         cvm::to_str(residue)+
00731         " (index "+cvm::to_str(aid)+
00732         ") for collective variables calculation.\n");
00733 
00734   int const index = add_atom_slot(aid);
00735   modifyRequestedAtoms().add(aid);
00736   update_atom_properties(index);
00737   return index;
00738 }

int colvarproxy_namd::init_atom ( int  atom_number  ) 

Definition at line 653 of file colvarproxy_namd.C.

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

00654 {
00655   // save time by checking first whether this atom has been requested before
00656   // (this is more common than a non-valid atom number)
00657   int aid = (atom_number-1);
00658 
00659   for (size_t i = 0; i < atoms_ids.size(); i++) {
00660     if (atoms_ids[i] == aid) {
00661       // this atom id was already recorded
00662       atoms_ncopies[i] += 1;
00663       return i;
00664     }
00665   }
00666 
00667   aid = check_atom_id(atom_number);
00668 
00669   if (aid < 0) {
00670     return INPUT_ERROR;
00671   }
00672 
00673   int const index = add_atom_slot(aid);
00674   modifyRequestedAtoms().add(aid);
00675   update_atom_properties(index);
00676   return index;
00677 }

int colvarproxy_namd::init_atom_group ( std::vector< int > const &  atoms_ids  ) 

Definition at line 1106 of file colvarproxy_namd.C.

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

01107 {
01108   if (cvm::debug())
01109     cvm::log("Reguesting from NAMD a group of size "+cvm::to_str(atoms_ids.size())+
01110         " for collective variables calculation.\n");
01111 
01112   // Note: modifyRequestedGroups is supposed to be in sync with the colvarproxy arrays,
01113   // and to stay that way during a simulation
01114 
01115   // compare this new group to those already allocated inside GlobalMaster
01116   int ig;
01117   for (ig = 0; ig < modifyRequestedGroups().size(); ig++) {
01118     AtomIDList const &namd_group = modifyRequestedGroups()[ig];
01119     bool b_match = true;
01120 
01121     if (namd_group.size() != ((int) atoms_ids.size())) {
01122       b_match = false;
01123     } else {
01124       int ia;
01125       for (ia = 0; ia < namd_group.size(); ia++) {
01126         int const aid = atoms_ids[ia];
01127         if (namd_group[ia] != aid) {
01128           b_match = false;
01129           break;
01130         }
01131       }
01132     }
01133 
01134     if (b_match) {
01135       if (cvm::debug())
01136         cvm::log("Group was already added.\n");
01137       // this group already exists
01138       atom_groups_ncopies[ig] += 1;
01139       return ig;
01140     }
01141   }
01142 
01143   // add this group (note: the argument of add_atom_group_slot() is redundant for NAMD, and provided only for consistency)
01144   size_t const index = add_atom_group_slot(atom_groups_ids.size());
01145   modifyRequestedGroups().resize(atom_groups_ids.size());
01146   // the following is done in calculate()
01147   // modifyGroupForces().resize(atom_groups_ids.size());
01148   AtomIDList &namd_group = modifyRequestedGroups()[index];
01149   namd_group.resize(atoms_ids.size());
01150   int const n_all_atoms = Node::Object()->molecule->numAtoms;
01151   for (size_t ia = 0; ia < atoms_ids.size(); ia++) {
01152     int const aid = atoms_ids[ia];
01153     if (cvm::debug())
01154       cvm::log("Adding atom "+cvm::to_str(aid+1)+
01155           " for collective variables calculation.\n");
01156     if ( (aid < 0) || (aid >= n_all_atoms) ) {
01157       cvm::error("Error: invalid atom number specified, "+
01158                  cvm::to_str(aid+1)+"\n", INPUT_ERROR);
01159       return -1;
01160     }
01161     namd_group[ia] = aid;
01162   }
01163 
01164   update_group_properties(index);
01165 
01166   if (cvm::debug()) {
01167     cvm::log("Group has index "+cvm::to_str(index)+"\n");
01168     cvm::log("modifyRequestedGroups length = "+cvm::to_str(modifyRequestedGroups().size())+
01169         ", modifyGroupForces length = "+cvm::to_str(modifyGroupForces().size())+"\n");
01170   }
01171 
01172   return index;
01173 }

void colvarproxy_namd::init_tcl_pointers (  )  [protected]

Definition at line 547 of file colvarproxy_namd.C.

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

Referenced by colvarproxy_namd().

00548 {
00549 #ifdef NAMD_TCL
00550   // Store pointer to NAMD's Tcl interpreter
00551   tcl_interp_ = reinterpret_cast<void *>(Node::Object()->getScript()->interp);
00552 #endif
00553 }

int colvarproxy_namd::init_volmap ( const char *  volmap_name  ) 

Definition at line 1244 of file colvarproxy_namd.C.

References error(), 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.

01245 {
01246   if (volmap_name == NULL) {
01247     return cvm::error("Error: no grid object name provided.", INPUT_ERROR);
01248   }
01249   int volmap_id = simparams->mgridforcelist.index_for_key(volmap_name);
01250 
01251   int error_code = init_volmap(volmap_id);
01252 
01253   if (error_code == COLVARS_OK) {
01254     // Check that the scale factor is correctly set to zero
01255     Molecule *mol = Node::Object()->molecule;
01256     GridforceGrid const *grid = mol->get_gridfrc_grid(volmap_id);
01257     Vector const gfScale = grid->get_scale();
01258     if ((gfScale.x != 0.0) || (gfScale.y != 0.0) || (gfScale.z != 0.0)) {
01259       error_code |= cvm::error("Error: GridForce map \""+
01260                                std::string(volmap_name)+
01261                                "\" has non-zero scale factors.\n",
01262                                INPUT_ERROR);
01263     }
01264   }
01265 
01266   return error_code;
01267 }

int colvarproxy_namd::init_volmap ( int  volmap_id  ) 

Definition at line 1220 of file colvarproxy_namd.C.

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

Referenced by init_volmap().

01221 {
01222   for (size_t i = 0; i < volmaps_ids.size(); i++) {
01223     if (volmaps_ids[i] == volmap_id) {
01224       // this map has already been requested
01225       volmaps_ncopies[i] += 1;
01226       return i;
01227     }
01228   }
01229 
01230   Molecule *mol = Node::Object()->molecule;
01231 
01232   if ((volmap_id < 0) || (volmap_id >= mol->numGridforceGrids)) {
01233     return cvm::error("Error: invalid numeric ID ("+cvm::to_str(volmap_id)+
01234                       ") for map.\n", INPUT_ERROR);
01235   }
01236 
01237   int const index = add_volmap_slot(volmap_id);
01238   modifyRequestedGridObjects().add(volmap_id);
01239 
01240   return index;
01241 }

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 943 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().

00947 {
00948   if (pdb_field_str.size() == 0)
00949     cvm::error("Error: must define which PDB field to use "
00950                "in order to define atoms from a PDB file.\n", INPUT_ERROR);
00951 
00952   PDB *pdb = new PDB(pdb_filename);
00953   size_t const pdb_natoms = pdb->num_atoms();
00954 
00955   e_pdb_field pdb_field_index = pdb_field_str2enum(pdb_field_str);
00956 
00957   for (size_t ipdb = 0; ipdb < pdb_natoms; ipdb++) {
00958 
00959     double atom_pdb_field_value = 0.0;
00960 
00961     switch (pdb_field_index) {
00962     case e_pdb_occ:
00963       atom_pdb_field_value = (pdb->atom(ipdb))->occupancy();
00964       break;
00965     case e_pdb_beta:
00966       atom_pdb_field_value = (pdb->atom(ipdb))->temperaturefactor();
00967       break;
00968     case e_pdb_x:
00969       atom_pdb_field_value = (pdb->atom(ipdb))->xcoor();
00970       break;
00971     case e_pdb_y:
00972       atom_pdb_field_value = (pdb->atom(ipdb))->ycoor();
00973       break;
00974     case e_pdb_z:
00975       atom_pdb_field_value = (pdb->atom(ipdb))->zcoor();
00976       break;
00977     default:
00978       break;
00979     }
00980 
00981     if ( (pdb_field_value) &&
00982          (atom_pdb_field_value != pdb_field_value) ) {
00983       continue;
00984     } else if (atom_pdb_field_value == 0.0) {
00985       continue;
00986     }
00987 
00988     if (atoms.is_enabled(colvardeps::f_ag_scalable)) {
00989       atoms.add_atom_id(ipdb);
00990     } else {
00991       atoms.add_atom(cvm::atom(ipdb+1));
00992     }
00993   }
00994 
00995   delete pdb;
00996   return (cvm::get_error() ? COLVARS_ERROR : COLVARS_OK);
00997 }

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 831 of file colvarproxy_namd.C.

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

00836 {
00837   if (pdb_field_str.size() == 0 && indices.size() == 0) {
00838     cvm::error("Bug alert: either PDB field should be defined or list of "
00839                "atom IDs should be available when loading atom coordinates!\n", BUG_ERROR);
00840   }
00841 
00842   e_pdb_field pdb_field_index;
00843   bool const use_pdb_field = (pdb_field_str.size() > 0);
00844   if (use_pdb_field) {
00845     pdb_field_index = pdb_field_str2enum(pdb_field_str);
00846   }
00847 
00848   // next index to be looked up in PDB file (if list is supplied)
00849   std::vector<int>::const_iterator current_index = indices.begin();
00850 
00851   PDB *pdb = new PDB(pdb_filename);
00852   size_t const pdb_natoms = pdb->num_atoms();
00853 
00854   if (pos.size() != pdb_natoms) {
00855 
00856     bool const pos_allocated = (pos.size() > 0);
00857 
00858     size_t ipos = 0, ipdb = 0;
00859     for ( ; ipdb < pdb_natoms; ipdb++) {
00860 
00861       if (use_pdb_field) {
00862         // PDB field mode: skip atoms with wrong value in PDB field
00863         double atom_pdb_field_value = 0.0;
00864 
00865         switch (pdb_field_index) {
00866         case e_pdb_occ:
00867           atom_pdb_field_value = (pdb->atom(ipdb))->occupancy();
00868           break;
00869         case e_pdb_beta:
00870           atom_pdb_field_value = (pdb->atom(ipdb))->temperaturefactor();
00871           break;
00872         case e_pdb_x:
00873           atom_pdb_field_value = (pdb->atom(ipdb))->xcoor();
00874           break;
00875         case e_pdb_y:
00876           atom_pdb_field_value = (pdb->atom(ipdb))->ycoor();
00877           break;
00878         case e_pdb_z:
00879           atom_pdb_field_value = (pdb->atom(ipdb))->zcoor();
00880           break;
00881         default:
00882           break;
00883         }
00884 
00885         if ( (pdb_field_value) &&
00886              (atom_pdb_field_value != pdb_field_value) ) {
00887           continue;
00888         } else if (atom_pdb_field_value == 0.0) {
00889           continue;
00890         }
00891 
00892       } else {
00893         // Atom ID mode: use predefined atom IDs from the atom group
00894         if (((int) ipdb) != *current_index) {
00895           // Skip atoms not in the list
00896           continue;
00897         } else {
00898           current_index++;
00899         }
00900       }
00901 
00902       if (!pos_allocated) {
00903         pos.push_back(cvm::atom_pos(0.0, 0.0, 0.0));
00904       } else if (ipos >= pos.size()) {
00905         cvm::error("Error: the PDB file \""+
00906                    std::string(pdb_filename)+
00907                    "\" contains coordinates for "
00908                    "more atoms than needed.\n", BUG_ERROR);
00909       }
00910 
00911       pos[ipos] = cvm::atom_pos((pdb->atom(ipdb))->xcoor(),
00912                                 (pdb->atom(ipdb))->ycoor(),
00913                                 (pdb->atom(ipdb))->zcoor());
00914       ipos++;
00915       if (!use_pdb_field && current_index == indices.end())
00916         break;
00917     }
00918 
00919     if (ipos < pos.size() || (!use_pdb_field && current_index != indices.end())) {
00920       size_t n_requested = use_pdb_field ? pos.size() : indices.size();
00921       cvm::error("Error: number of matching records in the PDB file \""+
00922                  std::string(pdb_filename)+"\" ("+cvm::to_str(ipos)+
00923                  ") does not match the number of requested coordinates ("+
00924                  cvm::to_str(n_requested)+").\n", INPUT_ERROR);
00925       return COLVARS_ERROR;
00926     }
00927   } else {
00928 
00929     // when the PDB contains exactly the number of atoms of the array,
00930     // ignore the fields and just read coordinates
00931     for (size_t ia = 0; ia < pos.size(); ia++) {
00932       pos[ia] = cvm::atom_pos((pdb->atom(ia))->xcoor(),
00933                               (pdb->atom(ia))->ycoor(),
00934                               (pdb->atom(ia))->zcoor());
00935     }
00936   }
00937 
00938   delete pdb;
00939   return COLVARS_OK;
00940 }

void colvarproxy_namd::log ( std::string const &  message  ) 

Definition at line 596 of file colvarproxy_namd.C.

References endi(), and iout.

Referenced by calculate(), colvarproxy_namd(), error(), exit(), init_atom_group(), output_stream(), request_total_force(), script_obj_to_str_vector(), setup(), update_atom_properties(), update_atoms_map(), and update_group_properties().

00597 {
00598   std::istringstream is(message);
00599   std::string line;
00600   while (std::getline(is, line))
00601     iout << "colvars: " << line << "\n";
00602   iout << endi;
00603 }

int colvarproxy_namd::num_replicas (  )  [virtual]

Definition at line 1377 of file colvarproxy_namd.C.

01377                                    {
01378   return CmiNumPartitions();
01379 }

std::ostream * colvarproxy_namd::output_stream ( std::string const &  output_name,
std::ios_base::openmode  mode 
)

Definition at line 1000 of file colvarproxy_namd.C.

References backup_file(), debug, error(), ofstream_namd::is_open(), and log().

01002 {
01003   if (cvm::debug()) {
01004     cvm::log("Using colvarproxy_namd::output_stream()\n");
01005   }
01006 
01007   std::ostream *os = get_output_stream(output_name);
01008   if (os != NULL) return os;
01009 
01010   if (!(mode & (std::ios_base::app | std::ios_base::ate))) {
01011     colvarproxy::backup_file(output_name);
01012   }
01013   ofstream_namd *osf = new ofstream_namd(output_name.c_str(), mode);
01014   if (!osf->is_open()) {
01015     cvm::error("Error: cannot write to file \""+output_name+"\".\n",
01016                FILE_ERROR);
01017     return NULL;
01018   }
01019   output_stream_names.push_back(output_name);
01020   output_files.push_back(osf);
01021   return osf;
01022 }

cvm::rvector colvarproxy_namd::position_distance ( cvm::atom_pos const &  pos1,
cvm::atom_pos const &  pos2 
) const

Definition at line 764 of file colvarproxy_namd.C.

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

00767 {
00768   Position const p1(pos1.x, pos1.y, pos1.z);
00769   Position const p2(pos2.x, pos2.y, pos2.z);
00770   // return p2 - p1
00771   if (this->lattice != NULL) {
00772     Vector const d = this->lattice->delta(p2, p1);
00773     return cvm::rvector(d.x, d.y, d.z);
00774   }
00775   else {
00776     return colvarproxy_system::position_distance(pos1, pos2);
00777   }
00778 }

cvm::real colvarproxy_namd::rand_gaussian (  )  [inline]

Definition at line 105 of file colvarproxy_namd.h.

References Random::gaussian(), and random.

00106   {
00107     return random.gaussian();
00108   }

void colvarproxy_namd::replica_comm_barrier (  )  [virtual]

Definition at line 1382 of file colvarproxy_namd.C.

References replica_barrier().

01382                                             {
01383   replica_barrier();
01384 }

int colvarproxy_namd::replica_comm_recv ( char *  msg_data,
int  buf_len,
int  src_rep 
) [virtual]

Definition at line 1387 of file colvarproxy_namd.C.

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

01388                                                      {
01389   DataMessage *recvMsg = NULL;
01390   replica_recv(&recvMsg, src_rep, CkMyPe());
01391   CmiAssert(recvMsg != NULL);
01392   int retval = recvMsg->size;
01393   if (buf_len >= retval) {
01394     memcpy(msg_data,recvMsg->data,retval);
01395   } else {
01396     retval = 0;
01397   }
01398   CmiFree(recvMsg);
01399   return retval;
01400 }

int colvarproxy_namd::replica_comm_send ( char *  msg_data,
int  msg_len,
int  dest_rep 
) [virtual]

Definition at line 1403 of file colvarproxy_namd.C.

References replica_send().

01404                                                       {
01405   replica_send(msg_data, msg_len, dest_rep, CkMyPe());
01406   return msg_len;
01407 }

int colvarproxy_namd::replica_enabled (  )  [virtual]

Definition at line 1363 of file colvarproxy_namd.C.

01363                                       {
01364 #if CMK_HAS_PARTITION
01365   return COLVARS_OK;
01366 #else
01367   return COLVARS_NOT_IMPLEMENTED;
01368 #endif
01369 }

int colvarproxy_namd::replica_index (  )  [virtual]

Definition at line 1372 of file colvarproxy_namd.C.

01372                                     {
01373   return CmiMyPartition();
01374 }

void colvarproxy_namd::request_total_force ( bool  yesno  ) 

Definition at line 583 of file colvarproxy_namd.C.

References debug, log(), and GlobalMaster::requestTotalForce().

00584 {
00585   if (cvm::debug()) {
00586     cvm::log("colvarproxy_namd::request_total_force()\n");
00587   }
00588   total_force_requested = yesno;
00589   requestTotalForce(total_force_requested);
00590   if (cvm::debug()) {
00591     cvm::log("colvarproxy_namd::request_total_force() end\n");
00592   }
00593 }

int colvarproxy_namd::reset (  ) 

Definition at line 256 of file colvarproxy_namd.C.

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

00257 {
00258   int error_code = COLVARS_OK;
00259 
00260   // Unrequest all atoms and group from NAMD
00261   modifyRequestedAtoms().clear();
00262   modifyRequestedGroups().clear();
00263   modifyRequestedGridObjects().clear();
00264 
00265   atoms_map.clear();
00266 
00267   // Clear internal Proxy records
00268   error_code |= colvarproxy::reset();
00269 
00270   return error_code;
00271 }

int colvarproxy_namd::run_colvar_callback ( std::string const &  name,
std::vector< const colvarvalue * > const &  cvcs,
colvarvalue &  value 
)

Definition at line 560 of file colvarproxy_namd.C.

00564 {
00565   return colvarproxy::tcl_run_colvar_callback(name, cvc_values, value);
00566 }

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 568 of file colvarproxy_namd.C.

00572 {
00573   return colvarproxy::tcl_run_colvar_gradient_callback(name, cvc_values,
00574                                                        gradient);
00575 }

int colvarproxy_namd::run_force_callback (  ) 

Definition at line 555 of file colvarproxy_namd.C.

00556 {
00557   return colvarproxy::tcl_run_force_callback();
00558 }

int colvarproxy_namd::scalable_group_coms (  )  [inline]

Definition at line 203 of file colvarproxy_namd.h.

00204   {
00205     return COLVARS_OK;
00206   }

char const * colvarproxy_namd::script_obj_to_str ( unsigned char *  obj  ) 

Definition at line 1069 of file colvarproxy_namd.C.

01070 {
01071 #ifdef NAMD_TCL
01072   return colvarproxy_tcl::tcl_get_str(obj);
01073 #else
01074   // This is most likely not going to be executed
01075   return colvarproxy::script_obj_to_str(obj);
01076 #endif
01077 }

std::vector< std::string > colvarproxy_namd::script_obj_to_str_vector ( unsigned char *  obj  ) 

Definition at line 1080 of file colvarproxy_namd.C.

References debug, and log().

01081 {
01082   if (cvm::debug()) {
01083     cvm::log("Called colvarproxy_namd::script_obj_to_str_vector().\n");
01084   }
01085   std::vector<std::string> result;
01086 #ifdef NAMD_TCL
01087   Tcl_Interp *interp = reinterpret_cast<Tcl_Interp *>(tcl_interp_);
01088   Tcl_Obj *tcl_obj = reinterpret_cast<Tcl_Obj *>(obj);
01089   Tcl_Obj **tcl_list_elems = NULL;
01090   int count = 0;
01091   if (Tcl_ListObjGetElements(interp, tcl_obj, &count, &tcl_list_elems) ==
01092       TCL_OK) {
01093     result.reserve(count);
01094     for (int i = 0; i < count; i++) {
01095       result.push_back(Tcl_GetString(tcl_list_elems[i]));
01096     }
01097   } else {
01098     Tcl_SetResult(interp,
01099                   const_cast<char *>("Cannot parse Tcl list."), TCL_STATIC);
01100   }
01101 #endif
01102   return result;
01103 }

int colvarproxy_namd::set_unit_system ( std::string const &  units_in,
bool  check_only 
)

Definition at line 1210 of file colvarproxy_namd.C.

References error().

01211 {
01212   if (units_in != "real") {
01213     cvm::error("Error: Specified unit system \"" + units_in + "\" is unsupported in NAMD. Supported units are \"real\" (A, kcal/mol).\n");
01214     return COLVARS_ERROR;
01215   }
01216   return COLVARS_OK;
01217 }

int colvarproxy_namd::setup (  ) 

Definition at line 200 of file colvarproxy_namd.C.

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

Referenced by calculate().

00201 {
00202   if (colvars->size() == 0) return COLVARS_OK;
00203 
00204   cvm::log("Updating NAMD interface:\n");
00205 
00206   errno = 0;
00207 
00208   if (simparams->wrapAll) {
00209     cvm::log("Warning: enabling wrapAll can lead to inconsistent results "
00210              "for Colvars calculations: please disable wrapAll, "
00211              "as is the default option in NAMD.\n");
00212   }
00213 
00214   cvm::log("updating atomic data ("+cvm::to_str(atoms_ids.size())+" atoms).\n");
00215 
00216   size_t i;
00217   for (i = 0; i < atoms_ids.size(); i++) {
00218     update_atom_properties(i);
00219 
00220     // zero out mutable arrays
00221     atoms_positions[i] = cvm::rvector(0.0, 0.0, 0.0);
00222     atoms_total_forces[i] = cvm::rvector(0.0, 0.0, 0.0);
00223     atoms_new_colvar_forces[i] = cvm::rvector(0.0, 0.0, 0.0);
00224   }
00225 
00226   size_t n_group_atoms = 0;
00227   for (int ig = 0; ig < modifyRequestedGroups().size(); ig++) {
00228     n_group_atoms += modifyRequestedGroups()[ig].size();
00229   }
00230 
00231   cvm::log("updating group data ("+cvm::to_str(atom_groups_ids.size())+
00232            " scalable groups, "+
00233            cvm::to_str(n_group_atoms)+" atoms in total).\n");
00234 
00235   // Note: groupMassBegin, groupMassEnd may be used here, but they won't work for charges
00236   for (int ig = 0; ig < modifyRequestedGroups().size(); ig++) {
00237 
00238     // update mass and charge
00239     update_group_properties(ig);
00240 
00241     atom_groups_coms[ig] = cvm::rvector(0.0, 0.0, 0.0);
00242     atom_groups_total_forces[ig] = cvm::rvector(0.0, 0.0, 0.0);
00243     atom_groups_new_colvar_forces[ig] = cvm::rvector(0.0, 0.0, 0.0);
00244   }
00245 
00246   log("updating grid object data ("+cvm::to_str(volmaps_ids.size())+
00247       " grid objects in total).\n");
00248   for (int imap = 0; imap < modifyGridObjForces().size(); imap++) {
00249     volmaps_new_colvar_forces[imap] = 0.0;
00250   }
00251 
00252   return COLVARS_OK;
00253 }

cvm::real colvarproxy_namd::temperature (  )  [inline]

Definition at line 100 of file colvarproxy_namd.h.

References thermostat_temperature.

00101   {
00102     return thermostat_temperature;
00103   }

bool colvarproxy_namd::total_forces_enabled (  )  const [inline]

Definition at line 77 of file colvarproxy_namd.h.

00078   {
00079     return total_force_requested;
00080   }

void colvarproxy_namd::update_atom_properties ( int  index  ) 

Definition at line 748 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().

00749 {
00750   Molecule *mol = Node::Object()->molecule;
00751   // update mass
00752   double const mass = mol->atommass(atoms_ids[index]);
00753   if (mass <= 0.001) {
00754     this->log("Warning: near-zero mass for atom "+
00755               cvm::to_str(atoms_ids[index]+1)+
00756               "; expect unstable dynamics if you apply forces to it.\n");
00757   }
00758   atoms_masses[index] = mass;
00759   // update charge
00760   atoms_charges[index] = mol->atomcharge(atoms_ids[index]);
00761 }

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, log(), and update_atom_properties().

Referenced by calculate().

00167 {
00168   for (AtomIDList::const_iterator a_i = begin; a_i != end; a_i++) {
00169 
00170     if (cvm::debug()) {
00171       cvm::log("Updating atoms_map for atom ID "+cvm::to_str(*a_i)+"\n");
00172     }
00173 
00174     if (atoms_map[*a_i] >= 0) continue;
00175 
00176     for (size_t i = 0; i < atoms_ids.size(); i++) {
00177       if (atoms_ids[i] == *a_i) {
00178         atoms_map[*a_i] = i;
00179         break;
00180       }
00181     }
00182 
00183     if (atoms_map[*a_i] < 0) {
00184       // this atom is probably managed by another GlobalMaster:
00185       // add it here anyway to avoid having to test for array boundaries at each step
00186       int const index = add_atom_slot(*a_i);
00187       atoms_map[*a_i] = index;
00188       update_atom_properties(index);
00189     }
00190   }
00191 
00192   if (cvm::debug()) {
00193     cvm::log("atoms_map = "+cvm::to_str(atoms_map)+".\n");
00194   }
00195 
00196   return COLVARS_OK;
00197 }

int colvarproxy_namd::update_group_properties ( int  index  ) 

Definition at line 1183 of file colvarproxy_namd.C.

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

Referenced by init_atom_group(), and setup().

01184 {
01185   AtomIDList const &namd_group = modifyRequestedGroups()[index];
01186   if (cvm::debug()) {
01187     cvm::log("Re-calculating total mass and charge for scalable group no. "+cvm::to_str(index+1)+" ("+
01188              cvm::to_str(namd_group.size())+" atoms).\n");
01189   }
01190 
01191   cvm::real total_mass = 0.0;
01192   cvm::real total_charge = 0.0;
01193   for (int i = 0; i < namd_group.size(); i++) {
01194     total_mass += Node::Object()->molecule->atommass(namd_group[i]);
01195     total_charge += Node::Object()->molecule->atomcharge(namd_group[i]);
01196   }
01197   atom_groups_masses[index] = total_mass;
01198   atom_groups_charges[index] = total_charge;
01199 
01200   if (cvm::debug()) {
01201     cvm::log("total mass = "+cvm::to_str(total_mass)+
01202              ", total charge = "+cvm::to_str(total_charge)+"\n");
01203   }
01204 
01205   return COLVARS_OK;
01206 }


Friends And Related Function Documentation

friend class cvm::atom [friend]

Definition at line 56 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 35 of file colvarproxy_namd.h.

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

Definition at line 46 of file colvarproxy_namd.h.

Referenced by calculate(), and colvarproxy_namd().

Definition at line 47 of file colvarproxy_namd.h.

Referenced by calculate().

NAMD-style PRNG object.

Definition at line 44 of file colvarproxy_namd.h.

Referenced by colvarproxy_namd(), and rand_gaussian().

Used to submit restraint energy as MISC.

Definition at line 50 of file colvarproxy_namd.h.

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

Pointer to the NAMD simulation input object.

Definition at line 38 of file colvarproxy_namd.h.

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

Self-explained.

Definition at line 41 of file colvarproxy_namd.h.

Referenced by colvarproxy_namd(), and temperature().


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

Generated on 4 Jun 2020 for NAMD by  doxygen 1.6.1