ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/Molecule.cpp
(Generate patch)

Comparing trunk/OOPSE/libmdtools/Molecule.cpp (file contents):
Revision 427 by mmeineke, Thu Mar 27 20:48:37 2003 UTC vs.
Revision 1097 by gezelter, Mon Apr 12 20:32:20 2004 UTC

# Line 1 | Line 1
1 < #include <cstdlib>
1 > #include <stdlib.h>
2  
3  
4   #include "Molecule.hpp"
5 + #include "simError.h"
6  
7  
8  
# Line 11 | Line 12 | Molecule::Molecule( void ){
12    myBonds = NULL;
13    myBends = NULL;
14    myTorsions = NULL;
15 +  myRigidBodies = NULL;
16  
17   }
18  
# Line 39 | Line 41 | Molecule::~Molecule( void ){
41      delete[] myTorsions;
42    }
43  
44 <  if( myExcludes != NULL ){
45 <    for(i=0; i<nExcludes; i++) if(myExcludes[i] != NULL ) delete myExcludes[i];
46 <    delete[] myExcludes;
44 >  if( myRigidBodies != NULL ){
45 >    for(i=0; i<nRigidBodies; i++) if(myRigidBodies[i] != NULL )
46 >      delete myRigidBodies[i];
47 >    delete[] myRigidBodies;
48    }
49 +  
50   }
51  
52  
# Line 53 | Line 57 | void Molecule::initialize( molInit &theInit ){
57    nBonds = theInit.nBonds;
58    nBends = theInit.nBends;
59    nTorsions = theInit.nTorsions;
60 <  nExcludes = theInit.nExcludes;
60 >  nRigidBodies = theInit.nRigidBodies;
61    nOriented = theInit.nOriented;
62  
63    myAtoms = theInit.myAtoms;
64    myBonds = theInit.myBonds;
65    myBends = theInit.myBends;
66    myTorsions = theInit.myTorsions;
67 <  myExcludes = theInit.myExcludses;
67 >  myRigidBodies = theInit.myRigidBodies;
68      
69   }
70  
# Line 68 | Line 72 | void Molecule::calcForces( void ){
72    
73    int i;
74  
75 +  for(i=0; i<nRigidBodies; i++) {
76 +    myRigidBodies[i]->updateAtoms();
77 +  }
78 +
79    for(i=0; i<nBonds; i++){
80      myBonds[i]->calc_forces();
81    }
# Line 79 | Line 87 | void Molecule::calcForces( void ){
87    for(i=0; i<nTorsions; i++){
88      myTorsions[i]->calc_forces();
89    }
90 +
91 +  // Rigid Body forces and torques are done after the fortran force loop
92 +
93   }
94  
95  
96 < void Molecule::getPotential( void ){
96 > double Molecule::getPotential( void ){
97    
98    int i;
99    double myPot = 0.0;
# Line 101 | Line 112 | void Molecule::getPotential( void ){
112  
113    return myPot;
114   }
115 +
116 + void Molecule::printMe( void ){
117 +  
118 +  int i;
119 +
120 +  for(i=0; i<nBonds; i++){
121 +    myBonds[i]->printMe();
122 +  }
123 +
124 +  for(i=0; i<nBends; i++){
125 +    myBends[i]->printMe();
126 +  }
127 +
128 +  for(i=0; i<nTorsions; i++){
129 +    myTorsions[i]->printMe();
130 +  }
131 +
132 + }
133 +
134 + void Molecule::moveCOM(double delta[3]){
135 +  double aPos[3];
136 +  int i, j;
137 +
138 +  for(i=0; i<nAtoms; i++) {
139 +    if(myAtoms[i] != NULL ) {
140 +      
141 +      myAtoms[i]->getPos( aPos );
142 +      
143 +      for (j=0; j< 3; j++)
144 +        aPos[j] += delta[j];
145 +
146 +      myAtoms[i]->setPos( aPos );
147 +    }
148 +  }
149 +
150 +  for(i=0; i<nRigidBodies; i++) {
151 +
152 +    if (myRigidBodies[i] != NULL) {
153 +      
154 +      myRigidBodies[i]->getPos( aPos );
155 +
156 +      for (j=0; j< 3; j++)
157 +        aPos[j] += delta[j];
158 +      
159 +      myRigidBodies[i]->setPos( aPos );
160 +    }
161 +  }  
162 + }
163 +
164 + void Molecule::atoms2rigidBodies( void ) {
165 +  int i;
166 +  for (i = 0; i < nRigidBodies; i++) {
167 +    if (myRigidBodies[i] != NULL) {
168 +      myRigidBodies[i]->calcForcesAndTorques();  
169 +    }
170 +  }
171 + }
172 +
173 + void Molecule::getCOM( double COM[3] ) {
174 +
175 +  double mass, mtot;
176 +  double aPos[3];
177 +  int i, j;
178 +
179 +  for (j=0; j<3; j++)
180 +    COM[j] = 0.0;
181 +
182 +  mtot   = 0.0;
183 +
184 +  for (i=0; i < nAtoms; i++) {
185 +    if (myAtoms[i] != NULL) {
186 +
187 +      mass = myAtoms[i]->getMass();
188 +      mtot   += mass;
189 +      
190 +      myAtoms[i]->getPos( aPos );
191 +
192 +      for( j = 0; j < 3; j++)
193 +        COM[j] += aPos[j] * mass;
194 +
195 +    }
196 +  }
197 +
198 +  for (j = 0; j < 3; j++)
199 +    COM[j] /= mtot;
200 + }
201 +
202 + double Molecule::getCOMvel( double COMvel[3] ) {
203 +
204 +  double mass, mtot;
205 +  double aVel[3];
206 +  int i, j;
207 +
208 +
209 +  for (j=0; j<3; j++)
210 +    COMvel[j] = 0.0;
211 +
212 +  mtot   = 0.0;
213 +
214 +  for (i=0; i < nAtoms; i++) {
215 +    if (myAtoms[i] != NULL) {
216 +
217 +      mass = myAtoms[i]->getMass();
218 +      mtot   += mass;
219 +
220 +      myAtoms[i]->getVel(aVel);
221 +
222 +      for (j=0; j<3; j++)
223 +        COMvel[j] += aVel[j]*mass;
224 +
225 +    }
226 +  }
227 +
228 +  for (j=0; j<3; j++)
229 +    COMvel[j] /= mtot;
230 +
231 +  return mtot;
232 +
233 + }
234 +
235 + double Molecule::getTotalMass()
236 + {
237 +  int natoms;
238 +  Atom** atoms;
239 +  double totalMass;
240 +  
241 +  natoms = getNAtoms();
242 +  atoms = getMyAtoms();
243 +  totalMass = 0;
244 +  for(int i =0; i < natoms; i++){
245 +    totalMass += atoms[i]->getMass();
246 +  }
247 +
248 +  return totalMass;
249 + }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines