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 11 by mmeineke, Tue Jul 9 18:40:59 2002 UTC vs.
Revision 124 by chuckv, Mon Sep 30 20:35:42 2002 UTC

# Line 6 | Line 6
6   #include "parse_me.h"
7   #include "LRI.hpp"
8   #include "Integrator.hpp"
9 + #include "mpiInterface.h"
10  
11   SimSetup::SimSetup(){
12    stamps = new MakeStamps();
# Line 21 | Line 22 | void SimSetup::parseFile( char* fileName ){
22  
23    inFileName = fileName;
24    set_interface_stamps( stamps, globals );
25 + #ifdef MPI
26 +  mpiEventInit();
27 + #endif
28    yacc_BASS( fileName );
29 + #ifdef MPI
30 +  throwMPIEvent(NULL);
31 + #endif
32 +
33   }
34  
35 + #ifdef MPI
36 + void SimSetup::receiveParse(void){
37 +
38 +    set_interface_stamps( stamps, globals );
39 +    mpiEventInit();
40 +    mpiEventLoop();
41 +
42 + }
43 + #endif
44 +
45   void SimSetup::createSim( void ){
46  
47    MakeStamps *the_stamps;
# Line 43 | Line 61 | void SimSetup::createSim( void ){
61    n_components = the_globals->getNComponents();
62    strcpy( force_field, the_globals->getForceField() );
63    strcpy( ensemble, the_globals->getEnsemble() );
64 <  
64 >
65    if( !strcmp( force_field, "TraPPE" ) ) the_ff = new TraPPEFF();
66    else if( !strcmp( force_field, "DipoleTest" ) ) the_ff = new DipoleTestFF();
67    else if( !strcmp( force_field, "TraPPE_Ex" ) ) the_ff = new TraPPE_ExFF();
68    else{
69 <    std::cerr<< "SimSetup Error. Unrecognized force field -> "
69 >    std::cerr<< "SimSetup Error. Unrecognized force field -> "
70               << force_field << "\n";
71      exit(8);
72    }
# Line 57 | Line 75 | void SimSetup::createSim( void ){
75    the_components = the_globals->getComponents();
76    components_nmol = new int[n_components];
77    comp_stamps = new MoleculeStamp*[n_components];
78 <  
78 >
79    if( !the_globals->haveNMol() ){
80 <    // we don't have the total number of molecules, so we assume it is
80 >    // we don't have the total number of molecules, so we assume it is
81      // given in each component
82  
83      tot_nmol = 0;
84      for( i=0; i<n_components; i++ ){
85 <      
85 >
86        if( !the_components[i]->haveNMol() ){
87          // we have a problem
88          std::cerr << "SimSetup Error. No global NMol or component NMol"
# Line 78 | Line 96 | void SimSetup::createSim( void ){
96    }
97    else{
98      std::cerr << "NOT A SUPPORTED FEATURE\n";
99 <    
99 >
100   //     tot_nmol = the_globals->getNMol();
101 <    
101 >
102   //     //we have the total number of molecules, now we check for molfractions
103   //     for( i=0; i<n_components; i++ ){
104 <      
104 >
105   //       if( !the_components[i]->haveMolFraction() ){
106 <        
106 >
107   //      if( !the_components[i]->haveNMol() ){
108   //        //we have a problem
109   //        std::cerr << "SimSetup error. Neither molFraction nor "
# Line 97 | Line 115 | void SimSetup::createSim( void ){
115  
116    for( i=0; i<n_components; i++ ){
117  
118 <    comp_stamps[i] =
118 >    comp_stamps[i] =
119        the_stamps->getMolecule( the_components[i]->getType() );
120    }
121  
104 –  
122  
123 +
124    // caclulate the number of atoms, bonds, bends and torsions
125  
126    tot_atoms = 0;
# Line 110 | Line 128 | void SimSetup::createSim( void ){
128    tot_bends = 0;
129    tot_torsions = 0;
130    for( i=0; i<n_components; i++ ){
131 <    
131 >
132      tot_atoms += components_nmol[i] * comp_stamps[i]->getNAtoms();
133      tot_bonds += components_nmol[i] * comp_stamps[i]->getNBonds();
134      tot_bends += components_nmol[i] * comp_stamps[i]->getNBends();
135      tot_torsions += components_nmol[i] * comp_stamps[i]->getNTorsions();
136    }
137 <  
137 >
138    tot_SRI = tot_bonds + tot_bends + tot_torsions;
139 <  
139 >
140    simnfo->n_atoms = tot_atoms;
141    simnfo->n_bonds = tot_bonds;
142    simnfo->n_bends = tot_bends;
# Line 126 | Line 144 | void SimSetup::createSim( void ){
144    simnfo->n_SRI = tot_SRI;
145  
146    // create the atom and short range interaction arrays
147 <  
147 >
148    the_atoms = new Atom*[tot_atoms];
149 <  //  the_molecules = new Molecule[tot_nmol];
150 <  
151 <  
149 >  the_molecules = new Molecule[tot_nmol];
150 >
151 >
152    if( tot_SRI ){
153      the_sris = new SRI*[tot_SRI];
154      the_excludes = new ex_pair[tot_SRI];
# Line 143 | Line 161 | void SimSetup::createSim( void ){
161    simnfo->n_exclude = tot_SRI;
162    simnfo->excludes = the_excludes;
163  
164 +
165    // initialize the arrays
166 <  
166 >
167    the_ff->setSimInfo( simnfo );
168 <    
168 >
169    makeAtoms();
170  
171    if( tot_bonds ){
# Line 187 | Line 206 | void SimSetup::createSim( void ){
206      simnfo->box_z = the_globals->getBox();
207    }
208    else if( the_globals->haveDensity() ){
209 <    
209 >
210      double vol;
211      vol = (double)tot_nmol / the_globals->getDensity();
212      simnfo->box_x = pow( vol, ( 1.0 / 3.0 ) );
# Line 213 | Line 232 | void SimSetup::createSim( void ){
232      }
233      simnfo->box_z = the_globals->getBoxZ();
234    }
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
235  
236 <    delete fileInit;    
237 <  }
238 <  else{
236 >
237 > //   if( the_globals->haveInitialConfig() ){
238 > //        InitializeFromFile* fileInit;
239 > //     fileInit = new InitializeFromFile( the_globals->getInitialConfig() );
240 >
241 > //     fileInit->read_xyz( simnfo ); // default velocities on
242 >
243 > //     delete fileInit;
244 > //   }
245 > //   else{
246 >
247      initFromBass();
227 –  }
248  
229 –  if( the_globals->haveFinalConfig() ){
230 –    strcpy( simnfo->finalName, the_globals->getFinalConfig() );
231 –  }
232 –  else{
233 –    strcpy( simnfo->finalName, inFileName );
234 –    char* endTest;
235 –    int nameLength = strlen( simnfo->finalName );
236 –    endTest = &(simnfo->finalName[nameLength - 5]);
237 –    if( !strcmp( endTest, ".bass" ) ){
238 –      strcpy( endTest, ".eor" );
239 –    }
240 –    else if( !strcmp( endTest, ".BASS" ) ){
241 –      strcpy( endTest, ".eor" );
242 –    }
243 –    else{
244 –      endTest = &(simnfo->finalName[nameLength - 4]);
245 –      if( !strcmp( endTest, ".bss" ) ){
246 –        strcpy( endTest, ".eor" );
247 –      }
248 –      else if( !strcmp( endTest, ".mdl" ) ){
249 –        strcpy( endTest, ".eor" );
250 –      }
251 –      else{
252 –        strcat( simnfo->finalName, ".eor" );
253 –      }
254 –    }
255 –  }
256 –    
257 –  // make the sample and status out names
249  
250 <  strcpy( simnfo->sampleName, inFileName );
251 <  char* endTest;
252 <  int nameLength = strlen( simnfo->sampleName );
253 <  endTest = &(simnfo->sampleName[nameLength - 5]);
254 <  if( !strcmp( endTest, ".bass" ) ){
255 <    strcpy( endTest, ".dump" );
256 <  }
257 <  else if( !strcmp( endTest, ".BASS" ) ){
258 <    strcpy( endTest, ".dump" );
259 <  }
260 <  else{
261 <    endTest = &(simnfo->sampleName[nameLength - 4]);
262 <    if( !strcmp( endTest, ".bss" ) ){
263 <      strcpy( endTest, ".dump" );
264 <    }
265 <    else if( !strcmp( endTest, ".mdl" ) ){
266 <      strcpy( endTest, ".dump" );
267 <    }
268 <    else{
269 <      strcat( simnfo->sampleName, ".dump" );
270 <    }
271 <  }
272 <  
273 <  strcpy( simnfo->statusName, inFileName );
274 <  nameLength = strlen( simnfo->statusName );
275 <  endTest = &(simnfo->statusName[nameLength - 5]);
276 <  if( !strcmp( endTest, ".bass" ) ){
277 <    strcpy( endTest, ".stat" );
278 <  }
279 <  else if( !strcmp( endTest, ".BASS" ) ){
280 <    strcpy( endTest, ".stat" );
281 <  }
282 <  else{
283 <    endTest = &(simnfo->statusName[nameLength - 4]);
284 <    if( !strcmp( endTest, ".bss" ) ){
285 <      strcpy( endTest, ".stat" );
286 <    }
287 <    else if( !strcmp( endTest, ".mdl" ) ){
288 <      strcpy( endTest, ".stat" );
289 <    }
290 <    else{
291 <      strcat( simnfo->statusName, ".stat" );
292 <    }
293 <  }
294 <  
250 > //   }
251 >
252 > //   if( the_globals->haveFinalConfig() ){
253 > //     strcpy( simnfo->finalName, the_globals->getFinalConfig() );
254 > //   }
255 > //   else{
256 > //     strcpy( simnfo->finalName, inFileName );
257 > //     char* endTest;
258 > //     int nameLength = strlen( simnfo->finalName );
259 > //     endTest = &(simnfo->finalName[nameLength - 5]);
260 > //     if( !strcmp( endTest, ".bass" ) ){
261 > //       strcpy( endTest, ".eor" );
262 > //     }
263 > //     else if( !strcmp( endTest, ".BASS" ) ){
264 > //       strcpy( endTest, ".eor" );
265 > //     }
266 > //     else{
267 > //       endTest = &(simnfo->finalName[nameLength - 4]);
268 > //       if( !strcmp( endTest, ".bss" ) ){
269 > //      strcpy( endTest, ".eor" );
270 > //       }
271 > //       else if( !strcmp( endTest, ".mdl" ) ){
272 > //      strcpy( endTest, ".eor" );
273 > //       }
274 > //       else{
275 > //      strcat( simnfo->finalName, ".eor" );
276 > //       }
277 > //     }
278 > //   }
279 >
280 > //   // make the sample and status out names
281 >
282 > //   strcpy( simnfo->sampleName, inFileName );
283 > //   char* endTest;
284 > //   int nameLength = strlen( simnfo->sampleName );
285 > //   endTest = &(simnfo->sampleName[nameLength - 5]);
286 > //   if( !strcmp( endTest, ".bass" ) ){
287 > //     strcpy( endTest, ".dump" );
288 > //   }
289 > //   else if( !strcmp( endTest, ".BASS" ) ){
290 > //     strcpy( endTest, ".dump" );
291 > //   }
292 > //   else{
293 > //     endTest = &(simnfo->sampleName[nameLength - 4]);
294 > //     if( !strcmp( endTest, ".bss" ) ){
295 > //       strcpy( endTest, ".dump" );
296 > //     }
297 > //     else if( !strcmp( endTest, ".mdl" ) ){
298 > //       strcpy( endTest, ".dump" );
299 > //     }
300 > //     else{
301 > //       strcat( simnfo->sampleName, ".dump" );
302 > //     }
303 > //   }
304 >
305 > //   strcpy( simnfo->statusName, inFileName );
306 > //   nameLength = strlen( simnfo->statusName );
307 > //   endTest = &(simnfo->statusName[nameLength - 5]);
308 > //   if( !strcmp( endTest, ".bass" ) ){
309 > //     strcpy( endTest, ".stat" );
310 > //   }
311 > //   else if( !strcmp( endTest, ".BASS" ) ){
312 > //     strcpy( endTest, ".stat" );
313 > //   }
314 > //   else{
315 > //     endTest = &(simnfo->statusName[nameLength - 4]);
316 > //     if( !strcmp( endTest, ".bss" ) ){
317 > //       strcpy( endTest, ".stat" );
318 > //     }
319 > //     else if( !strcmp( endTest, ".mdl" ) ){
320 > //       strcpy( endTest, ".stat" );
321 > //     }
322 > //     else{
323 > //       strcat( simnfo->statusName, ".stat" );
324 > //     }
325 > //   }
326 >
327 >
328    // set the status, sample, and themal kick times
329  
330    if( the_globals->haveSampleTime() ){
331 <    simnfo->sampleTime = the_globals->getSampleTime();
331 >    simnfo->sampleTime = the_globals->getSampleTime();
332      simnfo->statusTime = simnfo->sampleTime;
333      simnfo->thermalTime = simnfo->sampleTime;
334    }
335    else{
336 <    simnfo->sampleTime = the_globals->getRunTime();
336 >    simnfo->sampleTime = the_globals->getRunTime();
337      simnfo->statusTime = simnfo->sampleTime;
338      simnfo->thermalTime = simnfo->sampleTime;
339    }
# Line 325 | Line 349 | void SimSetup::createSim( void ){
349    // check for the temperature set flag
350  
351    if( the_globals->haveTempSet() ) simnfo->setTemp = the_globals->getTempSet();
352 <  
353 <    
352 >
353 >
354    // make the longe range forces and the integrator
355 <  
355 >
356    new AllLong( simnfo );
357 <      
357 >
358    if( !strcmp( force_field, "TraPPE" ) ) new Verlet( *simnfo );
359    if( !strcmp( force_field, "DipoleTest" ) ) new Symplectic( simnfo );
360    if( !strcmp( force_field, "TraPPE_Ex" ) ) new Symplectic( simnfo );
361   }
362  
363   void SimSetup::makeAtoms( void ){
364 <  
364 >
365    int i, j, k, index;
366    double ux, uy, uz, uSqr, u;
367    AtomStamp* current_atom;
368    DirectionalAtom* dAtom;
369 +  int molIndex, molStart, molEnd, nMemb;
370  
371 +
372 +  molIndex = 0;
373    index = 0;
374    for( i=0; i<n_components; i++ ){
375 <    
375 >
376      for( j=0; j<components_nmol[i]; j++ ){
377 <      
377 >
378 >      molStart = index;
379 >      nMemb = comp_stamps[i]->getNAtoms();
380        for( k=0; k<comp_stamps[i]->getNAtoms(); k++ ){
381 <        
381 >
382          current_atom = comp_stamps[i]->getAtom( k );
383 <        if( current_atom->haveOrientation() ){
383 >        if( current_atom->haveOrientation() ){
384  
385            dAtom = new DirectionalAtom;
386            simnfo->n_oriented++;
387            the_atoms[index] = dAtom;
388 <      
388 >
389            ux = current_atom->getOrntX();
390            uy = current_atom->getOrntY();
391            uz = current_atom->getOrntZ();
392 <          
392 >
393            uSqr = (ux * ux) + (uy * uy) + (uz * uz);
394 <          
394 >
395            u = sqrt( uSqr );
396            ux = ux / u;
397            uy = uy / u;
398            uz = uz / u;
399 <          
399 >
400            dAtom->setSUx( ux );
401            dAtom->setSUy( uy );
402            dAtom->setSUz( uz );
# Line 377 | Line 406 | void SimSetup::makeAtoms( void ){
406          }
407          the_atoms[index]->setType( current_atom->getType() );
408          the_atoms[index]->setIndex( index );
409 <        
409 >
410          // increment the index and repeat;
411          index++;
412        }
413 +
414 +      molEnd = index -1;
415 +      the_molecules[molIndex].setNMembers( nMemb );
416 +      the_molecules[molIndex].setStartAtom( molStart );
417 +      the_molecules[molIndex].setEndAtom( molEnd );
418 +      molIndex++;
419 +
420      }
421    }
422 <  
422 >
423    the_ff->initializeAtoms();
424   }
425  
# Line 397 | Line 433 | void SimSetup::makeBonds( void ){
433    index = 0;
434    offset = 0;
435    for( i=0; i<n_components; i++ ){
436 <    
436 >
437      for( j=0; j<components_nmol[i]; j++ ){
438 <      
438 >
439        for( k=0; k<comp_stamps[i]->getNBonds(); k++ ){
440 <        
440 >
441          current_bond = comp_stamps[i]->getBond( k );
442          the_bonds[index].a = current_bond->getA() + offset;
443          the_bonds[index].b = current_bond->getB() + offset;
# Line 415 | Line 451 | void SimSetup::makeBonds( void ){
451        offset += comp_stamps[i]->getNAtoms();
452      }
453    }
454 <  
454 >
455    the_ff->initializeBonds( the_bonds );
456   }
457  
# Line 429 | Line 465 | void SimSetup::makeBends( void ){
465    index = 0;
466    offset = 0;
467    for( i=0; i<n_components; i++ ){
468 <    
468 >
469      for( j=0; j<components_nmol[i]; j++ ){
470 <      
470 >
471        for( k=0; k<comp_stamps[i]->getNBends(); k++ ){
472 <        
472 >
473          current_bend = comp_stamps[i]->getBend( k );
474          the_bends[index].a = current_bend->getA() + offset;
475          the_bends[index].b = current_bend->getB() + offset;
# Line 448 | Line 484 | void SimSetup::makeBends( void ){
484        offset += comp_stamps[i]->getNAtoms();
485      }
486    }
487 <  
487 >
488    the_ff->initializeBends( the_bends );
489   }
490  
# Line 462 | Line 498 | void SimSetup::makeTorsions( void ){
498    index = 0;
499    offset = 0;
500    for( i=0; i<n_components; i++ ){
501 <    
501 >
502      for( j=0; j<components_nmol[i]; j++ ){
503 <      
503 >
504        for( k=0; k<comp_stamps[i]->getNTorsions(); k++ ){
505 <        
505 >
506          current_torsion = comp_stamps[i]->getTorsion( k );
507          the_torsions[index].a = current_torsion->getA() + offset;
508          the_torsions[index].b = current_torsion->getB() + offset;
# Line 482 | Line 518 | void SimSetup::makeTorsions( void ){
518        offset += comp_stamps[i]->getNAtoms();
519      }
520    }
521 <  
521 >
522    the_ff->initializeTorsions( the_torsions );
523   }
524  
489 – void SimSetup::makeMolecules( void ){
490 –
491 –  //empy for now
492 – }
493 –
525   void SimSetup::initFromBass( void ){
526  
527    int i, j, k;
# Line 528 | Line 559 | void SimSetup::initFromBass( void ){
559      celly = simnfo->box_y / temp3;
560      cellz = simnfo->box_z / temp3;
561    }
562 <  
562 >
563    current_mol = 0;
564    current_comp_mol = 0;
565    current_comp = 0;
566    current_atom_ndx = 0;
567 <  
567 >
568    for( i=0; i < n_cells ; i++ ){
569      for( j=0; j < n_cells; j++ ){
570        for( k=0; k < n_cells; k++ ){
571 <        
571 >
572          makeElement( i * cellx,
573                       j * celly,
574                       k * cellz );
575 <        
575 >
576          makeElement( i * cellx + 0.5 * cellx,
577                       j * celly + 0.5 * celly,
578                       k * cellz );
579 <        
579 >
580          makeElement( i * cellx,
581                       j * celly + 0.5 * celly,
582                       k * cellz + 0.5 * cellz );
583 <        
583 >
584          makeElement( i * cellx + 0.5 * cellx,
585                       j * celly,
586                       k * cellz + 0.5 * cellz );
# Line 559 | Line 590 | void SimSetup::initFromBass( void ){
590  
591    if( have_extra ){
592      done = 0;
593 <    
593 >
594      int start_ndx;
595      for( i=0; i < (n_cells+1) && !done; i++ ){
596        for( j=0; j < (n_cells+1) && !done; j++ ){
597 <        
597 >
598          if( i < n_cells ){
599 <          
599 >
600            if( j < n_cells ){
601              start_ndx = n_cells;
602            }
603            else start_ndx = 0;
604          }
605          else start_ndx = 0;
606 <        
606 >
607          for( k=start_ndx; k < (n_cells+1) && !done; k++ ){
608 <          
608 >
609            makeElement( i * cellx,
610                         j * celly,
611                         k * cellz );
612            done = ( current_mol >= tot_nmol );
613 <          
613 >
614            if( !done && n_per_extra > 1 ){
615              makeElement( i * cellx + 0.5 * cellx,
616                           j * celly + 0.5 * celly,
617                           k * cellz );
618              done = ( current_mol >= tot_nmol );
619            }
620 <          
620 >
621            if( !done && n_per_extra > 2){
622              makeElement( i * cellx,
623                           j * celly + 0.5 * celly,
624                           k * cellz + 0.5 * cellz );
625              done = ( current_mol >= tot_nmol );
626            }
627 <          
627 >
628            if( !done && n_per_extra > 3){
629              makeElement( i * cellx + 0.5 * cellx,
630                           j * celly,
# Line 604 | Line 635 | void SimSetup::initFromBass( void ){
635        }
636      }
637    }
638 <  
639 <  
638 >
639 >
640    for( i=0; i<simnfo->n_atoms; i++ ){
641      simnfo->atoms[i]->set_vx( 0.0 );
642      simnfo->atoms[i]->set_vy( 0.0 );
# Line 621 | Line 652 | void SimSetup::makeElement( double x, double y, double
652    double rotMat[3][3];
653  
654    for( k=0; k<comp_stamps[current_comp]->getNAtoms(); k++ ){
655 <    
655 >
656      current_atom = comp_stamps[current_comp]->getAtom( k );
657      if( !current_atom->havePosition() ){
658        std::cerr << "Component " << comp_stamps[current_comp]->getID()
# Line 631 | Line 662 | void SimSetup::makeElement( double x, double y, double
662                  << " position.\n";
663        exit(8);
664      }
665 <    
665 >
666      the_atoms[current_atom_ndx]->setX( x + current_atom->getPosX() );
667      the_atoms[current_atom_ndx]->setY( y + current_atom->getPosY() );
668      the_atoms[current_atom_ndx]->setZ( z + current_atom->getPosZ() );
669 <    
669 >
670      if( the_atoms[current_atom_ndx]->isDirectional() ){
671 <      
671 >
672        dAtom = (DirectionalAtom *)the_atoms[current_atom_ndx];
673 <      
673 >
674        rotMat[0][0] = 1.0;
675        rotMat[0][1] = 0.0;
676        rotMat[0][2] = 0.0;
# Line 657 | Line 688 | void SimSetup::makeElement( double x, double y, double
688  
689      current_atom_ndx++;
690    }
691 <  
691 >
692    current_mol++;
693    current_comp_mol++;
694  
695    if( current_comp_mol >= components_nmol[current_comp] ){
696 <    
696 >
697      current_comp_mol = 0;
698      current_comp++;
699    }

Diff Legend

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