ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE-4/src/applications/simpleBuilder/simpleBuilder.cpp
(Generate patch)

Comparing trunk/OOPSE-4/src/applications/simpleBuilder/simpleBuilder.cpp (file contents):
Revision 2185 by tim, Tue Apr 12 22:42:13 2005 UTC vs.
Revision 2204 by gezelter, Fri Apr 15 22:04:00 2005 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 67 | Line 67 | int main(int argc, char *argv []) {
67  
68   int main(int argc, char *argv []) {
69  
70 <    //register force fields
71 <    registerForceFields();
72 <    registerLattice();
70 >  //register force fields
71 >  registerForceFields();
72 >  registerLattice();
73      
74 <    gengetopt_args_info args_info;
75 <    std::string latticeType;
76 <    std::string inputFileName;
77 <    std::string outPrefix;
78 <    std::string outMdFileName;
79 <    std::string outInitFileName;
80 <    Lattice *simpleLat;
81 <    int numMol;
82 <    double latticeConstant;
83 <    std::vector<double> lc;
84 <    double mass;
85 <    const double rhoConvertConst = 1.661;
86 <    double density;
87 <    int nx,
74 >  gengetopt_args_info args_info;
75 >  std::string latticeType;
76 >  std::string inputFileName;
77 >  std::string outPrefix;
78 >  std::string outMdFileName;
79 >  std::string outInitFileName;
80 >  Lattice *simpleLat;
81 >  int numMol;
82 >  double latticeConstant;
83 >  std::vector<double> lc;
84 >  double mass;
85 >  const double rhoConvertConst = 1.661;
86 >  double density;
87 >  int nx,
88      ny,
89      nz;
90 <    Mat3x3d hmat;
91 <    MoLocator *locator;
92 <    std::vector<Vector3d> latticePos;
93 <    std::vector<Vector3d> latticeOrt;
94 <    int numMolPerCell;
95 <    int curMolIndex;
96 <    DumpWriter *writer;
90 >  Mat3x3d hmat;
91 >  MoLocator *locator;
92 >  std::vector<Vector3d> latticePos;
93 >  std::vector<Vector3d> latticeOrt;
94 >  int numMolPerCell;
95 >  int curMolIndex;
96 >  DumpWriter *writer;
97  
98 <    // parse command line arguments
99 <    if (cmdline_parser(argc, argv, &args_info) != 0)
100 <        exit(1);
98 >  // parse command line arguments
99 >  if (cmdline_parser(argc, argv, &args_info) != 0)
100 >    exit(1);
101  
102 <    density = args_info.density_arg;
102 >  density = args_info.density_arg;
103  
104 <    //get lattice type
105 <    latticeType = UpperCase(args_info.latticetype_arg);
104 >  //get lattice type
105 >  latticeType = UpperCase(args_info.latticetype_arg);
106  
107 <    simpleLat = LatticeFactory::getInstance()->createLattice(latticeType);
107 >  simpleLat = LatticeFactory::getInstance()->createLattice(latticeType);
108      
109 <    if (simpleLat == NULL) {
110 <      sprintf(painCave.errMsg, "Lattice Factory can not create %s lattice\n",
111 <              latticeType.c_str());
112 <      painCave.isFatal = 1;
113 <      simError();
114 <    }
109 >  if (simpleLat == NULL) {
110 >    sprintf(painCave.errMsg, "Lattice Factory can not create %s lattice\n",
111 >            latticeType.c_str());
112 >    painCave.isFatal = 1;
113 >    simError();
114 >  }
115  
116 <    //get the number of unit cell
117 <    nx = args_info.nx_arg;
116 >  //get the number of unit cell
117 >  nx = args_info.nx_arg;
118  
119 <    if (nx <= 0) {
120 <        std::cerr << "The number of unit cell in h direction must be greater than 0" << std::endl;
121 <        exit(1);
122 <    }
119 >  if (nx <= 0) {
120 >    std::cerr << "The number of unit cell in h direction must be greater than 0" << std::endl;
121 >    exit(1);
122 >  }
123  
124 <    ny = args_info.ny_arg;
124 >  ny = args_info.ny_arg;
125  
126 <    if (ny <= 0) {
127 <        std::cerr << "The number of unit cell in l direction must be greater than 0" << std::endl;
128 <        exit(1);
129 <    }
126 >  if (ny <= 0) {
127 >    std::cerr << "The number of unit cell in l direction must be greater than 0" << std::endl;
128 >    exit(1);
129 >  }
130  
131 <    nz = args_info.nz_arg;
131 >  nz = args_info.nz_arg;
132  
133 <    if (nz <= 0) {
134 <        std::cerr << "The number of unit cell in k direction must be greater than 0" << std::endl;
135 <        exit(1);
136 <    }
133 >  if (nz <= 0) {
134 >    std::cerr << "The number of unit cell in k direction must be greater than 0" << std::endl;
135 >    exit(1);
136 >  }
137  
138 <    //get input file name
139 <    if (args_info.inputs_num)
140 <        inputFileName = args_info.inputs[0];
141 <    else {
142 <        std::cerr << "You must specify a input file name.\n" << std::endl;
143 <        cmdline_parser_print_help();
144 <        exit(1);
145 <    }
138 >  //get input file name
139 >  if (args_info.inputs_num)
140 >    inputFileName = args_info.inputs[0];
141 >  else {
142 >    std::cerr << "You must specify a input file name.\n" << std::endl;
143 >    cmdline_parser_print_help();
144 >    exit(1);
145 >  }
146  
147 <    //parse md file and set up the system
148 <    SimCreator oldCreator;
149 <    SimInfo* oldInfo = oldCreator.createSim(inputFileName, false);
147 >  //parse md file and set up the system
148 >  SimCreator oldCreator;
149 >  SimInfo* oldInfo = oldCreator.createSim(inputFileName, false);
150  
151 <    if (oldInfo->getNMoleculeStamp()>= 2) {
152 <        std::cerr << "can not build the system with more than two components"
153 <            << std::endl;
154 <        exit(1);
155 <    }
151 >  if (oldInfo->getNMoleculeStamp()>= 2) {
152 >    std::cerr << "can not build the system with more than two components"
153 >              << std::endl;
154 >    exit(1);
155 >  }
156  
157 <    //get mass of molecule.
157 >  //get mass of molecule.
158  
159 <    mass = getMolMass(oldInfo->getMoleculeStamp(0), oldInfo->getForceField());
159 >  mass = getMolMass(oldInfo->getMoleculeStamp(0), oldInfo->getForceField());
160  
161 <    //creat lattice
162 <    simpleLat = LatticeFactory::getInstance()->createLattice(latticeType);
161 >  //creat lattice
162 >  simpleLat = LatticeFactory::getInstance()->createLattice(latticeType);
163  
164 <    if (simpleLat == NULL) {
165 <        std::cerr << "Error in creating lattice" << std::endl;
166 <        exit(1);
167 <    }
164 >  if (simpleLat == NULL) {
165 >    std::cerr << "Error in creating lattice" << std::endl;
166 >    exit(1);
167 >  }
168  
169 <    numMolPerCell = simpleLat->getNumSitesPerCell();
169 >  numMolPerCell = simpleLat->getNumSitesPerCell();
170  
171 <    //calculate lattice constant (in Angstrom)
172 <    latticeConstant = pow(rhoConvertConst * numMolPerCell * mass / density,
173 <                          1.0 / 3.0);
171 >  //calculate lattice constant (in Angstrom)
172 >  latticeConstant = pow(rhoConvertConst * numMolPerCell * mass / density,
173 >                        1.0 / 3.0);
174  
175 <    //set lattice constant
176 <    lc.push_back(latticeConstant);
177 <    simpleLat->setLatticeConstant(lc);
175 >  //set lattice constant
176 >  lc.push_back(latticeConstant);
177 >  simpleLat->setLatticeConstant(lc);
178  
179 <    //calculate the total number of molecules
180 <    numMol = nx * ny * nz * numMolPerCell;
179 >  //calculate the total number of molecules
180 >  numMol = nx * ny * nz * numMolPerCell;
181  
182 <    if (oldInfo->getNGlobalMolecules() != numMol) {
183 <        outPrefix = getPrefix(inputFileName.c_str()) + "_" + latticeType;
184 <        outMdFileName = outPrefix + ".md";
182 >  if (oldInfo->getNGlobalMolecules() != numMol) {
183 >    outPrefix = getPrefix(inputFileName.c_str()) + "_" + latticeType;
184 >    outMdFileName = outPrefix + ".md";
185  
186 <        //creat new .md file on fly which corrects the number of molecule    
187 <        createMdFile(inputFileName, outMdFileName, numMol);
188 <        std::cerr
189 <            << "SimpleBuilder Error: the number of molecule and the density are not matched"
190 <            << std::endl;
191 <        std::cerr << "A new .md file: " << outMdFileName
192 <            << " is generated, use it to rerun the simpleBuilder" << std::endl;
193 <        exit(1);
194 <    }
195 <
196 <    //determine the output file names  
197 <    if (args_info.output_given)
198 <        outInitFileName = args_info.output_arg;
199 <    else
200 <        outInitFileName = getPrefix(inputFileName.c_str()) + ".in";
186 >    //creat new .md file on fly which corrects the number of molecule    
187 >    createMdFile(inputFileName, outMdFileName, numMol);
188 >    std::cerr
189 >      << "SimpleBuilder Error: the number of molecule and the density are not matched"
190 >      << std::endl;
191 >    std::cerr << "A new .md file: " << outMdFileName
192 >              << " is generated, use it to rerun the simpleBuilder" << std::endl;
193 >    exit(1);
194 >  }
195 >
196 >  //determine the output file names  
197 >  if (args_info.output_given)
198 >    outInitFileName = args_info.output_arg;
199 >  else
200 >    outInitFileName = getPrefix(inputFileName.c_str()) + ".in";
201      
202 <    //creat Molocator
203 <    locator = new MoLocator(oldInfo->getMoleculeStamp(0), oldInfo->getForceField());
202 >  //creat Molocator
203 >  locator = new MoLocator(oldInfo->getMoleculeStamp(0), oldInfo->getForceField());
204  
205 <    //fill Hmat
206 <    hmat(0, 0)= nx * latticeConstant;
207 <    hmat(0, 1) = 0.0;
208 <    hmat(0, 2) = 0.0;
205 >  //fill Hmat
206 >  hmat(0, 0)= nx * latticeConstant;
207 >  hmat(0, 1) = 0.0;
208 >  hmat(0, 2) = 0.0;
209  
210 <    hmat(1, 0) = 0.0;
211 <    hmat(1, 1) = ny * latticeConstant;
212 <    hmat(1, 2) = 0.0;
210 >  hmat(1, 0) = 0.0;
211 >  hmat(1, 1) = ny * latticeConstant;
212 >  hmat(1, 2) = 0.0;
213  
214 <    hmat(2, 0) = 0.0;
215 <    hmat(2, 1) = 0.0;
216 <    hmat(2, 2) = nz * latticeConstant;
214 >  hmat(2, 0) = 0.0;
215 >  hmat(2, 1) = 0.0;
216 >  hmat(2, 2) = nz * latticeConstant;
217  
218 <    //set Hmat
219 <    oldInfo->getSnapshotManager()->getCurrentSnapshot()->setHmat(hmat);
218 >  //set Hmat
219 >  oldInfo->getSnapshotManager()->getCurrentSnapshot()->setHmat(hmat);
220  
221 <    //place the molecules
221 >  //place the molecules
222  
223 <    curMolIndex = 0;
223 >  curMolIndex = 0;
224  
225 <    //get the orientation of the cell sites
226 <    //for the same type of molecule in same lattice, it will not change
227 <    latticeOrt = simpleLat->getLatticePointsOrt();
225 >  //get the orientation of the cell sites
226 >  //for the same type of molecule in same lattice, it will not change
227 >  latticeOrt = simpleLat->getLatticePointsOrt();
228  
229 <    Molecule* mol;
230 <    SimInfo::MoleculeIterator mi;
231 <    mol = oldInfo->beginMolecule(mi);
232 <    for(int i = 0; i < nx; i++) {
233 <        for(int j = 0; j < ny; j++) {
234 <            for(int k = 0; k < nz; k++) {
229 >  Molecule* mol;
230 >  SimInfo::MoleculeIterator mi;
231 >  mol = oldInfo->beginMolecule(mi);
232 >  for(int i = 0; i < nx; i++) {
233 >    for(int j = 0; j < ny; j++) {
234 >      for(int k = 0; k < nz; k++) {
235  
236 <                //get the position of the cell sites
237 <                simpleLat->getLatticePointsPos(latticePos, i, j, k);
236 >        //get the position of the cell sites
237 >        simpleLat->getLatticePointsPos(latticePos, i, j, k);
238  
239 <                for(int l = 0; l < numMolPerCell; l++) {
240 <                    if (mol != NULL) {
241 <                        locator->placeMol(latticePos[l], latticeOrt[l], mol);
242 <                    } else {
243 <                        std::cerr << std::endl;                    
244 <                    }
245 <                    mol = oldInfo->nextMolecule(mi);
246 <                }
247 <            }
248 <        }
239 >        for(int l = 0; l < numMolPerCell; l++) {
240 >          if (mol != NULL) {
241 >            locator->placeMol(latticePos[l], latticeOrt[l], mol);
242 >          } else {
243 >            std::cerr << std::endl;                    
244 >          }
245 >          mol = oldInfo->nextMolecule(mi);
246 >        }
247 >      }
248      }
249 +  }
250  
251 <    //create dumpwriter and write out the coordinates
252 <    oldInfo->setFinalConfigFileName(outInitFileName);
253 <    writer = new DumpWriter(oldInfo);
251 >  //create dumpwriter and write out the coordinates
252 >  oldInfo->setFinalConfigFileName(outInitFileName);
253 >  writer = new DumpWriter(oldInfo);
254  
255 <    if (writer == NULL) {
256 <        std::cerr << "error in creating DumpWriter" << std::endl;
257 <        exit(1);
258 <    }
255 >  if (writer == NULL) {
256 >    std::cerr << "error in creating DumpWriter" << std::endl;
257 >    exit(1);
258 >  }
259  
260 <    writer->writeDump();
261 <    std::cout << "new initial configuration file: " << outInitFileName
262 <        << " is generated." << std::endl;
260 >  writer->writeDump();
261 >  std::cout << "new initial configuration file: " << outInitFileName
262 >            << " is generated." << std::endl;
263  
264 <    //delete objects
264 >  //delete objects
265  
266 <    //delete oldInfo and oldSimSetup
267 <    if (oldInfo != NULL)
268 <        delete oldInfo;
266 >  //delete oldInfo and oldSimSetup
267 >  if (oldInfo != NULL)
268 >    delete oldInfo;
269  
270 <    if (writer != NULL)
271 <        delete writer;
270 >  if (writer != NULL)
271 >    delete writer;
272      
273 <    delete simpleLat;
273 >  delete simpleLat;
274  
275 <    return 0;
275 >  return 0;
276   }
277  
278   void createMdFile(const std::string&oldMdFileName, const std::string&newMdFileName,
279                    int numMol) {
280 <    ifstream oldMdFile;
281 <    ofstream newMdFile;
282 <    const int MAXLEN = 65535;
283 <    char buffer[MAXLEN];
280 >  ifstream oldMdFile;
281 >  ofstream newMdFile;
282 >  const int MAXLEN = 65535;
283 >  char buffer[MAXLEN];
284  
285 <    //create new .md file based on old .md file
286 <    oldMdFile.open(oldMdFileName.c_str());
287 <    newMdFile.open(newMdFileName.c_str());
285 >  //create new .md file based on old .md file
286 >  oldMdFile.open(oldMdFileName.c_str());
287 >  newMdFile.open(newMdFileName.c_str());
288  
289 <    oldMdFile.getline(buffer, MAXLEN);
289 >  oldMdFile.getline(buffer, MAXLEN);
290  
291 <    while (!oldMdFile.eof()) {
291 >  while (!oldMdFile.eof()) {
292  
293 <        //correct molecule number
294 <        if (strstr(buffer, "nMol") != NULL) {
295 <            sprintf(buffer, "\t\tnMol = %d;", numMol);
296 <            newMdFile << buffer << std::endl;
297 <        } else
298 <            newMdFile << buffer << std::endl;
293 >    //correct molecule number
294 >    if (strstr(buffer, "nMol") != NULL) {
295 >      sprintf(buffer, "\t\tnMol = %d;", numMol);
296 >      newMdFile << buffer << std::endl;
297 >    } else
298 >      newMdFile << buffer << std::endl;
299  
300 <        oldMdFile.getline(buffer, MAXLEN);
301 <    }
300 >    oldMdFile.getline(buffer, MAXLEN);
301 >  }
302  
303 <    oldMdFile.close();
304 <    newMdFile.close();
303 >  oldMdFile.close();
304 >  newMdFile.close();
305   }
306  

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines