ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/branches/development/src/primitives/Molecule.cpp
Revision: 1465
Committed: Fri Jul 9 23:08:25 2010 UTC (14 years, 9 months ago) by chuckv
File size: 10145 byte(s)
Log Message:
Creating busticated version of OpenMD

File Contents

# User Rev Content
1 gezelter 507 /*
2 gezelter 246 * Copyright (c) 2005 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 gezelter 1390 * 1. Redistributions of source code must retain the above copyright
10 gezelter 246 * notice, this list of conditions and the following disclaimer.
11     *
12 gezelter 1390 * 2. Redistributions in binary form must reproduce the above copyright
13 gezelter 246 * 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 gezelter 1390 *
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 gezelter 246 */
41    
42     /**
43     * @file Molecule.cpp
44     * @author tlin
45     * @date 10/28/2004
46     * @version 1.0
47     */
48 gezelter 2
49 gezelter 246 #include <algorithm>
50     #include <set>
51 gezelter 2
52 tim 3 #include "primitives/Molecule.hpp"
53 gezelter 246 #include "utils/MemoryUtils.hpp"
54 tim 3 #include "utils/simError.h"
55 gezelter 2
56 gezelter 1390 namespace OpenMD {
57 gezelter 507 Molecule::Molecule(int stampId, int globalIndex, const std::string& molName)
58 gezelter 246 : stampId_(stampId), globalIndex_(globalIndex), moleculeName_(molName) {
59 gezelter 1211 }
60    
61 gezelter 507 Molecule::~Molecule() {
62 gezelter 1211
63 tim 398 MemoryUtils::deletePointers(atoms_);
64     MemoryUtils::deletePointers(bonds_);
65     MemoryUtils::deletePointers(bends_);
66     MemoryUtils::deletePointers(torsions_);
67 gezelter 1277 MemoryUtils::deletePointers(inversions_);
68 tim 398 MemoryUtils::deletePointers(rigidBodies_);
69     MemoryUtils::deletePointers(cutoffGroups_);
70     MemoryUtils::deletePointers(constraintPairs_);
71     MemoryUtils::deletePointers(constraintElems_);
72 gezelter 1211 // integrableObjects_ don't own the objects
73 gezelter 246 integrableObjects_.clear();
74    
75 gezelter 507 }
76 gezelter 1211
77 gezelter 507 void Molecule::addAtom(Atom* atom) {
78 gezelter 246 if (std::find(atoms_.begin(), atoms_.end(), atom) == atoms_.end()) {
79 gezelter 507 atoms_.push_back(atom);
80 gezelter 246 }
81 gezelter 507 }
82 gezelter 1211
83 gezelter 507 void Molecule::addBond(Bond* bond) {
84 gezelter 246 if (std::find(bonds_.begin(), bonds_.end(), bond) == bonds_.end()) {
85 gezelter 507 bonds_.push_back(bond);
86 gezelter 246 }
87 gezelter 507 }
88 gezelter 1211
89 gezelter 507 void Molecule::addBend(Bend* bend) {
90 gezelter 246 if (std::find(bends_.begin(), bends_.end(), bend) == bends_.end()) {
91 gezelter 507 bends_.push_back(bend);
92 gezelter 246 }
93 gezelter 507 }
94 gezelter 1211
95 gezelter 507 void Molecule::addTorsion(Torsion* torsion) {
96 gezelter 1211 if (std::find(torsions_.begin(), torsions_.end(), torsion) ==
97     torsions_.end()) {
98 gezelter 507 torsions_.push_back(torsion);
99 gezelter 246 }
100 gezelter 507 }
101 gezelter 1277
102     void Molecule::addInversion(Inversion* inversion) {
103     if (std::find(inversions_.begin(), inversions_.end(), inversion) ==
104     inversions_.end()) {
105     inversions_.push_back(inversion);
106     }
107     }
108 gezelter 1211
109 gezelter 507 void Molecule::addRigidBody(RigidBody *rb) {
110 gezelter 1211 if (std::find(rigidBodies_.begin(), rigidBodies_.end(), rb) ==
111     rigidBodies_.end()) {
112 gezelter 507 rigidBodies_.push_back(rb);
113 gezelter 246 }
114 gezelter 507 }
115 gezelter 1211
116 gezelter 507 void Molecule::addCutoffGroup(CutoffGroup* cp) {
117 gezelter 1211 if (std::find(cutoffGroups_.begin(), cutoffGroups_.end(), cp) ==
118     cutoffGroups_.end()) {
119 gezelter 507 cutoffGroups_.push_back(cp);
120 gezelter 1211 }
121 gezelter 507 }
122 gezelter 1211
123 gezelter 507 void Molecule::addConstraintPair(ConstraintPair* cp) {
124 gezelter 1211 if (std::find(constraintPairs_.begin(), constraintPairs_.end(), cp) ==
125     constraintPairs_.end()) {
126 gezelter 507 constraintPairs_.push_back(cp);
127 gezelter 1211 }
128 gezelter 507 }
129 gezelter 1211
130 gezelter 507 void Molecule::addConstraintElem(ConstraintElem* cp) {
131 gezelter 1211 if (std::find(constraintElems_.begin(), constraintElems_.end(), cp) ==
132     constraintElems_.end()) {
133 gezelter 507 constraintElems_.push_back(cp);
134 gezelter 246 }
135 gezelter 507 }
136 gezelter 1211
137 gezelter 507 void Molecule::complete() {
138 gezelter 246
139     std::set<Atom*> rigidAtoms;
140     RigidBody* rb;
141     std::vector<RigidBody*>::iterator rbIter;
142 gezelter 2
143 gezelter 246
144     for (rb = beginRigidBody(rbIter); rb != NULL; rb = nextRigidBody(rbIter)) {
145 gezelter 507 rigidAtoms.insert(rb->getBeginAtomIter(), rb->getEndAtomIter());
146 gezelter 246 }
147 gezelter 1211
148 gezelter 246 Atom* atom;
149     AtomIterator ai;
150     for (atom = beginAtom(ai); atom != NULL; atom = nextAtom(ai)) {
151 gezelter 1211
152 gezelter 507 if (rigidAtoms.find(*ai) == rigidAtoms.end()) {
153 gezelter 1211
154     // If an atom does not belong to a rigid body, it is an
155     // integrable object
156    
157 gezelter 507 integrableObjects_.push_back(*ai);
158     }
159 gezelter 246 }
160 gezelter 1211
161 gezelter 246 //find all free atoms (which do not belong to rigid bodies)
162 gezelter 1211 // performs the "difference" operation from set theory, the output
163     // range contains a copy of every element that is contained in
164     // [allAtoms.begin(), allAtoms.end()) and not contained in
165     // [rigidAtoms.begin(), rigidAtoms.end()).
166 gezelter 246 //std::set_difference(allAtoms.begin(), allAtoms.end(), rigidAtoms.begin(), rigidAtoms.end(),
167     // std::back_inserter(integrableObjects_));
168 gezelter 2
169 gezelter 246 //if (integrableObjects_.size() != allAtoms.size() - rigidAtoms.size()) {
170     // //Some atoms in rigidAtoms are not in allAtoms, something must be wrong
171     // sprintf(painCave.errMsg, "Atoms in rigidbody are not in the atom list of the same molecule");
172     //
173     // painCave.isFatal = 1;
174     // simError();
175     //}
176 tim 372 for (rb = beginRigidBody(rbIter); rb != NULL; rb = nextRigidBody(rbIter)) {
177     integrableObjects_.push_back(rb);
178     }
179     //integrableObjects_.insert(integrableObjects_.end(), rigidBodies_.begin(), rigidBodies_.end());
180 gezelter 507 }
181 gezelter 2
182 tim 963 RealType Molecule::getMass() {
183 gezelter 246 StuntDouble* sd;
184     std::vector<StuntDouble*>::iterator i;
185 tim 963 RealType mass = 0.0;
186 gezelter 1211
187     for (sd = beginIntegrableObject(i); sd != NULL; sd =
188     nextIntegrableObject(i)){
189 gezelter 507 mass += sd->getMass();
190 gezelter 246 }
191 gezelter 1211
192     return mass;
193 gezelter 507 }
194 gezelter 2
195 gezelter 507 Vector3d Molecule::getCom() {
196 gezelter 246 StuntDouble* sd;
197     std::vector<StuntDouble*>::iterator i;
198     Vector3d com;
199 tim 963 RealType totalMass = 0;
200     RealType mass;
201 gezelter 246
202 gezelter 1211 for (sd = beginIntegrableObject(i); sd != NULL; sd =
203     nextIntegrableObject(i)){
204 gezelter 507 mass = sd->getMass();
205     totalMass += mass;
206     com += sd->getPos() * mass;
207 gezelter 2 }
208 gezelter 1211
209 gezelter 246 com /= totalMass;
210 gezelter 2
211 gezelter 246 return com;
212 gezelter 507 }
213 gezelter 2
214 gezelter 507 void Molecule::moveCom(const Vector3d& delta) {
215 gezelter 246 StuntDouble* sd;
216     std::vector<StuntDouble*>::iterator i;
217    
218 gezelter 1211 for (sd = beginIntegrableObject(i); sd != NULL; sd =
219     nextIntegrableObject(i)){
220 gezelter 507 sd->setPos(sd->getPos() + delta);
221 gezelter 1211 }
222 gezelter 507 }
223 gezelter 2
224 gezelter 507 Vector3d Molecule::getComVel() {
225 gezelter 246 StuntDouble* sd;
226     std::vector<StuntDouble*>::iterator i;
227     Vector3d velCom;
228 tim 963 RealType totalMass = 0;
229     RealType mass;
230 gezelter 246
231 gezelter 1211 for (sd = beginIntegrableObject(i); sd != NULL; sd =
232     nextIntegrableObject(i)){
233 gezelter 507 mass = sd->getMass();
234     totalMass += mass;
235     velCom += sd->getVel() * mass;
236 gezelter 246 }
237 gezelter 1211
238 gezelter 246 velCom /= totalMass;
239 gezelter 1211
240 gezelter 246 return velCom;
241 gezelter 507 }
242 gezelter 2
243 tim 963 RealType Molecule::getPotential() {
244 gezelter 2
245 gezelter 246 Bond* bond;
246     Bend* bend;
247     Torsion* torsion;
248 gezelter 1277 Inversion* inversion;
249 gezelter 246 Molecule::BondIterator bondIter;;
250     Molecule::BendIterator bendIter;
251     Molecule::TorsionIterator torsionIter;
252 gezelter 1277 Molecule::InversionIterator inversionIter;
253 gezelter 2
254 tim 963 RealType potential = 0.0;
255 gezelter 2
256 gezelter 246 for (bond = beginBond(bondIter); bond != NULL; bond = nextBond(bondIter)) {
257 gezelter 507 potential += bond->getPotential();
258 gezelter 246 }
259 gezelter 2
260 gezelter 246 for (bend = beginBend(bendIter); bend != NULL; bend = nextBend(bendIter)) {
261 gezelter 507 potential += bend->getPotential();
262 gezelter 2 }
263    
264 gezelter 1211 for (torsion = beginTorsion(torsionIter); torsion != NULL; torsion =
265     nextTorsion(torsionIter)) {
266 gezelter 507 potential += torsion->getPotential();
267 gezelter 2 }
268 gezelter 1211
269 gezelter 1277 for (inversion = beginInversion(inversionIter); torsion != NULL;
270     inversion = nextInversion(inversionIter)) {
271     potential += inversion->getPotential();
272     }
273    
274 gezelter 246 return potential;
275 gezelter 1211
276 gezelter 507 }
277 gezelter 1211
278 cli2 1360 void Molecule::addProperty(GenericData* genData) {
279     properties_.addProperty(genData);
280     }
281    
282     void Molecule::removeProperty(const std::string& propName) {
283     properties_.removeProperty(propName);
284     }
285    
286     void Molecule::clearProperties() {
287     properties_.clearProperties();
288     }
289    
290     std::vector<std::string> Molecule::getPropertyNames() {
291     return properties_.getPropertyNames();
292     }
293    
294     std::vector<GenericData*> Molecule::getProperties() {
295     return properties_.getProperties();
296     }
297    
298     GenericData* Molecule::getPropertyByName(const std::string& propName) {
299     return properties_.getPropertyByName(propName);
300     }
301    
302    
303    
304    
305 gezelter 507 std::ostream& operator <<(std::ostream& o, Molecule& mol) {
306 gezelter 246 o << std::endl;
307     o << "Molecule " << mol.getGlobalIndex() << "has: " << std::endl;
308     o << mol.getNAtoms() << " atoms" << std::endl;
309     o << mol.getNBonds() << " bonds" << std::endl;
310     o << mol.getNBends() << " bends" << std::endl;
311     o << mol.getNTorsions() << " torsions" << std::endl;
312 gezelter 1277 o << mol.getNInversions() << " inversions" << std::endl;
313 gezelter 246 o << mol.getNRigidBodies() << " rigid bodies" << std::endl;
314     o << mol.getNIntegrableObjects() << "integrable objects" << std::endl;
315     o << mol.getNCutoffGroups() << "cutoff groups" << std::endl;
316     o << mol.getNConstraintPairs() << "constraint pairs" << std::endl;
317     return o;
318 gezelter 507 }
319 gezelter 1211
320 gezelter 1390 }//end namespace OpenMD

Properties

Name Value
svn:keywords Author Id Revision Date