# | Line 6 | Line 6 | |
---|---|---|
6 | #include "parse_me.h" | |
7 | #include "LRI.hpp" | |
8 | #include "Integrator.hpp" | |
9 | + | #include "simError.h" |
10 | ||
11 | + | #ifdef IS_MPI |
12 | + | #include "mpiBASS.h" |
13 | + | #include "mpiSimulation.hpp" |
14 | + | #include "bassDiag.hpp" |
15 | + | #endif |
16 | + | |
17 | SimSetup::SimSetup(){ | |
18 | stamps = new MakeStamps(); | |
19 | globals = new Globals(); | |
20 | + | |
21 | + | #ifdef IS_MPI |
22 | + | strcpy( checkPointMsg, "SimSetup creation successful" ); |
23 | + | MPIcheckPoint(); |
24 | + | #endif // IS_MPI |
25 | } | |
26 | ||
27 | SimSetup::~SimSetup(){ | |
# | Line 19 | Line 31 | void SimSetup::parseFile( char* fileName ){ | |
31 | ||
32 | void SimSetup::parseFile( char* fileName ){ | |
33 | ||
34 | < | inFileName = fileName; |
35 | < | set_interface_stamps( stamps, globals ); |
36 | < | yacc_BASS( fileName ); |
34 | > | #ifdef IS_MPI |
35 | > | if( worldRank == 0 ){ |
36 | > | #endif // is_mpi |
37 | > | |
38 | > | inFileName = fileName; |
39 | > | set_interface_stamps( stamps, globals ); |
40 | > | |
41 | > | #ifdef IS_MPI |
42 | > | mpiEventInit(); |
43 | > | #endif |
44 | > | |
45 | > | yacc_BASS( fileName ); |
46 | > | |
47 | > | #ifdef IS_MPI |
48 | > | throwMPIEvent(NULL); |
49 | > | } |
50 | > | else receiveParse(); |
51 | > | #endif |
52 | > | |
53 | } | |
54 | ||
55 | + | #ifdef IS_MPI |
56 | + | void SimSetup::receiveParse(void){ |
57 | + | |
58 | + | set_interface_stamps( stamps, globals ); |
59 | + | mpiEventInit(); |
60 | + | MPIcheckPoint(); |
61 | + | mpiEventLoop(); |
62 | + | |
63 | + | } |
64 | + | |
65 | + | |
66 | + | void SimSetup::testMe(void){ |
67 | + | bassDiag* dumpMe = new bassDiag(globals,stamps); |
68 | + | dumpMe->dumpStamps(); |
69 | + | delete dumpMe; |
70 | + | } |
71 | + | #endif |
72 | + | |
73 | void SimSetup::createSim( void ){ | |
74 | ||
75 | MakeStamps *the_stamps; | |
76 | Globals* the_globals; | |
77 | < | int i; |
77 | > | int i, j; |
78 | ||
79 | // get the stamps and globals; | |
80 | the_stamps = stamps; | |
# | Line 43 | Line 89 | void SimSetup::createSim( void ){ | |
89 | n_components = the_globals->getNComponents(); | |
90 | strcpy( force_field, the_globals->getForceField() ); | |
91 | strcpy( ensemble, the_globals->getEnsemble() ); | |
92 | < | |
92 | > | |
93 | if( !strcmp( force_field, "TraPPE" ) ) the_ff = new TraPPEFF(); | |
94 | else if( !strcmp( force_field, "DipoleTest" ) ) the_ff = new DipoleTestFF(); | |
95 | else if( !strcmp( force_field, "TraPPE_Ex" ) ) the_ff = new TraPPE_ExFF(); | |
96 | + | else if( !strcmp( force_field, "LJ" ) ) the_ff = new LJ_FF(); |
97 | else{ | |
98 | < | std::cerr<< "SimSetup Error. Unrecognized force field -> " |
99 | < | << force_field << "\n"; |
100 | < | exit(8); |
98 | > | sprintf( painCave.errMsg, |
99 | > | "SimSetup Error. Unrecognized force field -> %s\n", |
100 | > | force_field ); |
101 | > | painCave.isFatal = 1; |
102 | > | simError(); |
103 | } | |
104 | ||
105 | + | #ifdef IS_MPI |
106 | + | strcpy( checkPointMsg, "ForceField creation successful" ); |
107 | + | MPIcheckPoint(); |
108 | + | #endif // is_mpi |
109 | + | |
110 | // get the components and calculate the tot_nMol and indvidual n_mol | |
111 | the_components = the_globals->getComponents(); | |
112 | components_nmol = new int[n_components]; | |
113 | comp_stamps = new MoleculeStamp*[n_components]; | |
114 | < | |
114 | > | |
115 | if( !the_globals->haveNMol() ){ | |
116 | < | // we don't have the total number of molecules, so we assume it is |
116 | > | // we don't have the total number of molecules, so we assume it is |
117 | // given in each component | |
118 | ||
119 | tot_nmol = 0; | |
120 | for( i=0; i<n_components; i++ ){ | |
121 | < | |
121 | > | |
122 | if( !the_components[i]->haveNMol() ){ | |
123 | // we have a problem | |
124 | < | std::cerr << "SimSetup Error. No global NMol or component NMol" |
125 | < | << " given. Cannot calculate the number of atoms.\n"; |
126 | < | exit( 8 ); |
124 | > | sprintf( painCave.errMsg, |
125 | > | "SimSetup Error. No global NMol or component NMol" |
126 | > | " given. Cannot calculate the number of atoms.\n" ); |
127 | > | painCave.isFatal = 1; |
128 | > | simError(); |
129 | } | |
130 | ||
131 | tot_nmol += the_components[i]->getNMol(); | |
# | Line 77 | Line 133 | void SimSetup::createSim( void ){ | |
133 | } | |
134 | } | |
135 | else{ | |
136 | < | std::cerr << "NOT A SUPPORTED FEATURE\n"; |
136 | > | sprintf( painCave.errMsg, |
137 | > | "SimSetup error.\n" |
138 | > | "\tSorry, the ability to specify total" |
139 | > | " nMols and then give molfractions in the components\n" |
140 | > | "\tis not currently supported." |
141 | > | " Please give nMol in the components.\n" ); |
142 | > | painCave.isFatal = 1; |
143 | > | simError(); |
144 | ||
82 | – | // tot_nmol = the_globals->getNMol(); |
145 | ||
146 | < | // //we have the total number of molecules, now we check for molfractions |
147 | < | // for( i=0; i<n_components; i++ ){ |
148 | < | |
149 | < | // if( !the_components[i]->haveMolFraction() ){ |
150 | < | |
151 | < | // if( !the_components[i]->haveNMol() ){ |
152 | < | // //we have a problem |
153 | < | // std::cerr << "SimSetup error. Neither molFraction nor " |
154 | < | // << " nMol was given in component |
155 | < | |
146 | > | // tot_nmol = the_globals->getNMol(); |
147 | > | |
148 | > | // //we have the total number of molecules, now we check for molfractions |
149 | > | // for( i=0; i<n_components; i++ ){ |
150 | > | |
151 | > | // if( !the_components[i]->haveMolFraction() ){ |
152 | > | |
153 | > | // if( !the_components[i]->haveNMol() ){ |
154 | > | // //we have a problem |
155 | > | // std::cerr << "SimSetup error. Neither molFraction nor " |
156 | > | // << " nMol was given in component |
157 | > | |
158 | } | |
159 | ||
160 | + | #ifdef IS_MPI |
161 | + | strcpy( checkPointMsg, "Have the number of components" ); |
162 | + | MPIcheckPoint(); |
163 | + | #endif // is_mpi |
164 | + | |
165 | // make an array of molecule stamps that match the components used. | |
166 | + | // also extract the used stamps out into a separate linked list |
167 | ||
168 | + | simnfo->nComponents = n_components; |
169 | + | simnfo->componentsNmol = components_nmol; |
170 | + | simnfo->compStamps = comp_stamps; |
171 | + | simnfo->headStamp = new LinkedMolStamp(); |
172 | + | |
173 | + | char* id; |
174 | + | LinkedMolStamp* headStamp = simnfo->headStamp; |
175 | + | LinkedMolStamp* currentStamp = NULL; |
176 | for( i=0; i<n_components; i++ ){ | |
177 | ||
178 | < | comp_stamps[i] = |
179 | < | the_stamps->getMolecule( the_components[i]->getType() ); |
178 | > | id = the_components[i]->getType(); |
179 | > | comp_stamps[i] = NULL; |
180 | > | |
181 | > | // check to make sure the component isn't already in the list |
182 | > | |
183 | > | comp_stamps[i] = headStamp->match( id ); |
184 | > | if( comp_stamps[i] == NULL ){ |
185 | > | |
186 | > | // extract the component from the list; |
187 | > | |
188 | > | currentStamp = the_stamps->extractMolStamp( id ); |
189 | > | if( currentStamp == NULL ){ |
190 | > | sprintf( painCave.errMsg, |
191 | > | "SimSetup error: Component \"%s\" was not found in the " |
192 | > | "list of declared molecules\n" |
193 | > | id ); |
194 | > | painCave.isFatal = 1; |
195 | > | simError(); |
196 | > | } |
197 | > | |
198 | > | headStamp->add( currentStamp ); |
199 | > | comp_stamps[i] = headStamp->match( id ); |
200 | > | } |
201 | } | |
202 | ||
203 | + | #ifdef IS_MPI |
204 | + | strcpy( checkPointMsg, "Component stamps successfully extracted\n" ); |
205 | + | MPIcheckPoint(); |
206 | + | #endif // is_mpi |
207 | ||
208 | ||
209 | + | |
210 | + | |
211 | // caclulate the number of atoms, bonds, bends and torsions | |
212 | ||
213 | tot_atoms = 0; | |
# | Line 111 | Line 216 | void SimSetup::createSim( void ){ | |
216 | tot_torsions = 0; | |
217 | for( i=0; i<n_components; i++ ){ | |
218 | ||
219 | < | tot_atoms += components_nmol[i] * comp_stamps[i]->getNAtoms(); |
220 | < | tot_bonds += components_nmol[i] * comp_stamps[i]->getNBonds(); |
221 | < | tot_bends += components_nmol[i] * comp_stamps[i]->getNBends(); |
219 | > | tot_atoms += components_nmol[i] * comp_stamps[i]->getNAtoms(); |
220 | > | tot_bonds += components_nmol[i] * comp_stamps[i]->getNBonds(); |
221 | > | tot_bends += components_nmol[i] * comp_stamps[i]->getNBends(); |
222 | tot_torsions += components_nmol[i] * comp_stamps[i]->getNTorsions(); | |
223 | } | |
224 | < | |
224 | > | |
225 | tot_SRI = tot_bonds + tot_bends + tot_torsions; | |
226 | < | |
226 | > | |
227 | simnfo->n_atoms = tot_atoms; | |
228 | simnfo->n_bonds = tot_bonds; | |
229 | simnfo->n_bends = tot_bends; | |
230 | simnfo->n_torsions = tot_torsions; | |
231 | simnfo->n_SRI = tot_SRI; | |
232 | + | simnfo->n_mol = tot_nmol; |
233 | ||
128 | – | // create the atom and short range interaction arrays |
234 | ||
235 | < | the_atoms = new Atom*[tot_atoms]; |
236 | < | // the_molecules = new Molecule[tot_nmol]; |
235 | > | #ifdef IS_MPI |
236 | > | |
237 | > | // divide the molecules among processors here. |
238 | ||
239 | + | mpiSim = new mpiSimulation( simnfo ); |
240 | ||
241 | < | if( tot_SRI ){ |
242 | < | the_sris = new SRI*[tot_SRI]; |
243 | < | the_excludes = new ex_pair[tot_SRI]; |
241 | > | mpiSim->divideLabor(); |
242 | > | |
243 | > | // set up the local variables |
244 | > | |
245 | > | int localMol, allMol; |
246 | > | int local_atoms, local_bonds, local_bends, local_torsions, local_SRI; |
247 | > | |
248 | > | allMol = 0; |
249 | > | localMol = 0; |
250 | > | local_atoms = 0; |
251 | > | local_bonds = 0; |
252 | > | local_bends = 0; |
253 | > | local_torsions = 0; |
254 | > | for( i=0; i<n_components; i++ ){ |
255 | > | |
256 | > | for( j=0; j<components_nmol[i]; j++ ){ |
257 | > | |
258 | > | if( mpiSim->getMyMolStart() <= allMol && |
259 | > | allMol <= mpiSim->getMyMolEnd() ){ |
260 | > | |
261 | > | local_atoms += comp_stamps[i]->getNAtoms(); |
262 | > | local_bonds += comp_stamps[i]->getNBonds(); |
263 | > | local_bends += comp_stamps[i]->getNBends(); |
264 | > | local_torsions += comp_stamps[i]->getNTorsions(); |
265 | > | localMol++; |
266 | > | } |
267 | > | allMol++; |
268 | > | } |
269 | } | |
270 | + | local_SRI = local_bonds + local_bends + local_torsions; |
271 | + | |
272 | ||
273 | + | simnfo->n_atoms = mpiSim->getMyNlocal(); |
274 | + | |
275 | + | if( local_atoms != simnfo->n_atoms ){ |
276 | + | sprintf( painCave.errMsg, |
277 | + | "SimSetup error: mpiSim's localAtom (%d) and SimSetup's" |
278 | + | " localAtom (%d) are note equal.\n", |
279 | + | simnfo->n_atoms, |
280 | + | local_atoms ); |
281 | + | painCave.isFatal = 1; |
282 | + | simError(); |
283 | + | } |
284 | + | |
285 | + | simnfo->n_bonds = local_bonds; |
286 | + | simnfo->n_bends = local_bends; |
287 | + | simnfo->n_torsions = local_torsions; |
288 | + | simnfo->n_SRI = local_SRI; |
289 | + | simnfo->n_mol = localMol; |
290 | + | |
291 | + | strcpy( checkPointMsg, "Passed nlocal consistency check." ); |
292 | + | MPIcheckPoint(); |
293 | + | |
294 | + | |
295 | + | #endif // is_mpi |
296 | + | |
297 | + | |
298 | + | // create the atom and short range interaction arrays |
299 | + | |
300 | + | Atom::createArrays(simnfo->n_atoms); |
301 | + | the_atoms = new Atom*[simnfo->n_atoms]; |
302 | + | the_molecules = new Molecule[simnfo->n_mol]; |
303 | + | |
304 | + | |
305 | + | if( simnfo->n_SRI ){ |
306 | + | the_sris = new SRI*[simnfo->n_SRI]; |
307 | + | the_excludes = new ex_pair[simnfo->n_SRI]; |
308 | + | } |
309 | + | |
310 | // set the arrays into the SimInfo object | |
311 | ||
312 | simnfo->atoms = the_atoms; | |
# | Line 143 | Line 314 | void SimSetup::createSim( void ){ | |
314 | simnfo->n_exclude = tot_SRI; | |
315 | simnfo->excludes = the_excludes; | |
316 | ||
317 | + | |
318 | // initialize the arrays | |
319 | < | |
319 | > | |
320 | the_ff->setSimInfo( simnfo ); | |
321 | < | |
321 | > | |
322 | makeAtoms(); | |
323 | ||
324 | if( tot_bonds ){ | |
# | Line 161 | Line 333 | void SimSetup::createSim( void ){ | |
333 | makeTorsions(); | |
334 | } | |
335 | ||
164 | – | // makeMolecules(); |
336 | ||
337 | // get some of the tricky things that may still be in the globals | |
338 | ||
339 | if( simnfo->n_dipoles ){ | |
340 | ||
341 | if( !the_globals->haveRRF() ){ | |
342 | < | std::cerr << "SimSetup Error, system has dipoles, but no rRF was set.\n"; |
343 | < | exit(8); |
342 | > | sprintf( painCave.errMsg, |
343 | > | "SimSetup Error, system has dipoles, but no rRF was set.\n"); |
344 | > | painCave.isFatal = 1; |
345 | > | simError(); |
346 | } | |
347 | if( !the_globals->haveDielectric() ){ | |
348 | < | std::cerr << "SimSetup Error, system has dipoles, but no" |
349 | < | << " dielectric was set.\n"; |
350 | < | exit(8); |
348 | > | sprintf( painCave.errMsg, |
349 | > | "SimSetup Error, system has dipoles, but no" |
350 | > | " dielectric was set.\n" ); |
351 | > | painCave.isFatal = 1; |
352 | > | simError(); |
353 | } | |
354 | ||
355 | simnfo->rRF = the_globals->getRRF(); | |
356 | simnfo->dielectric = the_globals->getDielectric(); | |
357 | } | |
358 | ||
359 | + | #ifdef IS_MPI |
360 | + | strcpy( checkPointMsg, "rRf and dielectric check out" ); |
361 | + | MPIcheckPoint(); |
362 | + | #endif // is_mpi |
363 | + | |
364 | if( the_globals->haveBox() ){ | |
365 | simnfo->box_x = the_globals->getBox(); | |
366 | simnfo->box_y = the_globals->getBox(); | |
367 | simnfo->box_z = the_globals->getBox(); | |
368 | } | |
369 | else if( the_globals->haveDensity() ){ | |
370 | < | |
370 | > | |
371 | double vol; | |
372 | vol = (double)tot_nmol / the_globals->getDensity(); | |
373 | simnfo->box_x = pow( vol, ( 1.0 / 3.0 ) ); | |
# | Line 196 | Line 376 | void SimSetup::createSim( void ){ | |
376 | } | |
377 | else{ | |
378 | if( !the_globals->haveBoxX() ){ | |
379 | < | std::cerr << "SimSetup error, no periodic BoxX size given.\n"; |
380 | < | exit(8); |
379 | > | sprintf( painCave.errMsg, |
380 | > | "SimSetup error, no periodic BoxX size given.\n" ); |
381 | > | painCave.isFatal = 1; |
382 | > | simError(); |
383 | } | |
384 | simnfo->box_x = the_globals->getBoxX(); | |
385 | ||
386 | if( !the_globals->haveBoxY() ){ | |
387 | < | std::cerr << "SimSetup error, no periodic BoxY size given.\n"; |
388 | < | exit(8); |
387 | > | sprintf( painCave.errMsg, |
388 | > | "SimSetup error, no periodic BoxY size given.\n" ); |
389 | > | painCave.isFatal = 1; |
390 | > | simError(); |
391 | } | |
392 | simnfo->box_y = the_globals->getBoxY(); | |
393 | ||
394 | if( !the_globals->haveBoxZ() ){ | |
395 | < | std::cerr << "SimSetup error, no periodic BoxZ size given.\n"; |
396 | < | exit(8); |
395 | > | sprintf( painCave.errMsg, |
396 | > | "SimSetup error, no periodic BoxZ size given.\n" ); |
397 | > | painCave.isFatal = 1; |
398 | > | simError(); |
399 | } | |
400 | simnfo->box_z = the_globals->getBoxZ(); | |
401 | } | |
216 | – | |
217 | – | if( the_globals->haveInitialConfig() ){ |
218 | – | InitializeFromFile* fileInit; |
219 | – | fileInit = new InitializeFromFile( the_globals->getInitialConfig() ); |
220 | – | |
221 | – | fileInit->read_xyz( simnfo ); // default velocities on |
402 | ||
403 | < | delete fileInit; |
404 | < | } |
405 | < | else{ |
406 | < | initFromBass(); |
227 | < | } |
403 | > | #ifdef IS_MPI |
404 | > | strcpy( checkPointMsg, "Box size set up" ); |
405 | > | MPIcheckPoint(); |
406 | > | #endif // is_mpi |
407 | ||
408 | < | if( the_globals->haveFinalConfig() ){ |
409 | < | strcpy( simnfo->finalName, the_globals->getFinalConfig() ); |
410 | < | } |
411 | < | else{ |
412 | < | strcpy( simnfo->finalName, inFileName ); |
413 | < | char* endTest; |
414 | < | int nameLength = strlen( simnfo->finalName ); |
415 | < | endTest = &(simnfo->finalName[nameLength - 5]); |
416 | < | if( !strcmp( endTest, ".bass" ) ){ |
417 | < | strcpy( endTest, ".eor" ); |
408 | > | |
409 | > | |
410 | > | if( the_globals->haveInitialConfig() ){ |
411 | > | |
412 | > | InitializeFromFile* fileInit; |
413 | > | #ifdef IS_MPI // is_mpi |
414 | > | if( worldRank == 0 ){ |
415 | > | #endif //is_mpi |
416 | > | fileInit = new InitializeFromFile( the_globals->getInitialConfig() ); |
417 | > | #ifdef IS_MPI |
418 | > | }else fileInit = new InitializeFromFile( NULL ); |
419 | > | #endif |
420 | > | fileInit->read_xyz( simnfo ); // default velocities on |
421 | > | |
422 | > | delete fileInit; |
423 | > | } |
424 | > | else{ |
425 | > | |
426 | > | #ifdef IS_MPI |
427 | > | |
428 | > | // no init from bass |
429 | > | |
430 | > | sprintf( painCave.errMsg, |
431 | > | "Cannot intialize a parallel simulation without an initial configuration file.\n" ); |
432 | > | painCave.isFatal; |
433 | > | simError(); |
434 | > | |
435 | > | #else |
436 | > | |
437 | > | initFromBass(); |
438 | > | |
439 | > | |
440 | > | #endif |
441 | > | } |
442 | > | |
443 | > | #ifdef IS_MPI |
444 | > | strcpy( checkPointMsg, "Successfully read in the initial configuration" ); |
445 | > | MPIcheckPoint(); |
446 | > | #endif // is_mpi |
447 | > | |
448 | > | |
449 | > | |
450 | > | |
451 | > | |
452 | > | |
453 | > | |
454 | > | #ifdef IS_MPI |
455 | > | if( worldRank == 0 ){ |
456 | > | #endif // is_mpi |
457 | > | |
458 | > | if( the_globals->haveFinalConfig() ){ |
459 | > | strcpy( simnfo->finalName, the_globals->getFinalConfig() ); |
460 | } | |
240 | – | else if( !strcmp( endTest, ".BASS" ) ){ |
241 | – | strcpy( endTest, ".eor" ); |
242 | – | } |
461 | else{ | |
462 | < | endTest = &(simnfo->finalName[nameLength - 4]); |
463 | < | if( !strcmp( endTest, ".bss" ) ){ |
462 | > | strcpy( simnfo->finalName, inFileName ); |
463 | > | char* endTest; |
464 | > | int nameLength = strlen( simnfo->finalName ); |
465 | > | endTest = &(simnfo->finalName[nameLength - 5]); |
466 | > | if( !strcmp( endTest, ".bass" ) ){ |
467 | strcpy( endTest, ".eor" ); | |
468 | } | |
469 | < | else if( !strcmp( endTest, ".mdl" ) ){ |
469 | > | else if( !strcmp( endTest, ".BASS" ) ){ |
470 | strcpy( endTest, ".eor" ); | |
471 | } | |
472 | else{ | |
473 | < | strcat( simnfo->finalName, ".eor" ); |
473 | > | endTest = &(simnfo->finalName[nameLength - 4]); |
474 | > | if( !strcmp( endTest, ".bss" ) ){ |
475 | > | strcpy( endTest, ".eor" ); |
476 | > | } |
477 | > | else if( !strcmp( endTest, ".mdl" ) ){ |
478 | > | strcpy( endTest, ".eor" ); |
479 | > | } |
480 | > | else{ |
481 | > | strcat( simnfo->finalName, ".eor" ); |
482 | > | } |
483 | } | |
484 | } | |
255 | – | } |
485 | ||
486 | < | // make the sample and status out names |
487 | < | |
488 | < | strcpy( simnfo->sampleName, inFileName ); |
489 | < | char* endTest; |
490 | < | int nameLength = strlen( simnfo->sampleName ); |
491 | < | endTest = &(simnfo->sampleName[nameLength - 5]); |
492 | < | if( !strcmp( endTest, ".bass" ) ){ |
264 | < | strcpy( endTest, ".dump" ); |
265 | < | } |
266 | < | else if( !strcmp( endTest, ".BASS" ) ){ |
267 | < | strcpy( endTest, ".dump" ); |
268 | < | } |
269 | < | else{ |
270 | < | endTest = &(simnfo->sampleName[nameLength - 4]); |
271 | < | if( !strcmp( endTest, ".bss" ) ){ |
486 | > | // make the sample and status out names |
487 | > | |
488 | > | strcpy( simnfo->sampleName, inFileName ); |
489 | > | char* endTest; |
490 | > | int nameLength = strlen( simnfo->sampleName ); |
491 | > | endTest = &(simnfo->sampleName[nameLength - 5]); |
492 | > | if( !strcmp( endTest, ".bass" ) ){ |
493 | strcpy( endTest, ".dump" ); | |
494 | } | |
495 | < | else if( !strcmp( endTest, ".mdl" ) ){ |
495 | > | else if( !strcmp( endTest, ".BASS" ) ){ |
496 | strcpy( endTest, ".dump" ); | |
497 | } | |
498 | else{ | |
499 | < | strcat( simnfo->sampleName, ".dump" ); |
499 | > | endTest = &(simnfo->sampleName[nameLength - 4]); |
500 | > | if( !strcmp( endTest, ".bss" ) ){ |
501 | > | strcpy( endTest, ".dump" ); |
502 | > | } |
503 | > | else if( !strcmp( endTest, ".mdl" ) ){ |
504 | > | strcpy( endTest, ".dump" ); |
505 | > | } |
506 | > | else{ |
507 | > | strcat( simnfo->sampleName, ".dump" ); |
508 | > | } |
509 | } | |
510 | < | } |
511 | < | |
512 | < | strcpy( simnfo->statusName, inFileName ); |
513 | < | nameLength = strlen( simnfo->statusName ); |
514 | < | endTest = &(simnfo->statusName[nameLength - 5]); |
285 | < | if( !strcmp( endTest, ".bass" ) ){ |
286 | < | strcpy( endTest, ".stat" ); |
287 | < | } |
288 | < | else if( !strcmp( endTest, ".BASS" ) ){ |
289 | < | strcpy( endTest, ".stat" ); |
290 | < | } |
291 | < | else{ |
292 | < | endTest = &(simnfo->statusName[nameLength - 4]); |
293 | < | if( !strcmp( endTest, ".bss" ) ){ |
510 | > | |
511 | > | strcpy( simnfo->statusName, inFileName ); |
512 | > | nameLength = strlen( simnfo->statusName ); |
513 | > | endTest = &(simnfo->statusName[nameLength - 5]); |
514 | > | if( !strcmp( endTest, ".bass" ) ){ |
515 | strcpy( endTest, ".stat" ); | |
516 | } | |
517 | < | else if( !strcmp( endTest, ".mdl" ) ){ |
517 | > | else if( !strcmp( endTest, ".BASS" ) ){ |
518 | strcpy( endTest, ".stat" ); | |
519 | } | |
520 | else{ | |
521 | < | strcat( simnfo->statusName, ".stat" ); |
521 | > | endTest = &(simnfo->statusName[nameLength - 4]); |
522 | > | if( !strcmp( endTest, ".bss" ) ){ |
523 | > | strcpy( endTest, ".stat" ); |
524 | > | } |
525 | > | else if( !strcmp( endTest, ".mdl" ) ){ |
526 | > | strcpy( endTest, ".stat" ); |
527 | > | } |
528 | > | else{ |
529 | > | strcat( simnfo->statusName, ".stat" ); |
530 | > | } |
531 | } | |
532 | + | |
533 | + | #ifdef IS_MPI |
534 | } | |
535 | + | #endif // is_mpi |
536 | ||
537 | // set the status, sample, and themal kick times | |
538 | < | |
538 | > | |
539 | if( the_globals->haveSampleTime() ){ | |
540 | < | simnfo->sampleTime = the_globals->getSampleTime(); |
540 | > | simnfo->sampleTime = the_globals->getSampleTime(); |
541 | simnfo->statusTime = simnfo->sampleTime; | |
542 | simnfo->thermalTime = simnfo->sampleTime; | |
543 | } | |
544 | else{ | |
545 | < | simnfo->sampleTime = the_globals->getRunTime(); |
545 | > | simnfo->sampleTime = the_globals->getRunTime(); |
546 | simnfo->statusTime = simnfo->sampleTime; | |
547 | simnfo->thermalTime = simnfo->sampleTime; | |
548 | } | |
# | Line 325 | Line 558 | void SimSetup::createSim( void ){ | |
558 | // check for the temperature set flag | |
559 | ||
560 | if( the_globals->haveTempSet() ) simnfo->setTemp = the_globals->getTempSet(); | |
561 | < | |
562 | < | |
561 | > | |
562 | > | |
563 | // make the longe range forces and the integrator | |
564 | < | |
564 | > | |
565 | new AllLong( simnfo ); | |
566 | < | |
566 | > | |
567 | if( !strcmp( force_field, "TraPPE" ) ) new Verlet( *simnfo ); | |
568 | if( !strcmp( force_field, "DipoleTest" ) ) new Symplectic( simnfo ); | |
569 | if( !strcmp( force_field, "TraPPE_Ex" ) ) new Symplectic( simnfo ); | |
570 | } | |
571 | ||
572 | void SimSetup::makeAtoms( void ){ | |
573 | < | |
573 | > | |
574 | int i, j, k, index; | |
575 | double ux, uy, uz, uSqr, u; | |
576 | AtomStamp* current_atom; | |
577 | DirectionalAtom* dAtom; | |
578 | + | int molIndex, molStart, molEnd, nMemb, lMolIndex; |
579 | ||
580 | + | lMolIndex = 0; |
581 | + | molIndex = 0; |
582 | index = 0; | |
583 | for( i=0; i<n_components; i++ ){ | |
584 | < | |
584 | > | |
585 | for( j=0; j<components_nmol[i]; j++ ){ | |
350 | – | |
351 | – | for( k=0; k<comp_stamps[i]->getNAtoms(); k++ ){ |
352 | – | |
353 | – | current_atom = comp_stamps[i]->getAtom( k ); |
354 | – | if( current_atom->haveOrientation() ){ |
586 | ||
587 | < | dAtom = new DirectionalAtom; |
588 | < | simnfo->n_oriented++; |
589 | < | the_atoms[index] = dAtom; |
590 | < | |
591 | < | ux = current_atom->getOrntX(); |
592 | < | uy = current_atom->getOrntY(); |
593 | < | uz = current_atom->getOrntZ(); |
587 | > | #ifdef IS_MPI |
588 | > | if( mpiSim->getMyMolStart() <= molIndex && |
589 | > | molIndex <= mpiSim->getMyMolEnd() ){ |
590 | > | #endif // is_mpi |
591 | > | |
592 | > | molStart = index; |
593 | > | nMemb = comp_stamps[i]->getNAtoms(); |
594 | > | for( k=0; k<comp_stamps[i]->getNAtoms(); k++ ){ |
595 | ||
596 | < | uSqr = (ux * ux) + (uy * uy) + (uz * uz); |
596 | > | current_atom = comp_stamps[i]->getAtom( k ); |
597 | > | if( current_atom->haveOrientation() ){ |
598 | > | |
599 | > | dAtom = new DirectionalAtom(index); |
600 | > | simnfo->n_oriented++; |
601 | > | the_atoms[index] = dAtom; |
602 | > | |
603 | > | ux = current_atom->getOrntX(); |
604 | > | uy = current_atom->getOrntY(); |
605 | > | uz = current_atom->getOrntZ(); |
606 | > | |
607 | > | uSqr = (ux * ux) + (uy * uy) + (uz * uz); |
608 | > | |
609 | > | u = sqrt( uSqr ); |
610 | > | ux = ux / u; |
611 | > | uy = uy / u; |
612 | > | uz = uz / u; |
613 | > | |
614 | > | dAtom->setSUx( ux ); |
615 | > | dAtom->setSUy( uy ); |
616 | > | dAtom->setSUz( uz ); |
617 | > | } |
618 | > | else{ |
619 | > | the_atoms[index] = new GeneralAtom(index); |
620 | > | } |
621 | > | the_atoms[index]->setType( current_atom->getType() ); |
622 | > | the_atoms[index]->setIndex( index ); |
623 | ||
624 | < | u = sqrt( uSqr ); |
625 | < | ux = ux / u; |
368 | < | uy = uy / u; |
369 | < | uz = uz / u; |
370 | < | |
371 | < | dAtom->setSUx( ux ); |
372 | < | dAtom->setSUy( uy ); |
373 | < | dAtom->setSUz( uz ); |
624 | > | // increment the index and repeat; |
625 | > | index++; |
626 | } | |
375 | – | else{ |
376 | – | the_atoms[index] = new GeneralAtom; |
377 | – | } |
378 | – | the_atoms[index]->setType( current_atom->getType() ); |
379 | – | the_atoms[index]->setIndex( index ); |
627 | ||
628 | < | // increment the index and repeat; |
629 | < | index++; |
628 | > | molEnd = index -1; |
629 | > | the_molecules[lMolIndex].setNMembers( nMemb ); |
630 | > | the_molecules[lMolIndex].setStartAtom( molStart ); |
631 | > | the_molecules[lMolIndex].setEndAtom( molEnd ); |
632 | > | the_molecules[lMolIndex].setStampID( i ); |
633 | > | lMolIndex++; |
634 | > | |
635 | > | #ifdef IS_MPI |
636 | } | |
637 | + | #endif //is_mpi |
638 | + | |
639 | + | molIndex++; |
640 | } | |
641 | } | |
642 | < | |
642 | > | |
643 | the_ff->initializeAtoms(); | |
644 | } | |
645 | ||
646 | void SimSetup::makeBonds( void ){ | |
647 | ||
648 | < | int i, j, k, index, offset; |
648 | > | int i, j, k, index, offset, molIndex; |
649 | bond_pair* the_bonds; | |
650 | BondStamp* current_bond; | |
651 | ||
652 | the_bonds = new bond_pair[tot_bonds]; | |
653 | index = 0; | |
654 | offset = 0; | |
655 | + | molIndex = 0;g1 |
656 | + | |
657 | for( i=0; i<n_components; i++ ){ | |
658 | < | |
658 | > | |
659 | for( j=0; j<components_nmol[i]; j++ ){ | |
402 | – | |
403 | – | for( k=0; k<comp_stamps[i]->getNBonds(); k++ ){ |
404 | – | |
405 | – | current_bond = comp_stamps[i]->getBond( k ); |
406 | – | the_bonds[index].a = current_bond->getA() + offset; |
407 | – | the_bonds[index].b = current_bond->getB() + offset; |
660 | ||
661 | < | the_excludes[index].i = the_bonds[index].a; |
662 | < | the_excludes[index].j = the_bonds[index].b; |
663 | < | |
664 | < | // increment the index and repeat; |
665 | < | index++; |
666 | < | } |
667 | < | offset += comp_stamps[i]->getNAtoms(); |
668 | < | } |
661 | > | #ifdef IS_MPI |
662 | > | if( mpiSim->getMyMolStart() <= molIndex && |
663 | > | molIndex <= mpiSim->getMyMolEnd() ){ |
664 | > | #endif // is_mpi |
665 | > | |
666 | > | for( k=0; k<comp_stamps[i]->getNBonds(); k++ ){ |
667 | > | |
668 | > | current_bond = comp_stamps[i]->getBond( k ); |
669 | > | the_bonds[index].a = current_bond->getA() + offset; |
670 | > | the_bonds[index].b = current_bond->getB() + offset; |
671 | > | |
672 | > | the_excludes[index].i = the_bonds[index].a; |
673 | > | the_excludes[index].j = the_bonds[index].b; |
674 | > | |
675 | > | // increment the index and repeat; |
676 | > | index++; |
677 | > | } |
678 | > | offset += comp_stamps[i]->getNAtoms(); |
679 | > | |
680 | > | #ifdef IS_MPI |
681 | > | } |
682 | > | #endif is_mpi |
683 | > | |
684 | > | molIndex++; |
685 | > | } |
686 | } | |
687 | < | |
687 | > | |
688 | the_ff->initializeBonds( the_bonds ); | |
689 | } | |
690 | ||
691 | void SimSetup::makeBends( void ){ | |
692 | ||
693 | < | int i, j, k, index, offset; |
693 | > | int i, j, k, index, offset, molIndex; |
694 | bend_set* the_bends; | |
695 | BendStamp* current_bend; | |
696 | ||
697 | the_bends = new bend_set[tot_bends]; | |
698 | index = 0; | |
699 | offset = 0; | |
700 | + | molIndex = 0; |
701 | for( i=0; i<n_components; i++ ){ | |
702 | < | |
702 | > | |
703 | for( j=0; j<components_nmol[i]; j++ ){ | |
434 | – | |
435 | – | for( k=0; k<comp_stamps[i]->getNBends(); k++ ){ |
436 | – | |
437 | – | current_bend = comp_stamps[i]->getBend( k ); |
438 | – | the_bends[index].a = current_bend->getA() + offset; |
439 | – | the_bends[index].b = current_bend->getB() + offset; |
440 | – | the_bends[index].c = current_bend->getC() + offset; |
704 | ||
705 | < | the_excludes[index + tot_bonds].i = the_bends[index].a; |
706 | < | the_excludes[index + tot_bonds].j = the_bends[index].c; |
705 | > | #ifdef IS_MPI |
706 | > | if( mpiSim->getMyMolStart() <= molIndex && |
707 | > | molIndex <= mpiSim->getMyMolEnd() ){ |
708 | > | #endif // is_mpi |
709 | ||
710 | < | // increment the index and repeat; |
711 | < | index++; |
710 | > | for( k=0; k<comp_stamps[i]->getNBends(); k++ ){ |
711 | > | |
712 | > | current_bend = comp_stamps[i]->getBend( k ); |
713 | > | the_bends[index].a = current_bend->getA() + offset; |
714 | > | the_bends[index].b = current_bend->getB() + offset; |
715 | > | the_bends[index].c = current_bend->getC() + offset; |
716 | > | |
717 | > | the_excludes[index + tot_bonds].i = the_bends[index].a; |
718 | > | the_excludes[index + tot_bonds].j = the_bends[index].c; |
719 | > | |
720 | > | // increment the index and repeat; |
721 | > | index++; |
722 | > | } |
723 | > | offset += comp_stamps[i]->getNAtoms(); |
724 | > | |
725 | > | #ifdef IS_MPI |
726 | } | |
727 | < | offset += comp_stamps[i]->getNAtoms(); |
727 | > | #endif //is_mpi |
728 | > | |
729 | > | molIndex++; |
730 | } | |
731 | } | |
732 | < | |
732 | > | |
733 | the_ff->initializeBends( the_bends ); | |
734 | } | |
735 | ||
736 | void SimSetup::makeTorsions( void ){ | |
737 | ||
738 | < | int i, j, k, index, offset; |
738 | > | int i, j, k, index, offset, molIndex; |
739 | torsion_set* the_torsions; | |
740 | TorsionStamp* current_torsion; | |
741 | ||
742 | the_torsions = new torsion_set[tot_torsions]; | |
743 | index = 0; | |
744 | offset = 0; | |
745 | + | molIndex = 0; |
746 | for( i=0; i<n_components; i++ ){ | |
747 | < | |
747 | > | |
748 | for( j=0; j<components_nmol[i]; j++ ){ | |
749 | < | |
749 | > | |
750 | > | #ifdef IS_MPI |
751 | > | if( mpiSim->getMyMolStart() <= molIndex && |
752 | > | molIndex <= mpiSim->getMyMolEnd() ){ |
753 | > | #endif // is_mpi |
754 | > | |
755 | for( k=0; k<comp_stamps[i]->getNTorsions(); k++ ){ | |
756 | < | |
756 | > | |
757 | current_torsion = comp_stamps[i]->getTorsion( k ); | |
758 | the_torsions[index].a = current_torsion->getA() + offset; | |
759 | the_torsions[index].b = current_torsion->getB() + offset; | |
# | Line 480 | Line 767 | void SimSetup::makeTorsions( void ){ | |
767 | index++; | |
768 | } | |
769 | offset += comp_stamps[i]->getNAtoms(); | |
770 | + | |
771 | + | #ifdef IS_MPI |
772 | + | } |
773 | + | #endif //is_mpi |
774 | + | |
775 | + | molIndex++; |
776 | } | |
777 | } | |
778 | < | |
778 | > | |
779 | the_ff->initializeTorsions( the_torsions ); | |
780 | } | |
781 | ||
489 | – | void SimSetup::makeMolecules( void ){ |
490 | – | |
491 | – | //empy for now |
492 | – | } |
493 | – | |
782 | void SimSetup::initFromBass( void ){ | |
783 | ||
784 | int i, j, k; | |
# | Line 518 | Line 806 | void SimSetup::initFromBass( void ){ | |
806 | n_per_extra = (int)ceil( temp1 ); | |
807 | ||
808 | if( n_per_extra > 4){ | |
809 | < | std::cerr << "THere has been an error in constructing the non-complete lattice.\n"; |
810 | < | exit(8); |
809 | > | sprintf( painCave.errMsg, |
810 | > | "SimSetup error. There has been an error in constructing" |
811 | > | " the non-complete lattice.\n" ); |
812 | > | painCave.isFatal = 1; |
813 | > | simError(); |
814 | } | |
815 | } | |
816 | else{ | |
# | Line 528 | Line 819 | void SimSetup::initFromBass( void ){ | |
819 | celly = simnfo->box_y / temp3; | |
820 | cellz = simnfo->box_z / temp3; | |
821 | } | |
822 | < | |
822 | > | |
823 | current_mol = 0; | |
824 | current_comp_mol = 0; | |
825 | current_comp = 0; | |
826 | current_atom_ndx = 0; | |
827 | < | |
827 | > | |
828 | for( i=0; i < n_cells ; i++ ){ | |
829 | for( j=0; j < n_cells; j++ ){ | |
830 | for( k=0; k < n_cells; k++ ){ | |
831 | < | |
831 | > | |
832 | makeElement( i * cellx, | |
833 | j * celly, | |
834 | k * cellz ); | |
835 | < | |
835 | > | |
836 | makeElement( i * cellx + 0.5 * cellx, | |
837 | j * celly + 0.5 * celly, | |
838 | k * cellz ); | |
839 | < | |
839 | > | |
840 | makeElement( i * cellx, | |
841 | j * celly + 0.5 * celly, | |
842 | k * cellz + 0.5 * cellz ); | |
843 | < | |
843 | > | |
844 | makeElement( i * cellx + 0.5 * cellx, | |
845 | j * celly, | |
846 | k * cellz + 0.5 * cellz ); | |
# | Line 559 | Line 850 | void SimSetup::initFromBass( void ){ | |
850 | ||
851 | if( have_extra ){ | |
852 | done = 0; | |
853 | < | |
853 | > | |
854 | int start_ndx; | |
855 | for( i=0; i < (n_cells+1) && !done; i++ ){ | |
856 | for( j=0; j < (n_cells+1) && !done; j++ ){ | |
857 | < | |
857 | > | |
858 | if( i < n_cells ){ | |
859 | < | |
859 | > | |
860 | if( j < n_cells ){ | |
861 | start_ndx = n_cells; | |
862 | } | |
863 | else start_ndx = 0; | |
864 | } | |
865 | else start_ndx = 0; | |
866 | < | |
866 | > | |
867 | for( k=start_ndx; k < (n_cells+1) && !done; k++ ){ | |
868 | < | |
868 | > | |
869 | makeElement( i * cellx, | |
870 | j * celly, | |
871 | k * cellz ); | |
872 | done = ( current_mol >= tot_nmol ); | |
873 | < | |
873 | > | |
874 | if( !done && n_per_extra > 1 ){ | |
875 | makeElement( i * cellx + 0.5 * cellx, | |
876 | j * celly + 0.5 * celly, | |
877 | k * cellz ); | |
878 | done = ( current_mol >= tot_nmol ); | |
879 | } | |
880 | < | |
880 | > | |
881 | if( !done && n_per_extra > 2){ | |
882 | makeElement( i * cellx, | |
883 | j * celly + 0.5 * celly, | |
884 | k * cellz + 0.5 * cellz ); | |
885 | done = ( current_mol >= tot_nmol ); | |
886 | } | |
887 | < | |
887 | > | |
888 | if( !done && n_per_extra > 3){ | |
889 | makeElement( i * cellx + 0.5 * cellx, | |
890 | j * celly, | |
# | Line 604 | Line 895 | void SimSetup::initFromBass( void ){ | |
895 | } | |
896 | } | |
897 | } | |
898 | < | |
899 | < | |
898 | > | |
899 | > | |
900 | for( i=0; i<simnfo->n_atoms; i++ ){ | |
901 | simnfo->atoms[i]->set_vx( 0.0 ); | |
902 | simnfo->atoms[i]->set_vy( 0.0 ); | |
# | Line 621 | Line 912 | void SimSetup::makeElement( double x, double y, double | |
912 | double rotMat[3][3]; | |
913 | ||
914 | for( k=0; k<comp_stamps[current_comp]->getNAtoms(); k++ ){ | |
915 | < | |
915 | > | |
916 | current_atom = comp_stamps[current_comp]->getAtom( k ); | |
917 | if( !current_atom->havePosition() ){ | |
918 | < | std::cerr << "Component " << comp_stamps[current_comp]->getID() |
919 | < | << ", atom " << current_atom->getType() |
920 | < | << " does not have a position specified.\n" |
921 | < | << "The initialization routine is unable to give a start" |
922 | < | << " position.\n"; |
923 | < | exit(8); |
918 | > | sprintf( painCave.errMsg, |
919 | > | "SimSetup:initFromBass error.\n" |
920 | > | "\tComponent %s, atom %s does not have a position specified.\n" |
921 | > | "\tThe initialization routine is unable to give a start" |
922 | > | " position.\n", |
923 | > | comp_stamps[current_comp]->getID(), |
924 | > | current_atom->getType() ); |
925 | > | painCave.isFatal = 1; |
926 | > | simError(); |
927 | } | |
928 | < | |
928 | > | |
929 | the_atoms[current_atom_ndx]->setX( x + current_atom->getPosX() ); | |
930 | the_atoms[current_atom_ndx]->setY( y + current_atom->getPosY() ); | |
931 | the_atoms[current_atom_ndx]->setZ( z + current_atom->getPosZ() ); | |
932 | < | |
932 | > | |
933 | if( the_atoms[current_atom_ndx]->isDirectional() ){ | |
934 | < | |
934 | > | |
935 | dAtom = (DirectionalAtom *)the_atoms[current_atom_ndx]; | |
936 | < | |
936 | > | |
937 | rotMat[0][0] = 1.0; | |
938 | rotMat[0][1] = 0.0; | |
939 | rotMat[0][2] = 0.0; | |
# | Line 657 | Line 951 | void SimSetup::makeElement( double x, double y, double | |
951 | ||
952 | current_atom_ndx++; | |
953 | } | |
954 | < | |
954 | > | |
955 | current_mol++; | |
956 | current_comp_mol++; | |
957 | ||
958 | if( current_comp_mol >= components_nmol[current_comp] ){ | |
959 | < | |
959 | > | |
960 | current_comp_mol = 0; | |
961 | current_comp++; | |
962 | } |
– | Removed lines |
+ | Added lines |
< | Changed lines |
> | Changed lines |