ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/src/restraints/RestraintForceManager.cpp
(Generate patch)

Comparing trunk/src/restraints/RestraintForceManager.cpp (file contents):
Revision 1360 by cli2, Mon Sep 7 16:31:51 2009 UTC vs.
Revision 1390 by gezelter, Wed Nov 25 20:02:06 2009 UTC

# Line 6 | Line 6
6   * redistribute this software in source and binary code form, provided
7   * that the following conditions are met:
8   *
9 < * 1. Acknowledgement of the program authors must be made in any
10 < *    publication of scientific results based in part on use of the
11 < *    program.  An acceptable form of acknowledgement is citation of
12 < *    the article in which the program was described (Matthew
13 < *    A. Meineke, Charles F. Vardeman II, Teng Lin, Christopher
14 < *    J. Fennell and J. Daniel Gezelter, "OOPSE: An Object-Oriented
15 < *    Parallel Simulation Engine for Molecular Dynamics,"
16 < *    J. Comput. Chem. 26, pp. 252-271 (2005))
17 < *
18 < * 2. Redistributions of source code must retain the above copyright
9 > * 1. Redistributions of source code must retain the above copyright
10   *    notice, this list of conditions and the following disclaimer.
11   *
12 < * 3. Redistributions in binary form must reproduce the above copyright
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.
# Line 37 | Line 28
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, 24107 (2008).          
39 + * [4]  Vardeman & Gezelter, in progress (2009).                        
40   */
41  
42   #include <cmath>
# Line 45 | Line 45
45   #include "restraints/ObjectRestraint.hpp"
46   #include "io/RestReader.hpp"
47   #include "utils/simError.h"
48 < #include "utils/OOPSEConstant.hpp"
48 > #include "utils/PhysicalConstants.hpp"
49   #include "utils/StringUtils.hpp"
50   #include "selection/SelectionEvaluator.hpp"
51   #include "selection/SelectionManager.hpp"
# Line 54 | Line 54
54   #endif
55  
56  
57 < namespace oopse {
57 > namespace OpenMD {
58  
59    RestraintForceManager::RestraintForceManager(SimInfo* info): ForceManager(info) {
60  
# Line 78 | Line 78 | namespace oopse {
78        restTime_ = simParam->getStatusTime();
79      } else {
80        sprintf(painCave.errMsg,
81 <              "Restraint warning: If you use restraints without setting\n",
82 <              "\tstatusTime, no restraint data will be written to the rest\n",
81 >              "Restraint warning: If you use restraints without setting\n"
82 >              "\tstatusTime, no restraint data will be written to the rest\n"
83                "\tfile.\n");
84        painCave.isFatal = 0;
85        simError();
# Line 109 | Line 109 | namespace oopse {
109            molIndex = stamp[i]->getMolIndex();
110          }
111          
112 <        Molecule* mol = info_->getMoleculeByGlobalIndex(molIndex);
113 <        
114 <        if (mol == NULL) {
112 >        if (molIndex < 0) {
113            sprintf(painCave.errMsg,
114 <                  "Restraint Error: A molecular restraint was specified, but\n"
115 <                  "\tno molecule was found with global index %d.\n",
118 <                  molIndex);
114 >                  "Restraint Error: A molecular restraint was specified\n"
115 >                  "\twith a molIndex that was less than 0\n");
116            painCave.isFatal = 1;
117            simError();      
118 +        }
119 +        if (molIndex >= info_->getNGlobalMolecules()) {
120 +          sprintf(painCave.errMsg,
121 +                  "Restraint Error: A molecular restraint was specified with\n"
122 +                  "\ta molIndex that was greater than the total number of molecules\n");
123 +          painCave.isFatal = 1;
124 +          simError();      
125 +        }
126 +      
127 +        Molecule* mol = info_->getMoleculeByGlobalIndex(molIndex);
128 +        
129 +        if (mol == NULL) {
130 + #ifdef IS_MPI
131 +          // getMoleculeByGlobalIndex returns a NULL in parallel if
132 +          // this proc doesn't have the molecule.  Do a quick check to
133 +          // make sure another processor is supposed to have it.
134 +        
135 +          int myrank = MPI::COMM_WORLD.Get_rank();
136 +          if (info_->getMolToProc(molIndex) == myrank) {
137 +            // If we were supposed to have it but got a null, then freak out.
138 + #endif
139 +
140 +            sprintf(painCave.errMsg,
141 +                    "Restraint Error: A molecular restraint was specified, but\n"
142 +                    "\tno molecule was found with global index %d.\n",
143 +                    molIndex);
144 +            painCave.isFatal = 1;
145 +            simError();    
146 +
147 + #ifdef IS_MPI
148 +          }
149 + #endif
150          }
151          
152          MolecularRestraint* rest = new MolecularRestraint();
# Line 148 | Line 177 | namespace oopse {
177          if (stamp[i]->haveRestrainedSwingXAngle()) {
178            rest->setRestrainedSwingXAngle(stamp[i]->getRestrainedSwingXAngle() * M_PI/180.0);
179          }
180 +        if (stamp[i]->havePrint()) {
181 +          rest->setPrintRestraint(stamp[i]->getPrint());
182 +        }
183  
184          restraints_.push_back(rest);
185          mol->addProperty(new RestraintData("Restraint", rest));
# Line 170 | Line 202 | namespace oopse {
202  
203          SelectionEvaluator evaluator(info);
204          SelectionManager seleMan(info);
205 <        
205 >                
206          evaluator.loadScriptString(objectSelection);
207          seleMan.setSelectionSet(evaluator.evaluate());        
208          int selectionCount = seleMan.getSelectionCount();
# Line 211 | Line 243 | namespace oopse {
243            }
244            if (stamp[i]->haveRestrainedSwingYAngle()) {
245              rest->setRestrainedSwingYAngle(stamp[i]->getRestrainedSwingYAngle());
246 <          }        
246 >          }    
247 >          if (stamp[i]->havePrint()) {
248 >            rest->setPrintRestraint(stamp[i]->getPrint());
249 >          }
250 >    
251            restraints_.push_back(rest);
252            sd->addProperty(new RestraintData("Restraint", rest));
253            restrainedObjs_.push_back(sd);                    
# Line 311 | Line 347 | namespace oopse {
347            if (mRest == NULL) {
348              sprintf( painCave.errMsg,
349                       "Can not cast RestraintData to MolecularRestraint\n");
350 <            painCave.severity = OOPSE_ERROR;
350 >            painCave.severity = OPENMD_ERROR;
351              painCave.isFatal = 1;
352              simError();                      
353            }
354          } else {
355            sprintf( painCave.errMsg,
356                     "Can not cast GenericData to RestraintData\n");
357 <          painCave.severity = OOPSE_ERROR;
357 >          painCave.severity = OPENMD_ERROR;
358            painCave.isFatal = 1;
359            simError();      
360          }
361        } else {
362          sprintf( painCave.errMsg, "Can not find Restraint for RestrainedObject\n");
363 <        painCave.severity = OOPSE_ERROR;
363 >        painCave.severity = OPENMD_ERROR;
364          painCave.isFatal = 1;
365          simError();          
366        }
# Line 357 | Line 393 | namespace oopse {
393        unscaledPotential_ += mRest->getUnscaledPotential();
394  
395        restInfo = mRest->getRestraintInfo();
396 <      restInfo_.push_back(restInfo);
396 >
397 >      // only collect data on restraints that we're going to print:
398 >      if (mRest->getPrintRestraint())
399 >        restInfo_.push_back(restInfo);
400      }
401  
402      for(ro=restrainedObjs_.begin(); ro != restrainedObjs_.end(); ++ro){
# Line 373 | Line 412 | namespace oopse {
412            if (oRest == NULL) {
413              sprintf( painCave.errMsg,
414                       "Can not cast RestraintData to ObjectRestraint\n");
415 <            painCave.severity = OOPSE_ERROR;
415 >            painCave.severity = OPENMD_ERROR;
416              painCave.isFatal = 1;
417              simError();                      
418            }
419          } else {
420            sprintf( painCave.errMsg,
421                     "Can not cast GenericData to RestraintData\n");
422 <          painCave.severity = OOPSE_ERROR;
422 >          painCave.severity = OPENMD_ERROR;
423            painCave.isFatal = 1;
424            simError();      
425          }
426        } else {
427          sprintf( painCave.errMsg, "Can not find Restraint for RestrainedObject\n");
428 <        painCave.severity = OOPSE_ERROR;
428 >        painCave.severity = OPENMD_ERROR;
429          painCave.isFatal = 1;
430          simError();          
431        }
# Line 418 | Line 457 | namespace oopse {
457        unscaledPotential_ += oRest->getUnscaledPotential();
458  
459        restInfo = oRest->getRestraintInfo();
460 <      restInfo_.push_back(restInfo);
460 >
461 >      // only collect data on restraints that we're going to print:
462 >      if (oRest->getPrintRestraint())
463 >        restInfo_.push_back(restInfo);
464      }
465  
466      return unscaledPotential_ * scalingFactor;

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines