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

Comparing trunk/src/applications/dump2Xyz/Dump2XYZ.cpp (file contents):
Revision 413 by tim, Wed Mar 9 17:30:29 2005 UTC vs.
Revision 1782 by gezelter, Wed Aug 22 02:28:28 2012 UTC

# Line 1 | Line 1
1 < /*
1 > /*
2   * Copyright (c) 2005 The University of Notre Dame. All Rights Reserved.
3   *
4   * The University of Notre Dame grants you ("Licensee") a
# 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]  Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010).
40 + * [5]  Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011).
41   */
42  
43   #include <iostream>
# Line 47 | Line 48
48   #include "brains/Register.hpp"
49   #include "brains/SimCreator.hpp"
50   #include "brains/SimInfo.hpp"
51 + #include "brains/ForceManager.hpp"
52   #include "io/DumpReader.hpp"
53   #include "utils/simError.h"
54   #include "visitors/AtomVisitor.hpp"
55 + #include "visitors/ReplacementVisitor.hpp"
56   #include "visitors/CompositeVisitor.hpp"
57   #include "visitors/RigidBodyVisitor.hpp"
58   #include "visitors/OtherVisitor.hpp"
# Line 57 | Line 60
60   #include "selection/SelectionEvaluator.hpp"
61   #include "selection/SelectionManager.hpp"
62   #include "visitors/LipidTransVisitor.hpp"
63 + #include "visitors/AtomNameVisitor.hpp"
64  
65 < using namespace oopse;
65 > using namespace OpenMD;
66  
67 + using namespace std;
68   int main(int argc, char* argv[]){
69    
65  //register force fields
66  registerForceFields();
67  
70    gengetopt_args_info args_info;
71 <  std::string dumpFileName;
72 <  std::string mdFileName;
73 <  std::string xyzFileName;
71 >  string dumpFileName;
72 >  string xyzFileName;
73 >
74 >  bool printVel(false);
75 >  bool printFrc(false);
76 >  bool printVec(false);
77 >  bool printChrg(false);
78    
79    //parse the command line option
80    if (cmdline_parser (argc, argv, &args_info) != 0) {
# Line 79 | Line 85 | int main(int argc, char* argv[]){
85    if (args_info.input_given){
86      dumpFileName = args_info.input_arg;
87    } else {
88 <    std::cerr << "Does not have input file name" << std::endl;
88 >    cerr << "Does not have input file name" << endl;
89      exit(1);
90    }
91    
86  mdFileName = dumpFileName;
87  mdFileName = mdFileName.substr(0, mdFileName.rfind(".")) + ".md";
88
92    if (args_info.output_given){
93      xyzFileName = args_info.output_arg;
94    } else {
# Line 95 | Line 98 | int main(int argc, char* argv[]){
98    
99    //parse md file and set up the system
100    SimCreator creator;
101 <  SimInfo* info = creator.createSim(mdFileName, false);
101 >  SimInfo* info = creator.createSim(dumpFileName, false);
102 >  ForceManager* forceMan = new ForceManager(info);
103    
104 <  
101 <  
102 <  //creat visitor list
104 >  //create visitor list
105    CompositeVisitor* compositeVisitor = new CompositeVisitor();
106 <    
107 <  //creat RigidBody Visitor
106 >  
107 >  //create RigidBody Visitor
108    if(args_info.rigidbody_flag){
109      RBCOMVisitor* rbCOMVisitor = new RBCOMVisitor(info);
110      compositeVisitor->addVisitor(rbCOMVisitor, 900);
111    }
112    
113 <  //creat SSD atom visitor
113 >  //create SSD atom visitor
114    SSDAtomVisitor* ssdVisitor = new SSDAtomVisitor(info);
115    compositeVisitor->addVisitor(ssdVisitor, 800);
116 +
117 +  //create GBtail atom visitor
118 +  GBtailVisitor* gbtVisitor = new GBtailVisitor(info);
119 +  compositeVisitor->addVisitor(gbtVisitor, 790);
120 +
121 +  //create GBhead atom visitor
122 +  GBheadVisitor* gbhVisitor = new GBheadVisitor(info);
123 +  compositeVisitor->addVisitor(gbhVisitor, 789);
124    
125 <  LinearAtomVisitor* linearVisitor = new LinearAtomVisitor(info);
116 <  compositeVisitor->addVisitor(linearVisitor, 750);
117 <  
118 <  //creat default atom visitor
125 >  //create default atom visitor
126    DefaultAtomVisitor* defaultAtomVisitor = new DefaultAtomVisitor(info);
127    compositeVisitor->addVisitor(defaultAtomVisitor, 700);
128    
129 <  //creat waterType visitor
130 <  if(args_info.watertype_flag){
131 <    WaterTypeVisitor* waterTypeVisitor = new WaterTypeVisitor;
132 <    compositeVisitor->addVisitor(waterTypeVisitor, 600);
129 >  // if we gave the -w option, we want to skip the waters:
130 >  if (!args_info.water_given) {
131 >    //create waterType visitor
132 >    if(args_info.watertype_flag){
133 >      WaterTypeVisitor* waterTypeVisitor = new WaterTypeVisitor;
134 >      compositeVisitor->addVisitor(waterTypeVisitor, 600);
135 >    }
136 >  }
137 >  
138 >  if (args_info.basetype_flag) {
139 >    AtomNameVisitor* atomNameVisitor = new AtomNameVisitor(info);
140 >    compositeVisitor->addVisitor(atomNameVisitor, 550);    
141 >    cout << compositeVisitor->toString();
142    }
143    
144    //create ZconsVisitor
# Line 137 | Line 153 | int main(int argc, char* argv[]){
153      }
154    }
155    
156 <  //creat wrapping visitor
156 >  //create wrapping visitor
157    
158 <  if(args_info.periodicBox_flag){
159 <    WrappingVisitor* wrappingVisitor = new WrappingVisitor(info);
160 <    compositeVisitor->addVisitor(wrappingVisitor, 400);
161 <  }
158 >  //if(args_info.periodicBox_flag){
159 >  //  WrappingVisitor* wrappingVisitor = new WrappingVisitor(info);
160 >  //  compositeVisitor->addVisitor(wrappingVisitor, 400);
161 >  //}
162  
163 <  //creat replicate visitor
164 <  if(args_info.repeatX_given > 0 || args_info.repeatY_given > 0 ||args_info.repeatY_given > 0){
165 <    Vector3i replicateOpt(args_info.repeatX_arg, args_info.repeatY_arg, args_info.repeatZ_arg);
166 <    ReplicateVisitor* replicateVisitor = new ReplicateVisitor(info, replicateOpt);
163 >  //create replicate visitor
164 >  if(args_info.repeatX_given > 0 ||
165 >     args_info.repeatY_given > 0 ||
166 >     args_info.repeatY_given > 0) {
167 >    Vector3i replicateOpt(args_info.repeatX_arg,
168 >                          args_info.repeatY_arg,
169 >                          args_info.repeatZ_arg);
170 >    ReplicateVisitor* replicateVisitor = new ReplicateVisitor(info,
171 >                                                              replicateOpt);
172      compositeVisitor->addVisitor(replicateVisitor, 300);
173    }
174  
175  
176    //create rotation visitor
177    if (args_info.refsele_given&& args_info.originsele_given) {
178 <    compositeVisitor->addVisitor(new LipidTransVisitor(info, args_info.originsele_arg, args_info.refsele_arg), 250);
178 >    compositeVisitor->addVisitor(new LipidTransVisitor(info,
179 >                                                       args_info.originsele_arg,
180 >                                                       args_info.refsele_arg),
181 >                                 250);
182    } else if (args_info.refsele_given || args_info.originsele_given) {
183 <    std::cerr << "Both of --refsele and --originsele should appear by pair" << std::endl;
183 >    cerr << "Both of --refsele and --originsele should appear by pair"
184 >         << endl;
185      exit(1);
186    }
187 <    
188 <  //creat xyzVisitor
187 >  
188 >  //create xyzVisitor
189    XYZVisitor* xyzVisitor;
190 +
191    if (args_info.selection_given) {
192      xyzVisitor = new XYZVisitor(info, args_info.selection_arg);
193    } else {
194      xyzVisitor = new XYZVisitor(info);
195    }
196 <  compositeVisitor->addVisitor(xyzVisitor, 200);
196 >
197 >  if(args_info.velocities_flag){
198 >    printVel = true;
199 >    xyzVisitor->doVelocities(printVel);
200 >  }
201 >  if(args_info.forces_flag){
202 >    printFrc = true;
203 >    xyzVisitor->doForces(printFrc);
204 >  }
205 >  if(args_info.vectors_flag){
206 >    printVec = true;
207 >    xyzVisitor->doVectors(printVec);
208 >  }
209 >  if(args_info.charges_flag){
210 >    printChrg = true;
211 >    xyzVisitor->doCharges(printChrg);
212 >  }
213    
214 <  std::cout << compositeVisitor->toString();
214 >  compositeVisitor->addVisitor(xyzVisitor, 200);
215    
216 <  //creat prepareVisitor
216 >  //create prepareVisitor
217    PrepareVisitor* prepareVisitor = new PrepareVisitor();
218    
219    //open dump file
220    DumpReader* dumpReader = new DumpReader(info, dumpFileName);
221    int nframes = dumpReader->getNFrames();
222    
223 +  ofstream xyzStream(xyzFileName.c_str());
224    
182  std::ofstream xyzStream;    
183  xyzStream .open(xyzFileName.c_str());
184  
185  
225    SimInfo::MoleculeIterator miter;
226    Molecule::IntegrableObjectIterator  iiter;
227    Molecule::RigidBodyIterator rbIter;
228    Molecule* mol;
229 <  StuntDouble* integrableObject;
229 >  StuntDouble* sd;
230    RigidBody* rb;
231 +  Vector3d molCom;
232 +  Vector3d newMolCom;
233 +  Vector3d displacement;
234 +  Mat3x3d hmat;
235 +  Snapshot* currentSnapshot;
236        
237    for (int i = 0; i < nframes; i += args_info.frame_arg){
238      dumpReader->readFrame(i);
239      
240 +    if (printFrc) forceMan->calcForces();
241 +    
242 +    //wrapping the molecule
243 +    if(args_info.periodicBox_flag) {
244 +      currentSnapshot = info->getSnapshotManager()->getCurrentSnapshot();    
245 +      for (mol = info->beginMolecule(miter); mol != NULL;
246 +           mol = info->nextMolecule(miter)) {
247 +        
248 +        molCom = mol->getCom();
249 +        newMolCom = molCom;
250 +        currentSnapshot->wrapVector(newMolCom);
251 +        displacement = newMolCom - molCom;
252 +
253 +        for (sd = mol->beginIntegrableObject(iiter); sd != NULL;
254 +             sd = mol->nextIntegrableObject(iiter)) {  
255 +
256 +          sd->setPos(sd->getPos() + displacement);
257 +          
258 +        }
259 +      }    
260 +    }
261 +
262      //update atoms of rigidbody
263 <    for (mol = info->beginMolecule(miter); mol != NULL; mol = info->nextMolecule(miter)) {
263 >    for (mol = info->beginMolecule(miter); mol != NULL;
264 >         mol = info->nextMolecule(miter)) {
265        
266        //change the positions of atoms which belong to the rigidbodies
267 <      for (rb = mol->beginRigidBody(rbIter); rb != NULL; rb = mol->nextRigidBody(rbIter)) {
267 >      for (rb = mol->beginRigidBody(rbIter); rb != NULL;
268 >           rb = mol->nextRigidBody(rbIter)) {
269 >
270          rb->updateAtoms();
271 +        if (printVel) rb->updateAtomVel();
272 +
273        }
274      }
275      
276      //prepare visit
277 <    for (mol = info->beginMolecule(miter); mol != NULL; mol = info->nextMolecule(miter)) {
278 <      for (integrableObject = mol->beginIntegrableObject(iiter); integrableObject != NULL;
279 <           integrableObject = mol->nextIntegrableObject(iiter)) {
280 <        integrableObject->accept(prepareVisitor);
277 >    for (mol = info->beginMolecule(miter); mol != NULL;
278 >         mol = info->nextMolecule(miter)) {
279 >
280 >      for (sd = mol->beginIntegrableObject(iiter); sd != NULL;
281 >           sd = mol->nextIntegrableObject(iiter)) {
282 >
283 >        sd->accept(prepareVisitor);
284 >
285        }
286      }
287      
# Line 215 | Line 290 | int main(int argc, char* argv[]){
290  
291  
292      //visit stuntdouble
293 <    for (mol = info->beginMolecule(miter); mol != NULL; mol = info->nextMolecule(miter)) {
294 <      for (integrableObject = mol->beginIntegrableObject(iiter); integrableObject != NULL;
295 <           integrableObject = mol->nextIntegrableObject(iiter)) {
296 <        integrableObject->accept(compositeVisitor);
293 >    for (mol = info->beginMolecule(miter); mol != NULL;
294 >         mol = info->nextMolecule(miter)) {
295 >
296 >      for (sd = mol->beginIntegrableObject(iiter); sd != NULL;
297 >           sd = mol->nextIntegrableObject(iiter)) {
298 >
299 >        sd->accept(compositeVisitor);
300 >
301        }
302      }
303      
# Line 226 | Line 305 | int main(int argc, char* argv[]){
305      xyzVisitor->clear();
306      
307    }//end for (int i = 0; i < nframes; i += args_info.frame_arg)
308 <  
308 >
309    xyzStream.close();
310 <  
232 <  
310 >  delete prepareVisitor;
311    delete compositeVisitor;
312    delete info;
235  
313   }

Comparing trunk/src/applications/dump2Xyz/Dump2XYZ.cpp (property svn:keywords):
Revision 413 by tim, Wed Mar 9 17:30:29 2005 UTC vs.
Revision 1782 by gezelter, Wed Aug 22 02:28:28 2012 UTC

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

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines