| 1 | /* | 
| 2 | * Copyright (c) 2009 The University of Notre Dame. All Rights Reserved. | 
| 3 | * | 
| 4 | * The University of Notre Dame grants you ("Licensee") a | 
| 5 | * non-exclusive, royalty free, license to use, modify and | 
| 6 | * redistribute this software in source and binary code form, provided | 
| 7 | * that the following conditions are met: | 
| 8 | * | 
| 9 | * 1. Redistributions of source code must retain the above copyright | 
| 10 | *    notice, this list of conditions and the following disclaimer. | 
| 11 | * | 
| 12 | * 2. Redistributions in binary form must reproduce the above copyright | 
| 13 | *    notice, this list of conditions and the following disclaimer in the | 
| 14 | *    documentation and/or other materials provided with the | 
| 15 | *    distribution. | 
| 16 | * | 
| 17 | * This software is provided "AS IS," without a warranty of any | 
| 18 | * kind. All express or implied conditions, representations and | 
| 19 | * warranties, including any implied warranty of merchantability, | 
| 20 | * fitness for a particular purpose or non-infringement, are hereby | 
| 21 | * excluded.  The University of Notre Dame and its licensors shall not | 
| 22 | * be liable for any damages suffered by licensee as a result of | 
| 23 | * using, modifying or distributing the software or its | 
| 24 | * derivatives. In no event will the University of Notre Dame or its | 
| 25 | * licensors be liable for any lost revenue, profit or data, or for | 
| 26 | * direct, indirect, special, consequential, incidental or punitive | 
| 27 | * damages, however caused and regardless of the theory of liability, | 
| 28 | * arising out of the use of or inability to use software, even if the | 
| 29 | * University of Notre Dame has been advised of the possibility of | 
| 30 | * such damages. | 
| 31 | * | 
| 32 | * SUPPORT OPEN SCIENCE!  If you use OpenMD or its source code in your | 
| 33 | * research, please cite the appropriate papers when you publish your | 
| 34 | * work.  Good starting points are: | 
| 35 | * | 
| 36 | * [1]  Meineke, et al., J. Comp. Chem. 26, 252-271 (2005). | 
| 37 | * [2]  Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006). | 
| 38 | * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008). | 
| 39 | * [4]  Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010). | 
| 40 | * [5]  Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011). | 
| 41 | */ | 
| 42 |  | 
| 43 | #ifdef IS_MPI | 
| 44 | #include <mpi.h> | 
| 45 | #endif | 
| 46 |  | 
| 47 | #include "config.h" | 
| 48 | #include <cmath> | 
| 49 |  | 
| 50 | #include "restraints/RestraintForceManager.hpp" | 
| 51 | #include "restraints/MolecularRestraint.hpp" | 
| 52 | #include "restraints/ObjectRestraint.hpp" | 
| 53 | #include "io/RestReader.hpp" | 
| 54 | #include "utils/simError.h" | 
| 55 | #include "utils/PhysicalConstants.hpp" | 
| 56 | #include "utils/StringUtils.hpp" | 
| 57 | #include "selection/SelectionEvaluator.hpp" | 
| 58 | #include "selection/SelectionManager.hpp" | 
| 59 |  | 
| 60 | namespace OpenMD { | 
| 61 |  | 
| 62 | RestraintForceManager::RestraintForceManager(SimInfo* info): ForceManager(info) { | 
| 63 |  | 
| 64 | // order of affairs: | 
| 65 | // | 
| 66 | // 1) create restraints from the restraintStamps found in the MD | 
| 67 | // file. | 
| 68 | // | 
| 69 | // 2) Create RestraintReader to parse the input files for the ideal | 
| 70 | // structures.  This reader will set reference structures, and will | 
| 71 | // calculate molecular centers of mass, etc. | 
| 72 | // | 
| 73 | // 3) sit around and wait for calcForces to be called.  When it comes, | 
| 74 | // call the normal force manager calcForces, then loop through the | 
| 75 | // restrained objects and do their restraint forces. | 
| 76 |  | 
| 77 | currSnapshot_ = info_->getSnapshotManager()->getCurrentSnapshot(); | 
| 78 | Globals* simParam = info_->getSimParams(); | 
| 79 |  | 
| 80 | if (simParam->haveStatusTime()){ | 
| 81 | restTime_ = simParam->getStatusTime(); | 
| 82 | } else { | 
| 83 | sprintf(painCave.errMsg, | 
| 84 | "Restraint warning: If you use restraints without setting\n" | 
| 85 | "\tstatusTime, no restraint data will be written to the rest\n" | 
| 86 | "\tfile.\n"); | 
| 87 | painCave.isFatal = 0; | 
| 88 | simError(); | 
| 89 | restTime_ = simParam->getRunTime(); | 
| 90 | } | 
| 91 |  | 
| 92 | int nRestraintStamps = simParam->getNRestraintStamps(); | 
| 93 | std::vector<RestraintStamp*> stamp = simParam->getRestraintStamps(); | 
| 94 |  | 
| 95 | std::vector<int> stuntDoubleIndex; | 
| 96 |  | 
| 97 | for (int i = 0; i < nRestraintStamps; i++){ | 
| 98 |  | 
| 99 | std::string myType = toUpperCopy(stamp[i]->getType()); | 
| 100 |  | 
| 101 | if (myType.compare("MOLECULAR")==0){ | 
| 102 |  | 
| 103 | int molIndex; | 
| 104 | Vector3d refCom; | 
| 105 |  | 
| 106 | if (!stamp[i]->haveMolIndex()) { | 
| 107 | sprintf(painCave.errMsg, | 
| 108 | "Restraint Error: A molecular restraint was specified\n" | 
| 109 | "\twithout providing a value for molIndex.\n"); | 
| 110 | painCave.isFatal = 1; | 
| 111 | simError(); | 
| 112 | } else { | 
| 113 | molIndex = stamp[i]->getMolIndex(); | 
| 114 | } | 
| 115 |  | 
| 116 | if (molIndex < 0) { | 
| 117 | sprintf(painCave.errMsg, | 
| 118 | "Restraint Error: A molecular restraint was specified\n" | 
| 119 | "\twith a molIndex that was less than 0\n"); | 
| 120 | painCave.isFatal = 1; | 
| 121 | simError(); | 
| 122 | } | 
| 123 | if (molIndex >= info_->getNGlobalMolecules()) { | 
| 124 | sprintf(painCave.errMsg, | 
| 125 | "Restraint Error: A molecular restraint was specified with\n" | 
| 126 | "\ta molIndex that was greater than the total number of molecules\n"); | 
| 127 | painCave.isFatal = 1; | 
| 128 | simError(); | 
| 129 | } | 
| 130 |  | 
| 131 | Molecule* mol = info_->getMoleculeByGlobalIndex(molIndex); | 
| 132 |  | 
| 133 | if (mol == NULL) { | 
| 134 | #ifdef IS_MPI | 
| 135 | // getMoleculeByGlobalIndex returns a NULL in parallel if | 
| 136 | // this proc doesn't have the molecule.  Do a quick check to | 
| 137 | // make sure another processor is supposed to have it. | 
| 138 |  | 
| 139 | int myrank; | 
| 140 | MPI_Comm_rank( MPI_COMM_WORLD, &myrank); | 
| 141 |  | 
| 142 | if (info_->getMolToProc(molIndex) == myrank) { | 
| 143 |  | 
| 144 | // If we were supposed to have it but got a null, then freak out. | 
| 145 | #endif | 
| 146 |  | 
| 147 | sprintf(painCave.errMsg, | 
| 148 | "Restraint Error: A molecular restraint was specified, but\n" | 
| 149 | "\tno molecule was found with global index %d.\n", | 
| 150 | molIndex); | 
| 151 | painCave.isFatal = 1; | 
| 152 | simError(); | 
| 153 |  | 
| 154 | #ifdef IS_MPI | 
| 155 | } | 
| 156 | #endif | 
| 157 | } | 
| 158 |  | 
| 159 |  | 
| 160 | #ifdef IS_MPI | 
| 161 | // only handle this molecular restraint if this processor owns the | 
| 162 | // molecule | 
| 163 | int myrank; | 
| 164 | MPI_Comm_rank( MPI_COMM_WORLD, &myrank); | 
| 165 | if (info_->getMolToProc(molIndex) == myrank) { | 
| 166 |  | 
| 167 | #endif | 
| 168 |  | 
| 169 | MolecularRestraint* rest = new MolecularRestraint(); | 
| 170 |  | 
| 171 | std::string restPre("mol_"); | 
| 172 | std::stringstream restName; | 
| 173 | restName << restPre << molIndex; | 
| 174 | rest->setRestraintName(restName.str()); | 
| 175 |  | 
| 176 | if (stamp[i]->haveDisplacementSpringConstant()) { | 
| 177 | rest->setDisplacementForceConstant(stamp[i]->getDisplacementSpringConstant()); | 
| 178 | } | 
| 179 | if (stamp[i]->haveTwistSpringConstant()) { | 
| 180 | rest->setTwistForceConstant(stamp[i]->getTwistSpringConstant()); | 
| 181 | } | 
| 182 | if (stamp[i]->haveSwingXSpringConstant()) { | 
| 183 | rest->setSwingXForceConstant(stamp[i]->getSwingXSpringConstant()); | 
| 184 | } | 
| 185 | if (stamp[i]->haveSwingYSpringConstant()) { | 
| 186 | rest->setSwingYForceConstant(stamp[i]->getSwingYSpringConstant()); | 
| 187 | } | 
| 188 | if (stamp[i]->haveRestrainedTwistAngle()) { | 
| 189 | rest->setRestrainedTwistAngle(stamp[i]->getRestrainedTwistAngle() * M_PI/180.0); | 
| 190 | } | 
| 191 | if (stamp[i]->haveRestrainedSwingYAngle()) { | 
| 192 | rest->setRestrainedSwingYAngle(stamp[i]->getRestrainedSwingYAngle() * M_PI/180.0); | 
| 193 | } | 
| 194 | if (stamp[i]->haveRestrainedSwingXAngle()) { | 
| 195 | rest->setRestrainedSwingXAngle(stamp[i]->getRestrainedSwingXAngle() * M_PI/180.0); | 
| 196 | } | 
| 197 | if (stamp[i]->havePrint()) { | 
| 198 | rest->setPrintRestraint(stamp[i]->getPrint()); | 
| 199 | } | 
| 200 |  | 
| 201 | restraints_.push_back(rest); | 
| 202 | mol->addProperty(new RestraintData("Restraint", rest)); | 
| 203 | restrainedMols_.push_back(mol); | 
| 204 | #ifdef IS_MPI | 
| 205 | } | 
| 206 | #endif | 
| 207 | } else if (myType.compare("OBJECT") == 0) { | 
| 208 |  | 
| 209 | std::string objectSelection; | 
| 210 |  | 
| 211 | if (!stamp[i]->haveObjectSelection()) { | 
| 212 | sprintf(painCave.errMsg, | 
| 213 | "Restraint Error: An object restraint was specified\n" | 
| 214 | "\twithout providing a selection script in the\n" | 
| 215 | "\tobjectSelection variable.\n"); | 
| 216 | painCave.isFatal = 1; | 
| 217 | simError(); | 
| 218 | } else { | 
| 219 | objectSelection = stamp[i]->getObjectSelection(); | 
| 220 | } | 
| 221 |  | 
| 222 | SelectionEvaluator evaluator(info); | 
| 223 | SelectionManager seleMan(info); | 
| 224 |  | 
| 225 | evaluator.loadScriptString(objectSelection); | 
| 226 | seleMan.setSelectionSet(evaluator.evaluate()); | 
| 227 | int selectionCount = seleMan.getSelectionCount(); | 
| 228 |  | 
| 229 | #ifdef IS_MPI | 
| 230 | MPI_Allreduce(MPI_IN_PLACE, &selectionCount, 1, MPI_INT, MPI_SUM, | 
| 231 | MPI_COMM_WORLD); | 
| 232 | #endif | 
| 233 |  | 
| 234 | sprintf(painCave.errMsg, | 
| 235 | "Restraint Info: The specified restraint objectSelection,\n" | 
| 236 | "\t\t%s\n" | 
| 237 | "\twill result in %d integrable objects being\n" | 
| 238 | "\trestrained.\n", objectSelection.c_str(), selectionCount); | 
| 239 | painCave.severity = OPENMD_INFO; | 
| 240 | painCave.isFatal = 0; | 
| 241 | simError(); | 
| 242 |  | 
| 243 | int selei; | 
| 244 | StuntDouble* sd; | 
| 245 |  | 
| 246 | for (sd = seleMan.beginSelected(selei); sd != NULL; | 
| 247 | sd = seleMan.nextSelected(selei)) { | 
| 248 | stuntDoubleIndex.push_back(sd->getGlobalIntegrableObjectIndex()); | 
| 249 |  | 
| 250 | ObjectRestraint* rest = new ObjectRestraint(); | 
| 251 |  | 
| 252 | if (stamp[i]->haveDisplacementSpringConstant()) { | 
| 253 | rest->setDisplacementForceConstant(stamp[i]->getDisplacementSpringConstant()); | 
| 254 | } | 
| 255 | if (stamp[i]->haveTwistSpringConstant()) { | 
| 256 | rest->setTwistForceConstant(stamp[i]->getTwistSpringConstant()); | 
| 257 | } | 
| 258 | if (stamp[i]->haveSwingXSpringConstant()) { | 
| 259 | rest->setSwingXForceConstant(stamp[i]->getSwingXSpringConstant()); | 
| 260 | } | 
| 261 | if (stamp[i]->haveSwingYSpringConstant()) { | 
| 262 | rest->setSwingYForceConstant(stamp[i]->getSwingYSpringConstant()); | 
| 263 | } | 
| 264 | if (stamp[i]->haveRestrainedTwistAngle()) { | 
| 265 | rest->setRestrainedTwistAngle(stamp[i]->getRestrainedTwistAngle()); | 
| 266 | } | 
| 267 | if (stamp[i]->haveRestrainedSwingXAngle()) { | 
| 268 | rest->setRestrainedSwingXAngle(stamp[i]->getRestrainedSwingXAngle()); | 
| 269 | } | 
| 270 | if (stamp[i]->haveRestrainedSwingYAngle()) { | 
| 271 | rest->setRestrainedSwingYAngle(stamp[i]->getRestrainedSwingYAngle()); | 
| 272 | } | 
| 273 | if (stamp[i]->havePrint()) { | 
| 274 | rest->setPrintRestraint(stamp[i]->getPrint()); | 
| 275 | } | 
| 276 |  | 
| 277 | restraints_.push_back(rest); | 
| 278 | sd->addProperty(new RestraintData("Restraint", rest)); | 
| 279 | restrainedObjs_.push_back(sd); | 
| 280 | } | 
| 281 |  | 
| 282 | } | 
| 283 | } | 
| 284 |  | 
| 285 | // ThermodynamicIntegration subclasses RestraintForceManager, and there | 
| 286 | // are times when it won't use restraints at all, so only open the | 
| 287 | // restraint file if we are actually using restraints: | 
| 288 |  | 
| 289 | if (simParam->getUseRestraints()) { | 
| 290 | std::string refFile = simParam->getRestraint_file(); | 
| 291 | RestReader* rr = new RestReader(info, refFile, stuntDoubleIndex); | 
| 292 | rr->readReferenceStructure(); | 
| 293 | } | 
| 294 |  | 
| 295 | restOutput_ = getPrefix(info_->getFinalConfigFileName()) + ".rest"; | 
| 296 | restOut = new RestWriter(info_, restOutput_.c_str(), restraints_); | 
| 297 | if(!restOut){ | 
| 298 | sprintf(painCave.errMsg, "Restraint error: Failed to create RestWriter\n"); | 
| 299 | painCave.isFatal = 1; | 
| 300 | simError(); | 
| 301 | } | 
| 302 |  | 
| 303 | // todo: figure out the scale factor.  Right now, just scale them all to 1 | 
| 304 | std::vector<Restraint*>::const_iterator resti; | 
| 305 | for(resti=restraints_.begin(); resti != restraints_.end(); ++resti){ | 
| 306 | (*resti)->setScaleFactor(1.0); | 
| 307 | } | 
| 308 | } | 
| 309 |  | 
| 310 | RestraintForceManager::~RestraintForceManager(){ | 
| 311 | if (restOut) | 
| 312 | delete restOut; | 
| 313 | } | 
| 314 |  | 
| 315 | void RestraintForceManager::init() { | 
| 316 | currRestTime_ = currSnapshot_->getTime(); | 
| 317 | } | 
| 318 |  | 
| 319 | void RestraintForceManager::calcForces(){ | 
| 320 |  | 
| 321 | ForceManager::calcForces(); | 
| 322 | RealType restPot(0.0); | 
| 323 |  | 
| 324 | restPot = doRestraints(1.0); | 
| 325 |  | 
| 326 | #ifdef IS_MPI | 
| 327 | MPI_Allreduce(MPI_IN_PLACE, &restPot, 1, MPI_REALTYPE, MPI_SUM, | 
| 328 | MPI_COMM_WORLD); | 
| 329 | #endif | 
| 330 |  | 
| 331 | currSnapshot_ = info_->getSnapshotManager()->getCurrentSnapshot(); | 
| 332 | RealType pot = currSnapshot_->getLongRangePotential(); | 
| 333 | pot += restPot; | 
| 334 | currSnapshot_->setLongRangePotential(pot); | 
| 335 | currSnapshot_->setRestraintPotential(restPot); | 
| 336 |  | 
| 337 | //write out forces and current positions of restrained molecules | 
| 338 | if (currSnapshot_->getTime() >= currRestTime_){ | 
| 339 | restOut->writeRest(restInfo_); | 
| 340 | currRestTime_ += restTime_; | 
| 341 | } | 
| 342 | } | 
| 343 |  | 
| 344 | RealType RestraintForceManager::doRestraints(RealType scalingFactor){ | 
| 345 | std::vector<Molecule*>::const_iterator rm; | 
| 346 | GenericData* data; | 
| 347 | Molecule::IntegrableObjectIterator ioi; | 
| 348 | MolecularRestraint* mRest; | 
| 349 | StuntDouble* sd; | 
| 350 |  | 
| 351 | std::vector<StuntDouble*>::const_iterator ro; | 
| 352 | ObjectRestraint* oRest; | 
| 353 |  | 
| 354 | std::map<int, Restraint::RealPair> restInfo; | 
| 355 |  | 
| 356 | unscaledPotential_ = 0.0; | 
| 357 |  | 
| 358 | restInfo_.clear(); | 
| 359 |  | 
| 360 | for(rm=restrainedMols_.begin(); rm != restrainedMols_.end(); ++rm){ | 
| 361 |  | 
| 362 | // make sure this molecule (*rm) has a generic data for restraints: | 
| 363 | data = (*rm)->getPropertyByName("Restraint"); | 
| 364 | if (data != NULL) { | 
| 365 | // make sure we can reinterpret the generic data as restraint data: | 
| 366 | RestraintData* restData= dynamic_cast<RestraintData*>(data); | 
| 367 | if (restData != NULL) { | 
| 368 | // make sure we can reinterpet the restraint data as a pointer to | 
| 369 | // an MolecularRestraint: | 
| 370 | mRest = dynamic_cast<MolecularRestraint*>(restData->getData()); | 
| 371 | if (mRest == NULL) { | 
| 372 | sprintf( painCave.errMsg, | 
| 373 | "Can not cast RestraintData to MolecularRestraint\n"); | 
| 374 | painCave.severity = OPENMD_ERROR; | 
| 375 | painCave.isFatal = 1; | 
| 376 | simError(); | 
| 377 | } | 
| 378 | } else { | 
| 379 | sprintf( painCave.errMsg, | 
| 380 | "Can not cast GenericData to RestraintData\n"); | 
| 381 | painCave.severity = OPENMD_ERROR; | 
| 382 | painCave.isFatal = 1; | 
| 383 | simError(); | 
| 384 | } | 
| 385 | } else { | 
| 386 | sprintf( painCave.errMsg, "Can not find Restraint for RestrainedObject\n"); | 
| 387 | painCave.severity = OPENMD_ERROR; | 
| 388 | painCave.isFatal = 1; | 
| 389 | simError(); | 
| 390 | } | 
| 391 |  | 
| 392 | // phew.  At this point, we should have the pointer to the | 
| 393 | // correct MolecularRestraint in the variable mRest. | 
| 394 |  | 
| 395 | Vector3d molCom = (*rm)->getCom(); | 
| 396 |  | 
| 397 | std::vector<Vector3d> struc; | 
| 398 | std::vector<Vector3d> forces; | 
| 399 |  | 
| 400 |  | 
| 401 | for(sd = (*rm)->beginIntegrableObject(ioi); sd != NULL; | 
| 402 | sd = (*rm)->nextIntegrableObject(ioi)) { | 
| 403 | struc.push_back(sd->getPos()); | 
| 404 | } | 
| 405 |  | 
| 406 | mRest->setScaleFactor(scalingFactor); | 
| 407 | mRest->calcForce(struc, molCom); | 
| 408 | forces = mRest->getRestraintForces(); | 
| 409 | int index = 0; | 
| 410 |  | 
| 411 | for(sd = (*rm)->beginIntegrableObject(ioi); sd != NULL; | 
| 412 | sd = (*rm)->nextIntegrableObject(ioi)) { | 
| 413 | sd->addFrc(forces[index]); | 
| 414 | struc.push_back(sd->getPos()); | 
| 415 | index++; | 
| 416 | } | 
| 417 |  | 
| 418 | unscaledPotential_ += mRest->getUnscaledPotential(); | 
| 419 |  | 
| 420 | // only collect data on restraints that we're going to print: | 
| 421 | if (mRest->getPrintRestraint()) | 
| 422 | restInfo = mRest->getRestraintInfo(); | 
| 423 | restInfo_.push_back(restInfo); | 
| 424 | } | 
| 425 |  | 
| 426 | for(ro=restrainedObjs_.begin(); ro != restrainedObjs_.end(); ++ro){ | 
| 427 | // make sure this object (*ro) has a generic data for restraints: | 
| 428 | data = (*ro)->getPropertyByName("Restraint"); | 
| 429 | if (data != NULL) { | 
| 430 | // make sure we can reinterpret the generic data as restraint data: | 
| 431 | RestraintData* restData= dynamic_cast<RestraintData*>(data); | 
| 432 | if (restData != NULL) { | 
| 433 | // make sure we can reinterpet the restraint data as a pointer to | 
| 434 | // an ObjectRestraint: | 
| 435 | oRest = dynamic_cast<ObjectRestraint*>(restData->getData()); | 
| 436 | if (oRest == NULL) { | 
| 437 | sprintf( painCave.errMsg, | 
| 438 | "Can not cast RestraintData to ObjectRestraint\n"); | 
| 439 | painCave.severity = OPENMD_ERROR; | 
| 440 | painCave.isFatal = 1; | 
| 441 | simError(); | 
| 442 | } | 
| 443 | } else { | 
| 444 | sprintf( painCave.errMsg, | 
| 445 | "Can not cast GenericData to RestraintData\n"); | 
| 446 | painCave.severity = OPENMD_ERROR; | 
| 447 | painCave.isFatal = 1; | 
| 448 | simError(); | 
| 449 | } | 
| 450 | } else { | 
| 451 | sprintf( painCave.errMsg, "Can not find Restraint for RestrainedObject\n"); | 
| 452 | painCave.severity = OPENMD_ERROR; | 
| 453 | painCave.isFatal = 1; | 
| 454 | simError(); | 
| 455 | } | 
| 456 |  | 
| 457 | // phew.  At this point, we should have the pointer to the | 
| 458 | // correct Object restraint in the variable oRest. | 
| 459 | oRest->setScaleFactor(scalingFactor); | 
| 460 |  | 
| 461 | Vector3d pos = (*ro)->getPos(); | 
| 462 |  | 
| 463 | if ( (*ro)->isDirectional() ) { | 
| 464 |  | 
| 465 | // directional objects may have orientational restraints as well | 
| 466 | // as positional, so get the rotation matrix first: | 
| 467 |  | 
| 468 | RotMat3x3d A = (*ro)->getA(); | 
| 469 | oRest->calcForce(pos, A); | 
| 470 | (*ro)->addFrc(oRest->getRestraintForce()); | 
| 471 | (*ro)->addTrq(oRest->getRestraintTorque()); | 
| 472 | } else { | 
| 473 |  | 
| 474 | // plain vanilla positional restraints: | 
| 475 |  | 
| 476 | oRest->calcForce(pos); | 
| 477 | (*ro)->addFrc(oRest->getRestraintForce()); | 
| 478 | } | 
| 479 |  | 
| 480 | unscaledPotential_ += oRest->getUnscaledPotential(); | 
| 481 |  | 
| 482 | // only collect data on restraints that we're going to print: | 
| 483 | if (oRest->getPrintRestraint()) | 
| 484 | restInfo = oRest->getRestraintInfo(); | 
| 485 | restInfo_.push_back(restInfo); | 
| 486 | } | 
| 487 |  | 
| 488 | return unscaledPotential_ * scalingFactor; | 
| 489 | } | 
| 490 | } |