ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/branches/new_design/OOPSE-4/src/UseTheForce/DUFF.cpp
(Generate patch)

Comparing branches/new_design/OOPSE-4/src/UseTheForce/DUFF.cpp (file contents):
Revision 1757 by tim, Fri Nov 19 00:11:33 2004 UTC vs.
Revision 1758 by tim, Fri Nov 19 17:56:32 2004 UTC

# Line 25 | Line 25
25  
26   #include "UseTheForce/DUFF.hpp"
27   #include "UseTheForce/DarkSide/lj_interface.h"
28 + #include "UseTheForce/DarkSide/charge_interface.h"
29 + #include "UseTheForce/DarkSide/dipole_interface.h"
30 + #include "UseTheForce/DarkSide/sticky_interface.h"
31  
32   namespace oopse {
33  
# Line 33 | Line 36 | ForceField* createDUFF() {
36      return new DUFF();
37   }
38  
39 < //register createDUFF to ForceFieldCreator
39 > //register createDUFF to ForceFieldFactory
40   ForceFieldFactory::getInstance()->registerForceField("DUFF", createDUFF);
41  
42  
# Line 44 | Line 47 | ParseState DUFF::getSection(const std::string& section
47          case "AtomTypes" :
48              result = DUFF::AtomTypeSection;
49              break;
50 +        case "DirectionalAtomTypes" :
51 +            result = DUFF::DirectionalAtomTypeSection;
52 +            break;
53  
54          case "BondTypes" :
55              result = DUFF::BondTypeSection;
# Line 136 | Line 142 | void DUFF::parseAtomType(const std::string& line, int
142   }
143  
144  
145 < void DUFF::parseAtomType(const std::string& line, int lineNo){
145 > void DUFF::parseAtomType(const std::string& line, int lineNo, int& ident){
146      StringTokenizer tokenizer(line);
147 <    
142 <    if (tokenizer.countTokens() >= 4) {
143 <        atomTypeName = tokenizer.nextToken();
144 <        mass = tokenizer.nextTokenAsFloat();
145 <        epsilon = tokenizer.nextTokenAsFloat();
146 <        sigma = tokenizer.nextTokenAsFloat();
147 >    int nTokens = tokenizer.countTokens();    
148  
149 <        atomType = new AtomType();
149 >    //in AtomTypeSection, a line at least contains 5 tokens
150 >    //atomTypeName, is Directional, isLJ, isCharge and mass
151 >    if (nTokens < 5)  {
152 >        sprintf( painCave.errMsg,
153 >               "Not enough tokens when parsing DUFF Force Field : %s\n"
154 >               "in line %d : %s\n",
155 >               filename.c_str(), lineNo, line);
156 >        painCave.severity = OOPSE_ERROR;
157 >        painCave.isFatal = 1;
158 >        simError();                
159 >
160 >        
161 >    } else {
162 >
163 >        std::string atomTypeName = tokenizer.nextToken();
164 >        bool isDirectional = tokenizer.nextTokenAsBool();
165 >        bool isLJ = tokenizer.nextTokenAsBool();
166 >        bool isCharge = tokenizer.nextTokenAsBool();        
167 >        double mass = tokenizer.nextTokenAsDouble();
168 >        double epsilon;
169 >        double sigma;
170 >        double charge;
171 >        nTokens -= 5;
172 >
173 >        //parse epsilon and sigma
174 >        if (isLJ) {
175 >            if (nTokens >= 2) {
176 >                epsilon = tokenizer.nextTokenAsDouble();
177 >                sigma = tokenizer.nextTokenAsDouble();
178 >                nTokens -= 2;
179 >            } else {
180 >
181 >            }
182 >        }
183 >
184 >        //parse charge
185 >        if (isCharge) {
186 >            if (nTokens >= 1) {
187 >                charge = tokenizer.nextTokenAsDouble();
188 >                nTokens -= 1;
189 >            } else {
190 >
191 >            }
192 >        }
193 >        
194 >        AtomType* atomType;
195 >        if (isDirectional) {
196 >            atomType = new DirectionalAtomType();        
197 >        } else {
198 >            atomType = new AtomType();
199 >        }
200 >        
201          atomType->setName(atomTypeName);
202          atomType->setMass(mass);
203          
# Line 153 | Line 205 | void DUFF::parseAtomType(const std::string& line, int
205          //In  Lennard-Jones Force Field, we only need to set Lennard-Jones to true
206          atomType->setLennardJones();
207  
208 <        atomType->setIdent(ident);
209 <        atomType->addProperty(new DoubleGenericData("Epsilon", epsilon));                                
158 <        atomType->addProperty(new DoubleGenericData("Sigma", sigma));                                
208 >        atomType->setIdent(ident);
209 >        
210          atomType->complete();
160        //notify a new LJtype atom type is created
161        newLJtype(&ident, &sigma, &epsilon, &status);
211  
212 <        //add atom type to AtomTypeContainer
164 <        addAtomType(atomTypeName, atomType);
165 <        ++ident;
212 >        int setLJStatus;
213          
214 <    } else {
215 <        sprintf( painCave.errMsg,
216 <               "Not enough tokens when parsing Lennard-Jones Force Field : %s\n"
217 <               "in line %d : %s\n",
218 <               filename.c_str(), lineNo, line);
219 <        painCave.severity = OOPSE_ERROR;
220 <        painCave.isFatal = 1;
221 <        simError();                
214 >        //notify a new LJtype atom type is created
215 >        if (isLJ) {
216 >            newLJtype(&ident, &sigma, &epsilon, &setLJStatus);
217 >        }
218 >
219 >        int setChargeStatus;
220 >        if (isCharge) {
221 >            newChargeType(&ident, &charge, &setChargeStatus)
222 >        }
223 >
224 >        if (setLJStatus && setChargeStatus) {
225 >            //add atom type to AtomTypeContainer
226 >            addAtomType(atomTypeName, atomType);
227 >            ++ident;
228 >        } else {
229 >            //error in notifying fortran
230 >            delete atomType;
231 >        }
232      }    
233 +    
234   }
235  
236 +
237 + void DUFF::parseDirectionalAtomType(const std::string& line, int lineNo) {
238 +    StringTokenizer tokenizer(line);
239 +    int nTokens = tokenizer.countTokens();  
240 +    
241 +    //in irectionalAtomTypeSection, a line at least contains 6 tokens
242 +    //AtomTypeName, isDipole, isSticky, I_xx, I_yy and I_zz
243 +    if (nTokens < 6) {
244 +
245 +    } else {
246 +
247 +
248 +        std::string atomTypeName = tokenizer.nextToken();
249 +        bool isDipole = tokenizer.nextTokenAsBool();
250 +        bool isSticky = tokenizer.nextTokenAsBool();
251 +        double Ixx = tokenizer.nextTokenAsDouble();
252 +        double Iyy = tokenizer.nextTokenAsDouble();
253 +        double Izz = tokenizer.nextTokenAsDouble();
254 +        nTokens -= 6;
255 +
256 +        AtomType* atomType = getAtomType(atomTypeName);
257 +        if (atomType == NULL) {
258 +
259 +        }
260 +        
261 +        DirectionalAtomType* dAtomType = dynamic_cast<DirectionalAtomType*>(atomType);
262 +        if (dAtomType == NULL) {
263 +
264 +
265 +        }
266 +
267 +        dAtomType->setDipole(isDipole);
268 +        dAtomType->setSticky(isSticky);
269 +
270 +        Mat3x3d inertialMat;
271 +        inertialMat(0, 0) = Ixx;
272 +        inertialMat(1, 1) = Ixx;
273 +        inertialMat(2, 2) = Ixx;        
274 +        dAtomType->setI(inertialMat);
275 +
276 +        //read dipole moment
277 +        double dipole;
278 +        if (isDipole) {
279 +            if (nTokens >= 1) {
280 +                dipole = tokenizer.nextTokenAsDouble();
281 +                nTokens -= 1;
282 +            } else {
283 +
284 +            }
285 +        }
286 +
287 +        //read sticky parameters
288 +        double w0;
289 +        double v0;
290 +        double v0p;
291 +        double rl;
292 +        double ru;
293 +        double rlp;
294 +        double rup;        
295 +        if (isSticky) {
296 +            if (nTokens >= 7) {
297 +                w0 = tokenizer.nextTokenAsDouble();
298 +                v0 = tokenizer.nextTokenAsDouble();
299 +                v0p = tokenizer.nextTokenAsDouble();
300 +                rl = tokenizer.nextTokenAsDouble();
301 +                ru = tokenizer.nextTokenAsDouble();
302 +                rlp = tokenizer.nextTokenAsDouble();
303 +                rup = tokenizer.nextTokenAsDouble();    
304 +                nTokens -= 7;
305 +            } else {
306 +
307 +            }
308 +        }
309 +
310 +
311 +        //notify fotran a newDipoleType is created
312 +        int ident = dAtomType->getIdent();
313 +        int setDipoleStatus;
314 +        if (isDipole) {
315 +            newDipoleType(&ident, &dipole, &setDipoleStatus);
316 +        }              
317 +
318 +        //notify fotran a StickyType is created
319 +        int setStickyStatus;
320 +        if (isSticky) {
321 +            makeStickyType( &w0, &v0, &v0p, &rl, &ru, &rlp, &rup);
322 +        }
323 +
324 +
325 +        if (!setDipoleStatus || !setStickyStatus) {
326 +
327 +        }
328 +        
329 +    }
330 + }
331 +
332   void DUFF::parseBondType(const std::string& line, int lineNo){
333  
334      StringTokenizer tokenizer(line);
# Line 256 | Line 410 | void DUFF::parseBondType(const std::string& line, int
410                  for (int i = 0; i < nPairs; ++i) {
411                      power = tokenizer.nextTokenAsInt();
412                      coefficient = tokenizer.nextTokenAsDouble();
413 <
260 <                    if (power < 0) {
261 <
262 <                    } else {
263 <                        pbt->setCoefficient(power, coefficient);
264 <                    }
413 >                    pbt->setCoefficient(power, coefficient);
414                  }
415              }
416              
# Line 374 | Line 523 | void DUFF::parseBendType(const std::string& line, int
523                  for (int i = 0; i < nPairs; ++i) {
524                      power = tokenizer.nextTokenAsInt();
525                      coefficient = tokenizer.nextTokenAsDouble();
526 <
378 <                    if (power < 0) {
379 <
380 <                    } else {
381 <                        pbt->setCoefficient(power, coefficient);
382 <                    }
526 >                    pbt->setCoefficient(power, coefficient);
527                  }
528              }
529              
# Line 463 | Line 607 | void DUFF::parseTorsionType(const std::string& line, i
607                  for (int i = 0; i < nPairs; ++i) {
608                      power = tokenizer.nextTokenAsInt();
609                      coefficient = tokenizer.nextTokenAsDouble();
610 <
467 <                    if (power < 0) {
468 <
469 <                    } else {
470 <                        pbt->setCoefficient(power, coefficient);
471 <                    }
610 >                    pbt->setCoefficient(power, coefficient);
611                  }
612              }
613              
# Line 480 | Line 619 | void DUFF::parseTorsionType(const std::string& line, i
619              } else {
620                  int nSets = nTokens / 3;
621    
483                double kchi;
484                int n;
485                double delta;                
622                  CharmmTorsionType* ctt = new CharmmTorsionType();
623                  
624                  for (int i = 0; i < nSets; ++i) {
625 <                    kchi = tokenizer.nextTokenAsDouble();
626 <                    n = tokenizer.nextTokenAsInt();
627 <                    delta = tokenizer.nextTokenAsDouble();
625 >                    double kchi = tokenizer.nextTokenAsDouble();
626 >                    int n = tokenizer.nextTokenAsInt();
627 >                    double delta = tokenizer.nextTokenAsDouble();
628      
629 <                    if (n < 0 || n > 6) {
494 <
495 <                    } else {
496 <                        ctt->setCharmmTorsionParameter(kchi, n, delta);
497 <                    }
629 >                    ctt->setCharmmTorsionParameter(kchi, n, delta);
630                  }
631              }
632          default:
# Line 507 | Line 639 | void DUFF::parseTorsionType(const std::string& line, i
639   }
640  
641   } //end namespace oopse
510

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines