--- trunk/src/applications/dump2Xyz/Dump2XYZ.cpp 2006/05/23 18:46:11 968 +++ trunk/src/applications/dump2Xyz/Dump2XYZ.cpp 2010/07/09 19:29:05 1464 @@ -6,19 +6,10 @@ * redistribute this software in source and binary code form, provided * that the following conditions are met: * - * 1. Acknowledgement of the program authors must be made in any - * publication of scientific results based in part on use of the - * program. An acceptable form of acknowledgement is citation of - * the article in which the program was described (Matthew - * A. Meineke, Charles F. Vardeman II, Teng Lin, Christopher - * J. Fennell and J. Daniel Gezelter, "OOPSE: An Object-Oriented - * Parallel Simulation Engine for Molecular Dynamics," - * J. Comput. Chem. 26, pp. 252-271 (2005)) - * - * 2. Redistributions of source code must retain the above copyright + * 1. Redistributions of source code must retain the above copyright * notice, this list of conditions and the following disclaimer. * - * 3. Redistributions in binary form must reproduce the above copyright + * 2. Redistributions in binary form must reproduce the above copyright * notice, this list of conditions and the following disclaimer in the * documentation and/or other materials provided with the * distribution. @@ -37,6 +28,15 @@ * arising out of the use of or inability to use software, even if the * University of Notre Dame has been advised of the possibility of * such damages. + * + * SUPPORT OPEN SCIENCE! If you use OpenMD or its source code in your + * research, please cite the appropriate papers when you publish your + * work. Good starting points are: + * + * [1] Meineke, et al., J. Comp. Chem. 26, 252-271 (2005). + * [2] Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006). + * [3] Sun, Lin & Gezelter, J. Chem. Phys. 128, 24107 (2008). + * [4] Vardeman & Gezelter, in progress (2009). */ #include @@ -47,9 +47,11 @@ #include "brains/Register.hpp" #include "brains/SimCreator.hpp" #include "brains/SimInfo.hpp" +#include "brains/ForceManager.hpp" #include "io/DumpReader.hpp" #include "utils/simError.h" #include "visitors/AtomVisitor.hpp" +#include "visitors/ReplacementVisitor.hpp" #include "visitors/CompositeVisitor.hpp" #include "visitors/RigidBodyVisitor.hpp" #include "visitors/OtherVisitor.hpp" @@ -59,7 +61,7 @@ #include "visitors/LipidTransVisitor.hpp" #include "visitors/AtomNameVisitor.hpp" -using namespace oopse; +using namespace OpenMD; int main(int argc, char* argv[]){ @@ -68,8 +70,11 @@ int main(int argc, char* argv[]){ gengetopt_args_info args_info; std::string dumpFileName; - std::string mdFileName; std::string xyzFileName; + bool printVel; + bool printFrc; + bool printVec; + bool printChrg; //parse the command line option if (cmdline_parser (argc, argv, &args_info) != 0) { @@ -84,9 +89,6 @@ int main(int argc, char* argv[]){ exit(1); } - mdFileName = dumpFileName; - mdFileName = mdFileName.substr(0, mdFileName.rfind(".")) + ".md"; - if (args_info.output_given){ xyzFileName = args_info.output_arg; } else { @@ -96,13 +98,12 @@ int main(int argc, char* argv[]){ //parse md file and set up the system SimCreator creator; - SimInfo* info = creator.createSim(mdFileName, false); + SimInfo* info = creator.createSim(dumpFileName, false); + ForceManager* forceMan = new ForceManager(info); - - //create visitor list CompositeVisitor* compositeVisitor = new CompositeVisitor(); - + //create RigidBody Visitor if(args_info.rigidbody_flag){ RBCOMVisitor* rbCOMVisitor = new RBCOMVisitor(info); @@ -113,29 +114,23 @@ int main(int argc, char* argv[]){ SSDAtomVisitor* ssdVisitor = new SSDAtomVisitor(info); compositeVisitor->addVisitor(ssdVisitor, 800); - LinearAtomVisitor* linearVisitor = new LinearAtomVisitor(info); - compositeVisitor->addVisitor(linearVisitor, 750); - if (args_info.gb_given) { - linearVisitor->addGayBerneAtomType(args_info.gb_arg); - } - - GBLipidAtomVisitor* gbLipidVisitor = new GBLipidAtomVisitor(info); - compositeVisitor->addVisitor(gbLipidVisitor, 740); - //create default atom visitor DefaultAtomVisitor* defaultAtomVisitor = new DefaultAtomVisitor(info); compositeVisitor->addVisitor(defaultAtomVisitor, 700); - //create waterType visitor - if(args_info.watertype_flag){ - WaterTypeVisitor* waterTypeVisitor = new WaterTypeVisitor; - compositeVisitor->addVisitor(waterTypeVisitor, 600); - } - + // if we gave the -w option, we want to skip the waters: + if (!args_info.water_given) { + //create waterType visitor + if(args_info.watertype_flag){ + WaterTypeVisitor* waterTypeVisitor = new WaterTypeVisitor; + compositeVisitor->addVisitor(waterTypeVisitor, 600); + } + } + if (args_info.basetype_flag) { - AtomNameVisitor* atomNameVisitor = new AtomNameVisitor(info); - compositeVisitor->addVisitor(atomNameVisitor, 550); - + AtomNameVisitor* atomNameVisitor = new AtomNameVisitor(info); + compositeVisitor->addVisitor(atomNameVisitor, 550); + std::cout << compositeVisitor->toString(); } //create ZconsVisitor @@ -158,31 +153,57 @@ int main(int argc, char* argv[]){ //} //create replicate visitor - if(args_info.repeatX_given > 0 || args_info.repeatY_given > 0 ||args_info.repeatY_given > 0){ - Vector3i replicateOpt(args_info.repeatX_arg, args_info.repeatY_arg, args_info.repeatZ_arg); - ReplicateVisitor* replicateVisitor = new ReplicateVisitor(info, replicateOpt); + if(args_info.repeatX_given > 0 || + args_info.repeatY_given > 0 || + args_info.repeatY_given > 0) { + Vector3i replicateOpt(args_info.repeatX_arg, + args_info.repeatY_arg, + args_info.repeatZ_arg); + ReplicateVisitor* replicateVisitor = new ReplicateVisitor(info, + replicateOpt); compositeVisitor->addVisitor(replicateVisitor, 300); } //create rotation visitor if (args_info.refsele_given&& args_info.originsele_given) { - compositeVisitor->addVisitor(new LipidTransVisitor(info, args_info.originsele_arg, args_info.refsele_arg), 250); + compositeVisitor->addVisitor(new LipidTransVisitor(info, + args_info.originsele_arg, + args_info.refsele_arg), + 250); } else if (args_info.refsele_given || args_info.originsele_given) { - std::cerr << "Both of --refsele and --originsele should appear by pair" << std::endl; + std::cerr << "Both of --refsele and --originsele should appear by pair" + << std::endl; exit(1); } - + //create xyzVisitor XYZVisitor* xyzVisitor; + if (args_info.selection_given) { xyzVisitor = new XYZVisitor(info, args_info.selection_arg); } else { xyzVisitor = new XYZVisitor(info); } - compositeVisitor->addVisitor(xyzVisitor, 200); + + if(args_info.velocities_flag){ + printVel = true; + xyzVisitor->doVelocities(printVel); + } + if(args_info.forces_flag){ + printFrc = true; + xyzVisitor->doForces(printFrc); + } + if(args_info.vectors_flag){ + printVec = true; + xyzVisitor->doVectors(printVec); + } + if(args_info.charges_flag){ + printChrg = true; + xyzVisitor->doCharges(printChrg); + } - std::cout << compositeVisitor->toString(); + compositeVisitor->addVisitor(xyzVisitor, 200); //create prepareVisitor PrepareVisitor* prepareVisitor = new PrepareVisitor(); @@ -191,11 +212,8 @@ int main(int argc, char* argv[]){ DumpReader* dumpReader = new DumpReader(info, dumpFileName); int nframes = dumpReader->getNFrames(); + std::ofstream xyzStream(xyzFileName.c_str()); - std::ofstream xyzStream; - xyzStream .open(xyzFileName.c_str()); - - SimInfo::MoleculeIterator miter; Molecule::IntegrableObjectIterator iiter; Molecule::RigidBodyIterator rbIter; @@ -210,33 +228,42 @@ int main(int argc, char* argv[]){ for (int i = 0; i < nframes; i += args_info.frame_arg){ dumpReader->readFrame(i); - + + if (printFrc) forceMan->calcForces(); + //wrapping the molecule if(args_info.periodicBox_flag) { currentSnapshot = info->getSnapshotManager()->getCurrentSnapshot(); - for (mol = info->beginMolecule(miter); mol != NULL; mol = info->nextMolecule(miter)) { + for (mol = info->beginMolecule(miter); mol != NULL; + mol = info->nextMolecule(miter)) { molCom = mol->getCom(); newMolCom = molCom; currentSnapshot->wrapVector(newMolCom); displacement = newMolCom - molCom; - for (integrableObject = mol->beginIntegrableObject(iiter); integrableObject != NULL; + for (integrableObject = mol->beginIntegrableObject(iiter); + integrableObject != NULL; integrableObject = mol->nextIntegrableObject(iiter)) { integrableObject->setPos(integrableObject->getPos() + displacement); } } } //update atoms of rigidbody - for (mol = info->beginMolecule(miter); mol != NULL; mol = info->nextMolecule(miter)) { + for (mol = info->beginMolecule(miter); mol != NULL; + mol = info->nextMolecule(miter)) { //change the positions of atoms which belong to the rigidbodies - for (rb = mol->beginRigidBody(rbIter); rb != NULL; rb = mol->nextRigidBody(rbIter)) { + for (rb = mol->beginRigidBody(rbIter); rb != NULL; + rb = mol->nextRigidBody(rbIter)) { rb->updateAtoms(); + if (printVel) rb->updateAtomVel(); } } //prepare visit - for (mol = info->beginMolecule(miter); mol != NULL; mol = info->nextMolecule(miter)) { - for (integrableObject = mol->beginIntegrableObject(iiter); integrableObject != NULL; + for (mol = info->beginMolecule(miter); mol != NULL; + mol = info->nextMolecule(miter)) { + for (integrableObject = mol->beginIntegrableObject(iiter); + integrableObject != NULL; integrableObject = mol->nextIntegrableObject(iiter)) { integrableObject->accept(prepareVisitor); } @@ -247,8 +274,10 @@ int main(int argc, char* argv[]){ //visit stuntdouble - for (mol = info->beginMolecule(miter); mol != NULL; mol = info->nextMolecule(miter)) { - for (integrableObject = mol->beginIntegrableObject(iiter); integrableObject != NULL; + for (mol = info->beginMolecule(miter); mol != NULL; + mol = info->nextMolecule(miter)) { + for (integrableObject = mol->beginIntegrableObject(iiter); + integrableObject != NULL; integrableObject = mol->nextIntegrableObject(iiter)) { integrableObject->accept(compositeVisitor); } @@ -258,11 +287,9 @@ int main(int argc, char* argv[]){ xyzVisitor->clear(); }//end for (int i = 0; i < nframes; i += args_info.frame_arg) - + xyzStream.close(); - delete prepareVisitor; delete compositeVisitor; delete info; - }