ViewVC Help

View File | Revision Log | Show Annotations | View Changeset | Root Listing

root/OpenMD/trunk/src/restraints/ThermoIntegrationForceManager.cpp

Comparing trunk/src/restraints/ThermoIntegrationForceManager.cpp (file contents):
Revision 507 by gezelter, Fri Apr 15 22:04:00 2005 UTC vs.
Revision 1969 by gezelter, Wed Feb 26 14:14:50 2014 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, 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	<	#include <cmath>
43	>	#ifdef IS_MPI
44	>	#include <mpi.h>
45	>	#endif
46	>
47		#include "restraints/ThermoIntegrationForceManager.hpp"
44	–	#include "integrators/Integrator.hpp"
45	–	#include "math/SquareMatrix3.hpp"
46	–	#include "primitives/Molecule.hpp"
47	–	#include "utils/simError.h"
48	–	#include "utils/OOPSEConstant.hpp"
49	–	#include "utils/StringUtils.hpp"
48
49	<	namespace oopse {
49	>	namespace OpenMD {
50
51		ThermoIntegrationForceManager::ThermoIntegrationForceManager(SimInfo* info):
52	<	ForceManager(info){
53	<	currSnapshot_ = info_->getSnapshotManager()->getCurrentSnapshot();
54	<	simParam = info_->getSimParams();
52	>	RestraintForceManager(info){
53	>	currSnapshot_ = info_->getSnapshotManager()->getCurrentSnapshot();
54	>	simParam = info_->getSimParams();
55
56	<	if (simParam->haveThermIntLambda()){
57	<	tIntLambda_ = simParam->getThermIntLambda();
58	<	}
59	<	else{
60	<	tIntLambda_ = 1.0;
61	<	sprintf(painCave.errMsg,
62	<	"ThermoIntegration error: the transformation parameter (lambda) was\n"
63	<	"\tnot specified. OOPSE will use a default value of %f. To set\n"
64	<	"\tlambda, use the thermodynamicIntegrationLambda variable.\n",
65	<	tIntLambda_);
66	<	painCave.isFatal = 0;
67	<	simError();
68	<	}
56	>	if (simParam->haveThermodynamicIntegrationLambda()){
57	>	tIntLambda_ = simParam->getThermodynamicIntegrationLambda();
58	>	}
59	>	else{
60	>	tIntLambda_ = 1.0;
61	>	sprintf(painCave.errMsg,
62	>	"ThermoIntegration error: the transformation parameter\n"
63	>	"\t(lambda) was not specified. OpenMD will use a default\n"
64	>	"\tvalue of %f. To set lambda, use the \n"
65	>	"\tthermodynamicIntegrationLambda variable.\n",
66	>	tIntLambda_);
67	>	painCave.isFatal = 0;
68	>	simError();
69	>	}
70
71	<	if (simParam->haveThermIntK()){
72	<	tIntK_ = simParam->getThermIntK();
74	<	}
75	<	else{
76	<	tIntK_ = 1.0;
77	<	sprintf(painCave.errMsg,
78	<	"ThermoIntegration Warning: the tranformation parameter exponent\n"
79	<	"\t(k) was not specified. OOPSE will use a default value of %f.\n"
80	<	"\tTo set k, use the thermodynamicIntegrationK variable.\n",
81	<	tIntK_);
82	<	painCave.isFatal = 0;
83	<	simError();
84	<	}
85	<
86	<	if (simParam->getUseSolidThermInt()) {
87	<	// build a restraint object
88	<	restraint_ =  new Restraints(info_, tIntLambda_, tIntK_);
89	<
90	<	}
91	<
92	<	// build the scaling factor used to modulate the forces and torques
93	<	factor_ = pow(tIntLambda_, tIntK_);
94	<
71	>	if (simParam->haveThermodynamicIntegrationK()){
72	>	tIntK_ = simParam->getThermodynamicIntegrationK();
73		}
74	+	else{
75	+	tIntK_ = 1.0;
76	+	sprintf(painCave.errMsg,
77	+	"ThermoIntegration Warning: the tranformation parameter\n"
78	+	"\texponent (k) was not specified. OpenMD will use a default\n"
79	+	"\tvalue of %f. To set k, use the thermodynamicIntegrationK\n"
80	+	"\tvariable.\n",
81	+	tIntK_);
82	+	painCave.isFatal = 0;
83	+	simError();
84	+	}
85	+
86	+	// build the scaling factor used to modulate the forces and torques
87	+	factor_ = pow(tIntLambda_, tIntK_);
88	+	}
89
90		ThermoIntegrationForceManager::~ThermoIntegrationForceManager(){
91		}
92
93	<	void ThermoIntegrationForceManager::calcForces(bool needPotential,
101	<	bool needStress){
93	>	void ThermoIntegrationForceManager::calcForces(){
94		Snapshot* curSnapshot;
95		SimInfo::MoleculeIterator mi;
96		Molecule* mol;
97		Molecule::IntegrableObjectIterator ii;
98	<	StuntDouble* integrableObject;
98	>	StuntDouble* sd;
99		Vector3d frc;
100		Vector3d trq;
101		Mat3x3d tempTau;
102
103		// perform the standard calcForces first
104	<	ForceManager::calcForces(needPotential, needStress);
104	>	ForceManager::calcForces();
105
106		curSnapshot = info_->getSnapshotManager()->getCurrentSnapshot();
107
108	<	// now scale forces and torques of all the integrableObjects
108	>	// now scale forces and torques of all the sds
109
110		for (mol = info_->beginMolecule(mi); mol != NULL;
111		mol = info_->nextMolecule(mi)) {
112	<	for (integrableObject = mol->beginIntegrableObject(ii);
113	<	integrableObject != NULL;
114	<	integrableObject = mol->nextIntegrableObject(ii)) {
115	<	frc = integrableObject->getFrc();
112	>
113	>	for (sd = mol->beginIntegrableObject(ii); sd != NULL;
114	>	sd = mol->nextIntegrableObject(ii)) {
115	>
116	>	frc = sd->getFrc();
117		frc *= factor_;
118	<	integrableObject->setFrc(frc);
118	>	sd->setFrc(frc);
119
120	<	if (integrableObject->isDirectional()){
121	<	trq = integrableObject->getTrq();
120	>	if (sd->isDirectional()){
121	>	trq = sd->getTrq();
122		trq *= factor_;
123	<	integrableObject->setTrq(trq);
123	>	sd->setTrq(trq);
124		}
125		}
126		}
134	–
135	–	// set vraw to be the unmodulated potential
136	–	lrPot_ = curSnapshot->statData[Stats::LONG_RANGE_POTENTIAL];
137	–	curSnapshot->statData[Stats::VRAW] = lrPot_;
127
128	+	// set rawPotential to be the unmodulated potential
129	+	lrPot_ = curSnapshot->getLongRangePotential();
130	+	curSnapshot->setRawPotential(lrPot_);
131	+
132		// modulate the potential and update the snapshot
133		lrPot_ *= factor_;
134	<	curSnapshot->statData[Stats::LONG_RANGE_POTENTIAL] = lrPot_;
134	>	curSnapshot->setLongRangePotential(lrPot_);
135
136		// scale the pressure tensor
137	<	tempTau = curSnapshot->statData.getTau();
137	>	tempTau = curSnapshot->getStressTensor();
138		tempTau *= factor_;
139	<	curSnapshot->statData.setTau(tempTau);
139	>	curSnapshot->setStressTensor(tempTau);
140
141	<	// do crystal restraint forces for thermodynamic integration
142	<	if (simParam->getUseSolidThermInt()) {
143	<
151	<	lrPot_ += restraint_->Calc_Restraint_Forces();
152	<	curSnapshot->statData[Stats::LONG_RANGE_POTENTIAL] = lrPot_;
153	<
154	<	vHarm_ = restraint_->getVharm();
155	<	curSnapshot->statData[Stats::VHARM] = vHarm_;
156	<	}
141	>	// now, on to the applied restraining potentials (if needed):
142	>	RealType restPot_local = 0.0;
143	>	RealType vHarm_local = 0.0;
144
145	<	}
146	<
145	>	if (simParam->getUseRestraints()) {
146	>	// do restraints from RestraintForceManager:
147	>	restPot_local = doRestraints(1.0 - factor_);
148	>	vHarm_local = getUnscaledPotential();
149	>	}
150	>
151	>	#ifdef IS_MPI
152	>	RealType restPot;
153	>	MPI_Allreduce(&restPot_local, &restPot, 1,
154	>	MPI_REALTYPE, MPI_SUM, MPI_COMM_WORLD);
155	>	MPI_Allreduce(&vHarm_local, &vHarm_, 1,
156	>	MPI_REALTYPE, MPI_SUM, MPI_COMM_WORLD);
157	>	lrPot_ += restPot;
158	>	#else
159	>	lrPot_ += restPot_local;
160	>	vHarm_ = vHarm_local;
161	>	#endif
162	>
163	>	// give the final values to stats
164	>	curSnapshot->setLongRangePotential(lrPot_);
165	>	curSnapshot->setRestraintPotential(vHarm_);
166	>	}
167		}

Comparing trunk/src/restraints/ThermoIntegrationForceManager.cpp (property svn:keywords):
Revision 507 by gezelter, Fri Apr 15 22:04:00 2005 UTC vs.
Revision 1969 by gezelter, Wed Feb 26 14:14:50 2014 UTC

#	Line 0 | Line 1
1	+	Author Id Revision Date

Diff Legend

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