ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/mdtools/interface_implementation/SimSetup.cpp
(Generate patch)

Comparing trunk/mdtools/interface_implementation/SimSetup.cpp (file contents):
Revision 144 by mmeineke, Thu Oct 17 21:59:12 2002 UTC vs.
Revision 215 by chuckv, Thu Dec 19 21:59:51 2002 UTC

# 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 24 | Line 31 | void SimSetup::parseFile( char* fileName ){
31  
32   void SimSetup::parseFile( char* fileName ){
33  
27 –  inFileName = fileName;
28 –  set_interface_stamps( stamps, globals );
34   #ifdef IS_MPI
35 <  mpiEventInit();
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 <  yacc_BASS( fileName );
44 >
45 >    yacc_BASS( fileName );
46 >
47   #ifdef IS_MPI
48 <  throwMPIEvent(NULL);
48 >    throwMPIEvent(NULL);
49 >  }
50 >  else receiveParse();
51   #endif
52  
53   }
# Line 41 | Line 57 | void SimSetup::receiveParse(void){
57  
58      set_interface_stamps( stamps, globals );
59      mpiEventInit();
60 +    MPIcheckPoint();
61      mpiEventLoop();
62  
63   }
# Line 52 | Line 69 | void SimSetup::createSim( void ){
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 75 | Line 93 | void SimSetup::createSim( void ){
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];
# Line 95 | Line 121 | void SimSetup::createSim( void ){
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 105 | Line 133 | void SimSetup::createSim( void ){
133      }
134    }
135    else{
136 <    std::cerr << "NOT A SUPPORTED FEATURE\n";
137 <
138 < //     tot_nmol = the_globals->getNMol();
139 <
140 < //     //we have the total number of molecules, now we check for molfractions
141 < //     for( i=0; i<n_components; i++ ){
142 <
143 < //       if( !the_components[i]->haveMolFraction() ){
144 <
145 < //      if( !the_components[i]->haveNMol() ){
146 < //        //we have a problem
147 < //        std::cerr << "SimSetup error. Neither molFraction nor "
148 < //                  << " nMol was given in component
149 <
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 >    
145 >    
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 138 | Line 215 | void SimSetup::createSim( void ){
215    tot_bends = 0;
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();
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();
222      tot_torsions += components_nmol[i] * comp_stamps[i]->getNTorsions();
223    }
224  
# Line 152 | Line 229 | void SimSetup::createSim( void ){
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  
234 +  
235 + #ifdef IS_MPI
236 +
237 +  // divide the molecules among processors here.
238 +  
239 +  mpiSim = new mpiSimulation( simnfo );
240 +  
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 <  the_atoms = new Atom*[tot_atoms];
301 <  the_molecules = new Molecule[tot_nmol];
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( tot_SRI ){
306 <    the_sris = new SRI*[tot_SRI];
307 <    the_excludes = new ex_pair[tot_SRI];
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
# Line 190 | Line 333 | void SimSetup::createSim( void ){
333      makeTorsions();
334    }
335  
193 –  //  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();
# Line 225 | 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    }
402  
403 + #ifdef IS_MPI
404 +  strcpy( checkPointMsg, "Box size set up" );
405 +  MPIcheckPoint();
406 + #endif // is_mpi
407  
247 – //   if( the_globals->haveInitialConfig() ){
248 – //        InitializeFromFile* fileInit;
249 – //     fileInit = new InitializeFromFile( the_globals->getInitialConfig() );
408  
251 – //     fileInit->read_xyz( simnfo ); // default velocities on
409  
410 < //     delete fileInit;
411 < //   }
412 < //   else{
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 <    initFromBass();
422 >   delete fileInit;
423 > }
424 > else{
425  
426 + #ifdef IS_MPI
427  
428 < //   }
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 < //   if( the_globals->haveFinalConfig() ){
263 < //     strcpy( simnfo->finalName, the_globals->getFinalConfig() );
264 < //   }
265 < //   else{
266 < //     strcpy( simnfo->finalName, inFileName );
267 < //     char* endTest;
268 < //     int nameLength = strlen( simnfo->finalName );
269 < //     endTest = &(simnfo->finalName[nameLength - 5]);
270 < //     if( !strcmp( endTest, ".bass" ) ){
271 < //       strcpy( endTest, ".eor" );
272 < //     }
273 < //     else if( !strcmp( endTest, ".BASS" ) ){
274 < //       strcpy( endTest, ".eor" );
275 < //     }
276 < //     else{
277 < //       endTest = &(simnfo->finalName[nameLength - 4]);
278 < //       if( !strcmp( endTest, ".bss" ) ){
279 < //      strcpy( endTest, ".eor" );
280 < //       }
281 < //       else if( !strcmp( endTest, ".mdl" ) ){
282 < //      strcpy( endTest, ".eor" );
283 < //       }
284 < //       else{
285 < //      strcat( simnfo->finalName, ".eor" );
286 < //       }
287 < //     }
288 < //   }
437 >  initFromBass();
438  
290 – //   // make the sample and status out names
439  
440 < //   strcpy( simnfo->sampleName, inFileName );
441 < //   char* endTest;
294 < //   int nameLength = strlen( simnfo->sampleName );
295 < //   endTest = &(simnfo->sampleName[nameLength - 5]);
296 < //   if( !strcmp( endTest, ".bass" ) ){
297 < //     strcpy( endTest, ".dump" );
298 < //   }
299 < //   else if( !strcmp( endTest, ".BASS" ) ){
300 < //     strcpy( endTest, ".dump" );
301 < //   }
302 < //   else{
303 < //     endTest = &(simnfo->sampleName[nameLength - 4]);
304 < //     if( !strcmp( endTest, ".bss" ) ){
305 < //       strcpy( endTest, ".dump" );
306 < //     }
307 < //     else if( !strcmp( endTest, ".mdl" ) ){
308 < //       strcpy( endTest, ".dump" );
309 < //     }
310 < //     else{
311 < //       strcat( simnfo->sampleName, ".dump" );
312 < //     }
313 < //   }
440 > #endif
441 > }
442  
443 < //   strcpy( simnfo->statusName, inFileName );
444 < //   nameLength = strlen( simnfo->statusName );
445 < //   endTest = &(simnfo->statusName[nameLength - 5]);
446 < //   if( !strcmp( endTest, ".bass" ) ){
319 < //     strcpy( endTest, ".stat" );
320 < //   }
321 < //   else if( !strcmp( endTest, ".BASS" ) ){
322 < //     strcpy( endTest, ".stat" );
323 < //   }
324 < //   else{
325 < //     endTest = &(simnfo->statusName[nameLength - 4]);
326 < //     if( !strcmp( endTest, ".bss" ) ){
327 < //       strcpy( endTest, ".stat" );
328 < //     }
329 < //     else if( !strcmp( endTest, ".mdl" ) ){
330 < //       strcpy( endTest, ".stat" );
331 < //     }
332 < //     else{
333 < //       strcat( simnfo->statusName, ".stat" );
334 < //     }
335 < //   }
443 > #ifdef IS_MPI
444 >  strcpy( checkPointMsg, "Successfully read in the initial configuration" );
445 >  MPIcheckPoint();
446 > #endif // is_mpi
447  
448  
449 <  // set the status, sample, and themal kick times
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 +    }
461 +    else{
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, ".BASS" ) ){
470 +        strcpy( endTest, ".eor" );
471 +      }
472 +      else{
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 +    }
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" ) ){
493 +      strcpy( endTest, ".dump" );
494 +    }
495 +    else if( !strcmp( endTest, ".BASS" ) ){
496 +      strcpy( endTest, ".dump" );
497 +    }
498 +    else{
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 +    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, ".BASS" ) ){
518 +      strcpy( endTest, ".stat" );
519 +    }
520 +    else{
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 +  
539    if( the_globals->haveSampleTime() ){
540      simnfo->sampleTime = the_globals->getSampleTime();
541      simnfo->statusTime = simnfo->sampleTime;
# Line 376 | Line 575 | void SimSetup::makeAtoms( void ){
575    double ux, uy, uz, uSqr, u;
576    AtomStamp* current_atom;
577    DirectionalAtom* dAtom;
578 <  int molIndex, molStart, molEnd, nMemb;
578 >  int molIndex, molStart, molEnd, nMemb, lMolIndex;
579  
580 <
580 >  lMolIndex = 0;
581    molIndex = 0;
582    index = 0;
583    for( i=0; i<n_components; i++ ){
584  
585      for( j=0; j<components_nmol[i]; j++ ){
586  
587 <      molStart = index;
588 <      nMemb = comp_stamps[i]->getNAtoms();
589 <      for( k=0; k<comp_stamps[i]->getNAtoms(); k++ ){
587 > #ifdef IS_MPI
588 >      if( mpiSim->getMyMolStart() <= molIndex &&
589 >          molIndex <= mpiSim->getMyMolEnd() ){
590 > #endif // is_mpi        
591  
592 <        current_atom = comp_stamps[i]->getAtom( k );
593 <        if( current_atom->haveOrientation() ){
594 <
595 <          dAtom = new DirectionalAtom;
596 <          simnfo->n_oriented++;
597 <          the_atoms[index] = dAtom;
598 <
599 <          ux = current_atom->getOrntX();
600 <          uy = current_atom->getOrntY();
601 <          uz = current_atom->getOrntZ();
602 <
603 <          uSqr = (ux * ux) + (uy * uy) + (uz * uz);
604 <
605 <          u = sqrt( uSqr );
606 <          ux = ux / u;
607 <          uy = uy / u;
608 <          uz = uz / u;
609 <
610 <          dAtom->setSUx( ux );
611 <          dAtom->setSUy( uy );
612 <          dAtom->setSUz( uz );
613 <        }
614 <        else{
615 <          the_atoms[index] = new GeneralAtom;
616 <        }
617 <        the_atoms[index]->setType( current_atom->getType() );
618 <        the_atoms[index]->setIndex( index );
619 <
620 <        // increment the index and repeat;
621 <        index++;
622 <      }
592 >        molStart = index;
593 >        nMemb = comp_stamps[i]->getNAtoms();
594 >        for( k=0; k<comp_stamps[i]->getNAtoms(); k++ ){
595 >          
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 >          // increment the index and repeat;
625 >          index++;
626 >        }
627 >        
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 <      molEnd = index -1;
636 <      the_molecules[molIndex].setNMembers( nMemb );
637 <      the_molecules[molIndex].setStartAtom( molStart );
638 <      the_molecules[molIndex].setEndAtom( molEnd );
635 > #ifdef IS_MPI
636 >      }
637 > #endif //is_mpi
638 >      
639        molIndex++;
429 –
640      }
641    }
642  
# Line 435 | Line 645 | void SimSetup::makeBonds( void ){
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  
659      for( j=0; j<components_nmol[i]; j++ ){
660  
661 <      for( k=0; k<comp_stamps[i]->getNBonds(); k++ ){
662 <
663 <        current_bond = comp_stamps[i]->getBond( k );
664 <        the_bonds[index].a = current_bond->getA() + offset;
665 <        the_bonds[index].b = current_bond->getB() + offset;
666 <
667 <        the_excludes[index].i = the_bonds[index].a;
668 <        the_excludes[index].j = the_bonds[index].b;
669 <
670 <        // increment the index and repeat;
671 <        index++;
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 <      offset += comp_stamps[i]->getNAtoms();
683 <    }
682 > #endif is_mpi
683 >      
684 >      molIndex++;
685 >    }      
686    }
687  
688    the_ff->initializeBonds( the_bonds );
# Line 467 | Line 690 | void SimSetup::makeBends( void ){
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  
703      for( j=0; j<components_nmol[i]; j++ ){
704  
705 <      for( k=0; k<comp_stamps[i]->getNBends(); k++ ){
705 > #ifdef IS_MPI
706 >      if( mpiSim->getMyMolStart() <= molIndex &&
707 >          molIndex <= mpiSim->getMyMolEnd() ){
708 > #endif // is_mpi        
709  
710 <        current_bend = comp_stamps[i]->getBend( k );
711 <        the_bends[index].a = current_bend->getA() + offset;
712 <        the_bends[index].b = current_bend->getB() + offset;
713 <        the_bends[index].c = current_bend->getC() + offset;
714 <
715 <        the_excludes[index + tot_bonds].i = the_bends[index].a;
716 <        the_excludes[index + tot_bonds].j = the_bends[index].c;
717 <
718 <        // increment the index and repeat;
719 <        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  
# Line 500 | Line 735 | void SimSetup::makeTorsions( void ){
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  
748      for( j=0; j<components_nmol[i]; j++ ){
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  
757          current_torsion = comp_stamps[i]->getTorsion( k );
# Line 526 | 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  
# Line 559 | 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 665 | Line 915 | void SimSetup::makeElement( double x, double y, double
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  
929      the_atoms[current_atom_ndx]->setX( x + current_atom->getPosX() );

Diff Legend

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