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

Comparing trunk/OOPSE/libmdtools/SimSetup.cpp (file contents):
Revision 615 by mmeineke, Tue Jul 15 17:57:04 2003 UTC vs.
Revision 616 by mmeineke, Tue Jul 15 18:52:16 2003 UTC

# Line 78 | Line 78 | void SimSetup::createSim( void ){
78  
79   void SimSetup::createSim( void ){
80  
81  MakeStamps *the_stamps;
82  Globals* the_globals;
81    int i, j, k, globalAtomIndex;
82    
83    // gather all of the information from the Bass file
# Line 90 | Line 88 | void SimSetup::createSim( void ){
88  
89    sysObjectsCreation();
90  
93
94
95  // initialize the arrays
96
97
98
99  makeMolecules();
100  info->identArray = new int[info->n_atoms];
101  for(i=0; i<info->n_atoms; i++){
102    info->identArray[i] = the_atoms[i]->getIdent();
103  }
104  
105
91    // check on the post processing info
92    
93    finalInfoCheck();
94  
110
111
112
95    // initialize the system coordinates
96  
97    initSystemCoords();
# Line 119 | Line 101 | void SimSetup::createSim( void ){
101  
102    makeOutNames();
103    
122
123  
124
125  
126
127  
128
129
104    // make the integrator
105    
106 +  makeIntegrator();
107    
133  NVT*  myNVT = NULL;
134  NPTi* myNPTi = NULL;
135  NPTf* myNPTf = NULL;
136  NPTim* myNPTim = NULL;
137  NPTfm* myNPTfm = NULL;
138
139  switch( ensembleCase ){
140
141  case NVE_ENS:
142    new NVE( info, the_ff );
143    break;
144
145  case NVT_ENS:
146    myNVT = new NVT( info, the_ff );
147    myNVT->setTargetTemp(the_globals->getTargetTemp());
148
149    if (the_globals->haveTauThermostat())
150      myNVT->setTauThermostat(the_globals->getTauThermostat());
151
152    else {
153      sprintf( painCave.errMsg,
154               "SimSetup error: If you use the NVT\n"
155               "    ensemble, you must set tauThermostat.\n");
156      painCave.isFatal = 1;
157      simError();
158    }
159    break;
160
161  case NPTi_ENS:
162    myNPTi = new NPTi( info, the_ff );
163    myNPTi->setTargetTemp( the_globals->getTargetTemp() );
164
165    if (the_globals->haveTargetPressure())
166      myNPTi->setTargetPressure(the_globals->getTargetPressure());
167    else {
168      sprintf( painCave.errMsg,
169               "SimSetup error: If you use a constant pressure\n"
170               "    ensemble, you must set targetPressure in the BASS file.\n");
171      painCave.isFatal = 1;
172      simError();
173    }
174    
175    if( the_globals->haveTauThermostat() )
176      myNPTi->setTauThermostat( the_globals->getTauThermostat() );
177    else{
178      sprintf( painCave.errMsg,
179               "SimSetup error: If you use an NPT\n"
180               "    ensemble, you must set tauThermostat.\n");
181      painCave.isFatal = 1;
182      simError();
183    }
184
185    if( the_globals->haveTauBarostat() )
186      myNPTi->setTauBarostat( the_globals->getTauBarostat() );
187    else{
188      sprintf( painCave.errMsg,
189               "SimSetup error: If you use an NPT\n"
190               "    ensemble, you must set tauBarostat.\n");
191      painCave.isFatal = 1;
192      simError();
193    }
194    break;
195
196  case NPTf_ENS:
197    myNPTf = new NPTf( info, the_ff );
198    myNPTf->setTargetTemp( the_globals->getTargetTemp());
199
200    if (the_globals->haveTargetPressure())
201      myNPTf->setTargetPressure(the_globals->getTargetPressure());
202    else {
203      sprintf( painCave.errMsg,
204               "SimSetup error: If you use a constant pressure\n"
205               "    ensemble, you must set targetPressure in the BASS file.\n");
206      painCave.isFatal = 1;
207      simError();
208    }    
209
210    if( the_globals->haveTauThermostat() )
211      myNPTf->setTauThermostat( the_globals->getTauThermostat() );
212    else{
213      sprintf( painCave.errMsg,
214               "SimSetup error: If you use an NPT\n"
215               "    ensemble, you must set tauThermostat.\n");
216      painCave.isFatal = 1;
217      simError();
218    }
219
220    if( the_globals->haveTauBarostat() )
221      myNPTf->setTauBarostat( the_globals->getTauBarostat() );
222    else{
223      sprintf( painCave.errMsg,
224               "SimSetup error: If you use an NPT\n"
225               "    ensemble, you must set tauBarostat.\n");
226      painCave.isFatal = 1;
227      simError();
228    }
229    break;
230    
231  case NPTim_ENS:
232    myNPTim = new NPTim( info, the_ff );
233    myNPTim->setTargetTemp( the_globals->getTargetTemp());
234
235    if (the_globals->haveTargetPressure())
236      myNPTim->setTargetPressure(the_globals->getTargetPressure());
237    else {
238      sprintf( painCave.errMsg,
239               "SimSetup error: If you use a constant pressure\n"
240               "    ensemble, you must set targetPressure in the BASS file.\n");
241      painCave.isFatal = 1;
242      simError();
243    }
244    
245    if( the_globals->haveTauThermostat() )
246      myNPTim->setTauThermostat( the_globals->getTauThermostat() );
247    else{
248      sprintf( painCave.errMsg,
249               "SimSetup error: If you use an NPT\n"
250               "    ensemble, you must set tauThermostat.\n");
251      painCave.isFatal = 1;
252      simError();
253    }
254
255    if( the_globals->haveTauBarostat() )
256      myNPTim->setTauBarostat( the_globals->getTauBarostat() );
257    else{
258      sprintf( painCave.errMsg,
259               "SimSetup error: If you use an NPT\n"
260               "    ensemble, you must set tauBarostat.\n");
261      painCave.isFatal = 1;
262      simError();
263    }
264    break;
265
266  case NPTfm_ENS:
267    myNPTfm = new NPTfm( info, the_ff );
268    myNPTfm->setTargetTemp( the_globals->getTargetTemp());
269
270    if (the_globals->haveTargetPressure())
271      myNPTfm->setTargetPressure(the_globals->getTargetPressure());
272    else {
273      sprintf( painCave.errMsg,
274               "SimSetup error: If you use a constant pressure\n"
275               "    ensemble, you must set targetPressure in the BASS file.\n");
276      painCave.isFatal = 1;
277      simError();
278    }
279    
280    if( the_globals->haveTauThermostat() )
281      myNPTfm->setTauThermostat( the_globals->getTauThermostat() );
282    else{
283      sprintf( painCave.errMsg,
284               "SimSetup error: If you use an NPT\n"
285               "    ensemble, you must set tauThermostat.\n");
286      painCave.isFatal = 1;
287      simError();
288    }
289
290    if( the_globals->haveTauBarostat() )
291      myNPTfm->setTauBarostat( the_globals->getTauBarostat() );
292    else{
293      sprintf( painCave.errMsg,
294               "SimSetup error: If you use an NPT\n"
295               "    ensemble, you must set tauBarostat.\n");
296      painCave.isFatal = 1;
297      simError();
298    }
299    break;
300
301  default:
302    sprintf( painCave.errMsg,
303             "SimSetup Error. Unrecognized ensemble in case statement.\n");
304    painCave.isFatal = 1;
305    simError();
306  }
307
308
108   #ifdef IS_MPI
109    mpiSim->mpiRefresh();
110   #endif
111  
112    // initialize the Fortran
113  
114 +  initFortran();
115  
316  info->refreshSim();
317  
318  if( !strcmp( info->mixingRule, "standard") ){
319    the_ff->initForceField( LB_MIXING_RULE );
320  }
321  else if( !strcmp( info->mixingRule, "explicit") ){
322    the_ff->initForceField( EXPLICIT_MIXING_RULE );
323  }
324  else{
325    sprintf( painCave.errMsg,
326             "SimSetup Error: unknown mixing rule -> \"%s\"\n",
327             info->mixingRule );
328    painCave.isFatal = 1;
329    simError();
330  }
116  
117  
333 #ifdef IS_MPI
334  strcpy( checkPointMsg,
335          "Successfully intialized the mixingRule for Fortran." );
336  MPIcheckPoint();
337 #endif // is_mpi
118   }
119  
120  
121   void SimSetup::makeMolecules( void ){
122  
123    int i, j, exI, exJ, tempEx, stampID, atomOffset, excludeOffset;
124 <  molInit info;
124 >  molInit molInfo;
125    DirectionalAtom* dAtom;
126    LinkedAssign* extras;
127    LinkedAssign* current_extra;
# Line 370 | Line 150 | void SimSetup::makeMolecules( void ){
150      
151      stampID = the_molecules[i].getStampID();
152  
153 <    info.nAtoms    = comp_stamps[stampID]->getNAtoms();
154 <    info.nBonds    = comp_stamps[stampID]->getNBonds();
155 <    info.nBends    = comp_stamps[stampID]->getNBends();
156 <    info.nTorsions = comp_stamps[stampID]->getNTorsions();
157 <    info.nExcludes = info.nBonds + info.nBends + info.nTorsions;
153 >    molInfo.nAtoms    = comp_stamps[stampID]->getNAtoms();
154 >    molInfo.nBonds    = comp_stamps[stampID]->getNBonds();
155 >    molInfo.nBends    = comp_stamps[stampID]->getNBends();
156 >    molInfo.nTorsions = comp_stamps[stampID]->getNTorsions();
157 >    molInfo.nExcludes = molInfo.nBonds + molInfo.nBends + molInfo.nTorsions;
158  
159 <    info.myAtoms = &the_atoms[atomOffset];
160 <    info.myExcludes = &the_excludes[excludeOffset];
161 <    info.myBonds = new Bond*[info.nBonds];
162 <    info.myBends = new Bend*[info.nBends];
163 <    info.myTorsions = new Torsion*[info.nTorsions];
159 >    molInfo.myAtoms = &the_atoms[atomOffset];
160 >    molInfo.myExcludes = &the_excludes[excludeOffset];
161 >    molInfo.myBonds = new Bond*[molInfo.nBonds];
162 >    molInfo.myBends = new Bend*[molInfo.nBends];
163 >    molInfo.myTorsions = new Torsion*[molInfo.nTorsions];
164  
165 <    theBonds = new bond_pair[info.nBonds];
166 <    theBends = new bend_set[info.nBends];
167 <    theTorsions = new torsion_set[info.nTorsions];
165 >    theBonds = new bond_pair[molInfo.nBonds];
166 >    theBends = new bend_set[molInfo.nBends];
167 >    theTorsions = new torsion_set[molInfo.nTorsions];
168      
169      // make the Atoms
170      
171 <    for(j=0; j<info.nAtoms; j++){
171 >    for(j=0; j<molInfo.nAtoms; j++){
172        
173        currentAtom = comp_stamps[stampID]->getAtom( j );
174        if( currentAtom->haveOrientation() ){
175          
176          dAtom = new DirectionalAtom(j + atomOffset);
177          info->n_oriented++;
178 <        info.myAtoms[j] = dAtom;
178 >        molInfo.myAtoms[j] = dAtom;
179          
180          ux = currentAtom->getOrntX();
181          uy = currentAtom->getOrntY();
# Line 413 | Line 193 | void SimSetup::makeMolecules( void ){
193          dAtom->setSUz( uz );
194        }
195        else{
196 <        info.myAtoms[j] = new GeneralAtom(j + atomOffset);
196 >        molInfo.myAtoms[j] = new GeneralAtom(j + atomOffset);
197        }
198 <      info.myAtoms[j]->setType( currentAtom->getType() );
198 >      molInfo.myAtoms[j]->setType( currentAtom->getType() );
199      
200   #ifdef IS_MPI
201        
202 <      info.myAtoms[j]->setGlobalIndex( globalIndex[j+atomOffset] );
202 >      molInfo.myAtoms[j]->setGlobalIndex( globalIndex[j+atomOffset] );
203        
204   #endif // is_mpi
205      }
206      
207      // make the bonds
208 <    for(j=0; j<info.nBonds; j++){
208 >    for(j=0; j<molInfo.nBonds; j++){
209        
210        currentBond = comp_stamps[stampID]->getBond( j );
211        theBonds[j].a = currentBond->getA() + atomOffset;
# Line 452 | Line 232 | void SimSetup::makeMolecules( void ){
232        the_excludes[j+excludeOffset]->setPair( (exI+1), (exJ+1) );
233   #endif  //is_mpi
234      }
235 <    excludeOffset += info.nBonds;
235 >    excludeOffset += molInfo.nBonds;
236  
237      //make the bends
238 <    for(j=0; j<info.nBends; j++){
238 >    for(j=0; j<molInfo.nBends; j++){
239        
240        currentBend = comp_stamps[stampID]->getBend( j );
241        theBends[j].a = currentBend->getA() + atomOffset;
# Line 538 | Line 318 | void SimSetup::makeMolecules( void ){
318        the_excludes[j+excludeOffset]->setPair( (exI+1), (exJ+1) );
319   #endif  //is_mpi
320      }
321 <    excludeOffset += info.nBends;
321 >    excludeOffset += molInfo.nBends;
322  
323 <    for(j=0; j<info.nTorsions; j++){
323 >    for(j=0; j<molInfo.nTorsions; j++){
324        
325        currentTorsion = comp_stamps[stampID]->getTorsion( j );
326        theTorsions[j].a = currentTorsion->getA() + atomOffset;
# Line 568 | Line 348 | void SimSetup::makeMolecules( void ){
348        the_excludes[j+excludeOffset]->setPair( (exI+1), (exJ+1) );
349   #endif  //is_mpi
350      }
351 <    excludeOffset += info.nTorsions;
351 >    excludeOffset += molInfo.nTorsions;
352  
353      
354      // send the arrays off to the forceField for init.
355  
356 <    the_ff->initializeAtoms( info.nAtoms, info.myAtoms );
357 <    the_ff->initializeBonds( info.nBonds, info.myBonds, theBonds );
358 <    the_ff->initializeBends( info.nBends, info.myBends, theBends );
359 <    the_ff->initializeTorsions( info.nTorsions, info.myTorsions, theTorsions );
356 >    the_ff->initializeAtoms( molInfo.nAtoms, molInfo.myAtoms );
357 >    the_ff->initializeBonds( molInfo.nBonds, molInfo.myBonds, theBonds );
358 >    the_ff->initializeBends( molInfo.nBends, molInfo.myBends, theBends );
359 >    the_ff->initializeTorsions( molInfo.nTorsions, molInfo.myTorsions, theTorsions );
360  
361  
362 <    the_molecules[i].initialize( info );
362 >    the_molecules[i].initialize( molInfo );
363  
364  
365 <    atomOffset += info.nAtoms;
365 >    atomOffset += molInfo.nAtoms;
366      delete[] theBonds;
367      delete[] theBends;
368      delete[] theTorsions;
# Line 784 | Line 564 | void SimSetup::gatherInfo( void ){
564  
565  
566   void SimSetup::gatherInfo( void ){
567 +  int i,j,k;
568  
569    ensembleCase = -1;
570    ffCase = -1;
571  
572    // get the stamps and globals;
573 <  the_stamps = stamps;
574 <  the_globals = globals;
573 >  stamps = stamps;
574 >  globals = globals;
575  
576    // set the easy ones first
577 <  info->target_temp = the_globals->getTargetTemp();
578 <  info->dt = the_globals->getDt();
579 <  info->run_time = the_globals->getRunTime();
580 <  n_components = the_globals->getNComponents();
577 >  info->target_temp = globals->getTargetTemp();
578 >  info->dt = globals->getDt();
579 >  info->run_time = globals->getRunTime();
580 >  n_components = globals->getNComponents();
581  
582  
583    // get the forceField
584  
585 <  strcpy( force_field, the_globals->getForceField() );
585 >  strcpy( force_field, globals->getForceField() );
586  
587    if( !strcasecmp( force_field, "DUFF" )) ffCase = FF_DUFF;
588    else if( !strcasecmp( force_field, "LJ" )) ffCase = FF_LJ;
# Line 815 | Line 596 | void SimSetup::gatherInfo( void ){
596  
597    // get the ensemble
598  
599 <  strcpy( ensemble, the_globals->getEnsemble() );
599 >  strcpy( ensemble, globals->getEnsemble() );
600  
601    if( !strcasecmp( ensemble, "NVE" ))      ensembleCase = NVE_ENS;
602    else if( !strcasecmp( ensemble, "NVT" )) ensembleCase = NVT_ENS;
# Line 838 | Line 619 | void SimSetup::gatherInfo( void ){
619  
620    // get the mixing rule
621  
622 <  strcpy( info->mixingRule, the_globals->getMixingRule() );
623 <  info->usePBC = the_globals->getPBC();
622 >  strcpy( info->mixingRule, globals->getMixingRule() );
623 >  info->usePBC = globals->getPBC();
624          
625    
626    // get the components and calculate the tot_nMol and indvidual n_mol
627  
628 <  the_components = the_globals->getComponents();
628 >  the_components = globals->getComponents();
629    components_nmol = new int[n_components];
630  
631  
632 <  if( !the_globals->haveNMol() ){
632 >  if( !globals->haveNMol() ){
633      // we don't have the total number of molecules, so we assume it is
634      // given in each component
635  
# Line 881 | Line 662 | void SimSetup::gatherInfo( void ){
662  
663    // set the status, sample, and thermal kick times
664    
665 <  if( the_globals->haveSampleTime() ){
666 <    info->sampleTime = the_globals->getSampleTime();
665 >  if( globals->haveSampleTime() ){
666 >    info->sampleTime = globals->getSampleTime();
667      info->statusTime = info->sampleTime;
668      info->thermalTime = info->sampleTime;
669    }
670    else{
671 <    info->sampleTime = the_globals->getRunTime();
671 >    info->sampleTime = globals->getRunTime();
672      info->statusTime = info->sampleTime;
673      info->thermalTime = info->sampleTime;
674    }
675  
676 <  if( the_globals->haveStatusTime() ){
677 <    info->statusTime = the_globals->getStatusTime();
676 >  if( globals->haveStatusTime() ){
677 >    info->statusTime = globals->getStatusTime();
678    }
679  
680 <  if( the_globals->haveThermalTime() ){
681 <    info->thermalTime = the_globals->getThermalTime();
680 >  if( globals->haveThermalTime() ){
681 >    info->thermalTime = globals->getThermalTime();
682    }
683  
684    // check for the temperature set flag
685  
686 <  if( the_globals->haveTempSet() ) info->setTemp = the_globals->getTempSet();
686 >  if( globals->haveTempSet() ) info->setTemp = globals->getTempSet();
687  
688    // get some of the tricky things that may still be in the globals
689  
690    double boxVector[3];
691 <  if( the_globals->haveBox() ){
692 <    boxVector[0] = the_globals->getBox();
693 <    boxVector[1] = the_globals->getBox();
694 <    boxVector[2] = the_globals->getBox();
691 >  if( globals->haveBox() ){
692 >    boxVector[0] = globals->getBox();
693 >    boxVector[1] = globals->getBox();
694 >    boxVector[2] = globals->getBox();
695      
696      info->setBox( boxVector );
697    }
698 <  else if( the_globals->haveDensity() ){
698 >  else if( globals->haveDensity() ){
699  
700      double vol;
701 <    vol = (double)tot_nmol / the_globals->getDensity();
701 >    vol = (double)tot_nmol / globals->getDensity();
702       boxVector[0] = pow( vol, ( 1.0 / 3.0 ) );
703       boxVector[1] = boxVector[0];
704       boxVector[2] = boxVector[0];
# Line 925 | Line 706 | void SimSetup::gatherInfo( void ){
706      info->setBox( boxVector );
707    }
708    else{
709 <    if( !the_globals->haveBoxX() ){
709 >    if( !globals->haveBoxX() ){
710        sprintf( painCave.errMsg,
711                 "SimSetup error, no periodic BoxX size given.\n" );
712        painCave.isFatal = 1;
713        simError();
714      }
715 <    boxVector[0] = the_globals->getBoxX();
715 >    boxVector[0] = globals->getBoxX();
716  
717 <    if( !the_globals->haveBoxY() ){
717 >    if( !globals->haveBoxY() ){
718        sprintf( painCave.errMsg,
719                 "SimSetup error, no periodic BoxY size given.\n" );
720        painCave.isFatal = 1;
721        simError();
722      }
723 <    boxVector[1] = the_globals->getBoxY();
723 >    boxVector[1] = globals->getBoxY();
724  
725 <    if( !the_globals->haveBoxZ() ){
725 >    if( !globals->haveBoxZ() ){
726        sprintf( painCave.errMsg,
727                 "SimSetup error, no periodic BoxZ size given.\n" );
728        painCave.isFatal = 1;
729        simError();
730      }
731 <    boxVector[2] = the_globals->getBoxZ();
731 >    boxVector[2] = globals->getBoxZ();
732  
733      info->setBox( boxVector );
734    }
# Line 977 | Line 758 | void SimSetup::finalInfoCheck( void ){
758    }
759    
760   #ifdef IS_MPI
761 <  int myUse = usesDipoles
762 <  MPI_Allreduce( &myUse, &UsesDipoles, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD );
761 >  int myUse = usesDipoles;
762 >  MPI_Allreduce( &myUse, &usesDipoles, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD );
763   #endif //is_mpi
764  
765  
766 <  if (the_globals->getUseRF() ) {
766 >  if (globals->getUseRF() ) {
767      info->useReactionField = 1;
768      
769 <    if( !the_globals->haveECR() ){
769 >    if( !globals->haveECR() ){
770        sprintf( painCave.errMsg,
771                 "SimSetup Warning: using default value of 1/2 the smallest "
772                 "box length for the electrostaticCutoffRadius.\n"
# Line 998 | Line 779 | void SimSetup::finalInfoCheck( void ){
779        if (info->boxLz <= smallest) smallest = info->boxLz;
780        info->ecr = 0.5 * smallest;
781      } else {
782 <      info->ecr        = the_globals->getECR();
782 >      info->ecr        = globals->getECR();
783      }
784  
785 <    if( !the_globals->haveEST() ){
785 >    if( !globals->haveEST() ){
786        sprintf( painCave.errMsg,
787                 "SimSetup Warning: using default value of 0.05 * the "
788                 "electrostaticCutoffRadius for the electrostaticSkinThickness\n"
# Line 1010 | Line 791 | void SimSetup::finalInfoCheck( void ){
791        simError();
792        info->est = 0.05 * info->ecr;
793      } else {
794 <      info->est        = the_globals->getEST();
794 >      info->est        = globals->getEST();
795      }
796      
797 <    if(!the_globals->haveDielectric() ){
797 >    if(!globals->haveDielectric() ){
798        sprintf( painCave.errMsg,
799                 "SimSetup Error: You are trying to use Reaction Field without"
800                 "setting a dielectric constant!\n"
# Line 1021 | Line 802 | void SimSetup::finalInfoCheck( void ){
802        painCave.isFatal = 1;
803        simError();
804      }
805 <    info->dielectric = the_globals->getDielectric();  
805 >    info->dielectric = globals->getDielectric();  
806    }
807    else {
808      if (usesDipoles) {
809        
810 <      if( !the_globals->haveECR() ){
810 >      if( !globals->haveECR() ){
811          sprintf( painCave.errMsg,
812                   "SimSetup Warning: using default value of 1/2 the smallest "
813                   "box length for the electrostaticCutoffRadius.\n"
# Line 1039 | Line 820 | void SimSetup::finalInfoCheck( void ){
820          if (info->boxLz <= smallest) smallest = info->boxLz;
821          info->ecr = 0.5 * smallest;
822        } else {
823 <        info->ecr        = the_globals->getECR();
823 >        info->ecr        = globals->getECR();
824        }
825        
826 <      if( !the_globals->haveEST() ){
826 >      if( !globals->haveEST() ){
827          sprintf( painCave.errMsg,
828                   "SimSetup Warning: using default value of 5%% of the "
829                   "electrostaticCutoffRadius for the "
# Line 1052 | Line 833 | void SimSetup::finalInfoCheck( void ){
833          simError();
834          info->est = 0.05 * info->ecr;
835        } else {
836 <        info->est        = the_globals->getEST();
836 >        info->est        = globals->getEST();
837        }
838      }
839    }  
# Line 1066 | Line 847 | void SimSetup::initSystemCoords( void ){
847  
848   void SimSetup::initSystemCoords( void ){
849  
850 < if( the_globals->haveInitialConfig() ){
850 > if( globals->haveInitialConfig() ){
851  
852       InitializeFromFile* fileInit;
853   #ifdef IS_MPI // is_mpi
854       if( worldRank == 0 ){
855   #endif //is_mpi
856 <   fileInit = new InitializeFromFile( the_globals->getInitialConfig() );
856 >   fileInit = new InitializeFromFile( globals->getInitialConfig() );
857   #ifdef IS_MPI
858       }else fileInit = new InitializeFromFile( NULL );
859   #endif
# Line 1113 | Line 894 | void SimSetup::makeOutNames( void ){
894    if( worldRank == 0 ){
895   #endif // is_mpi
896      
897 <    if( the_globals->haveFinalConfig() ){
898 <      strcpy( info->finalName, the_globals->getFinalConfig() );
897 >    if( globals->haveFinalConfig() ){
898 >      strcpy( info->finalName, globals->getFinalConfig() );
899      }
900      else{
901        strcpy( info->finalName, inFileName );
# Line 1197 | Line 978 | void SimSetup::sysObjectsCreation( void ){
978  
979   void SimSetup::sysObjectsCreation( void ){
980  
981 +  int i;
982 +
983    // create the forceField
984  
985    createFF();
# Line 1219 | Line 1002 | void SimSetup::sysObjectsCreation( void ){
1002    
1003    makeSysArrays();
1004  
1005 +  // make and initialize the molecules (all but atomic coordinates)
1006    
1007 +  makeMolecules();
1008 +  info->identArray = new int[info->n_atoms];
1009 +  for(i=0; i<info->n_atoms; i++){
1010 +    info->identArray[i] = the_atoms[i]->getIdent();
1011 +  }
1012 +  
1013  
1014  
1015   }
# Line 1254 | Line 1044 | void SimSetup::compList( void ){
1044  
1045   void SimSetup::compList( void ){
1046  
1047 +  int i;
1048 +
1049    comp_stamps = new MoleculeStamp*[n_components];
1050  
1051    // make an array of molecule stamps that match the components used.
# Line 1279 | Line 1071 | void SimSetup::compList( void ){
1071        
1072        // extract the component from the list;
1073        
1074 <      currentStamp = the_stamps->extractMolStamp( id );
1074 >      currentStamp = stamps->extractMolStamp( id );
1075        if( currentStamp == NULL ){
1076          sprintf( painCave.errMsg,
1077                   "SimSetup error: Component \"%s\" was not found in the "
# Line 1303 | Line 1095 | void SimSetup::calcSysValues( void ){
1095   }
1096  
1097   void SimSetup::calcSysValues( void ){
1098 +  int i, j, k;
1099  
1100 +
1101    tot_atoms = 0;
1102    tot_bonds = 0;
1103    tot_bends = 0;
# Line 1333 | Line 1127 | void SimSetup::mpiMolDivide( void ){
1127  
1128   void SimSetup::mpiMolDivide( void ){
1129    
1130 +  int i, j, k;
1131    int localMol, allMol;
1132    int local_atoms, local_bonds, local_bends, local_torsions, local_SRI;
1133  
# Line 1402 | Line 1197 | void SimSetup::makeSysArrays( void ){
1197  
1198  
1199   void SimSetup::makeSysArrays( void ){
1200 +  int i, j, k;
1201  
1202 +
1203    // create the atom and short range interaction arrays
1204  
1205    Atom::createArrays(info->n_atoms);
# Line 1474 | Line 1271 | void SimSetup::makeSysArrays( void ){
1271    info->excludes = the_excludes;
1272  
1273    the_ff->setSimInfo( info );
1274 +
1275 + }
1276 +
1277 + void SimSetup::makeIntegrator( void ){
1278 +
1279 +  NVT*  myNVT = NULL;
1280 +  NPTi* myNPTi = NULL;
1281 +  NPTf* myNPTf = NULL;
1282 +  NPTim* myNPTim = NULL;
1283 +  NPTfm* myNPTfm = NULL;
1284 +
1285 +  switch( ensembleCase ){
1286 +
1287 +  case NVE_ENS:
1288 +    new NVE( info, the_ff );
1289 +    break;
1290 +
1291 +  case NVT_ENS:
1292 +    myNVT = new NVT( info, the_ff );
1293 +    myNVT->setTargetTemp(globals->getTargetTemp());
1294 +
1295 +    if (globals->haveTauThermostat())
1296 +      myNVT->setTauThermostat(globals->getTauThermostat());
1297 +
1298 +    else {
1299 +      sprintf( painCave.errMsg,
1300 +               "SimSetup error: If you use the NVT\n"
1301 +               "    ensemble, you must set tauThermostat.\n");
1302 +      painCave.isFatal = 1;
1303 +      simError();
1304 +    }
1305 +    break;
1306 +
1307 +  case NPTi_ENS:
1308 +    myNPTi = new NPTi( info, the_ff );
1309 +    myNPTi->setTargetTemp( globals->getTargetTemp() );
1310 +
1311 +    if (globals->haveTargetPressure())
1312 +      myNPTi->setTargetPressure(globals->getTargetPressure());
1313 +    else {
1314 +      sprintf( painCave.errMsg,
1315 +               "SimSetup error: If you use a constant pressure\n"
1316 +               "    ensemble, you must set targetPressure in the BASS file.\n");
1317 +      painCave.isFatal = 1;
1318 +      simError();
1319 +    }
1320 +    
1321 +    if( globals->haveTauThermostat() )
1322 +      myNPTi->setTauThermostat( globals->getTauThermostat() );
1323 +    else{
1324 +      sprintf( painCave.errMsg,
1325 +               "SimSetup error: If you use an NPT\n"
1326 +               "    ensemble, you must set tauThermostat.\n");
1327 +      painCave.isFatal = 1;
1328 +      simError();
1329 +    }
1330  
1331 +    if( globals->haveTauBarostat() )
1332 +      myNPTi->setTauBarostat( globals->getTauBarostat() );
1333 +    else{
1334 +      sprintf( painCave.errMsg,
1335 +               "SimSetup error: If you use an NPT\n"
1336 +               "    ensemble, you must set tauBarostat.\n");
1337 +      painCave.isFatal = 1;
1338 +      simError();
1339 +    }
1340 +    break;
1341 +
1342 +  case NPTf_ENS:
1343 +    myNPTf = new NPTf( info, the_ff );
1344 +    myNPTf->setTargetTemp( globals->getTargetTemp());
1345 +
1346 +    if (globals->haveTargetPressure())
1347 +      myNPTf->setTargetPressure(globals->getTargetPressure());
1348 +    else {
1349 +      sprintf( painCave.errMsg,
1350 +               "SimSetup error: If you use a constant pressure\n"
1351 +               "    ensemble, you must set targetPressure in the BASS file.\n");
1352 +      painCave.isFatal = 1;
1353 +      simError();
1354 +    }    
1355 +
1356 +    if( globals->haveTauThermostat() )
1357 +      myNPTf->setTauThermostat( globals->getTauThermostat() );
1358 +    else{
1359 +      sprintf( painCave.errMsg,
1360 +               "SimSetup error: If you use an NPT\n"
1361 +               "    ensemble, you must set tauThermostat.\n");
1362 +      painCave.isFatal = 1;
1363 +      simError();
1364 +    }
1365 +
1366 +    if( globals->haveTauBarostat() )
1367 +      myNPTf->setTauBarostat( globals->getTauBarostat() );
1368 +    else{
1369 +      sprintf( painCave.errMsg,
1370 +               "SimSetup error: If you use an NPT\n"
1371 +               "    ensemble, you must set tauBarostat.\n");
1372 +      painCave.isFatal = 1;
1373 +      simError();
1374 +    }
1375 +    break;
1376 +    
1377 +  case NPTim_ENS:
1378 +    myNPTim = new NPTim( info, the_ff );
1379 +    myNPTim->setTargetTemp( globals->getTargetTemp());
1380 +
1381 +    if (globals->haveTargetPressure())
1382 +      myNPTim->setTargetPressure(globals->getTargetPressure());
1383 +    else {
1384 +      sprintf( painCave.errMsg,
1385 +               "SimSetup error: If you use a constant pressure\n"
1386 +               "    ensemble, you must set targetPressure in the BASS file.\n");
1387 +      painCave.isFatal = 1;
1388 +      simError();
1389 +    }
1390 +    
1391 +    if( globals->haveTauThermostat() )
1392 +      myNPTim->setTauThermostat( globals->getTauThermostat() );
1393 +    else{
1394 +      sprintf( painCave.errMsg,
1395 +               "SimSetup error: If you use an NPT\n"
1396 +               "    ensemble, you must set tauThermostat.\n");
1397 +      painCave.isFatal = 1;
1398 +      simError();
1399 +    }
1400 +
1401 +    if( globals->haveTauBarostat() )
1402 +      myNPTim->setTauBarostat( globals->getTauBarostat() );
1403 +    else{
1404 +      sprintf( painCave.errMsg,
1405 +               "SimSetup error: If you use an NPT\n"
1406 +               "    ensemble, you must set tauBarostat.\n");
1407 +      painCave.isFatal = 1;
1408 +      simError();
1409 +    }
1410 +    break;
1411 +
1412 +  case NPTfm_ENS:
1413 +    myNPTfm = new NPTfm( info, the_ff );
1414 +    myNPTfm->setTargetTemp( globals->getTargetTemp());
1415 +
1416 +    if (globals->haveTargetPressure())
1417 +      myNPTfm->setTargetPressure(globals->getTargetPressure());
1418 +    else {
1419 +      sprintf( painCave.errMsg,
1420 +               "SimSetup error: If you use a constant pressure\n"
1421 +               "    ensemble, you must set targetPressure in the BASS file.\n");
1422 +      painCave.isFatal = 1;
1423 +      simError();
1424 +    }
1425 +    
1426 +    if( globals->haveTauThermostat() )
1427 +      myNPTfm->setTauThermostat( globals->getTauThermostat() );
1428 +    else{
1429 +      sprintf( painCave.errMsg,
1430 +               "SimSetup error: If you use an NPT\n"
1431 +               "    ensemble, you must set tauThermostat.\n");
1432 +      painCave.isFatal = 1;
1433 +      simError();
1434 +    }
1435 +
1436 +    if( globals->haveTauBarostat() )
1437 +      myNPTfm->setTauBarostat( globals->getTauBarostat() );
1438 +    else{
1439 +      sprintf( painCave.errMsg,
1440 +               "SimSetup error: If you use an NPT\n"
1441 +               "    ensemble, you must set tauBarostat.\n");
1442 +      painCave.isFatal = 1;
1443 +      simError();
1444 +    }
1445 +    break;
1446 +
1447 +  default:
1448 +    sprintf( painCave.errMsg,
1449 +             "SimSetup Error. Unrecognized ensemble in case statement.\n");
1450 +    painCave.isFatal = 1;
1451 +    simError();
1452 +  }
1453 +
1454   }
1455 +
1456 + void SimSetup::initFortran( void ){
1457 +
1458 +  info->refreshSim();
1459 +  
1460 +  if( !strcmp( info->mixingRule, "standard") ){
1461 +    the_ff->initForceField( LB_MIXING_RULE );
1462 +  }
1463 +  else if( !strcmp( info->mixingRule, "explicit") ){
1464 +    the_ff->initForceField( EXPLICIT_MIXING_RULE );
1465 +  }
1466 +  else{
1467 +    sprintf( painCave.errMsg,
1468 +             "SimSetup Error: unknown mixing rule -> \"%s\"\n",
1469 +             info->mixingRule );
1470 +    painCave.isFatal = 1;
1471 +    simError();
1472 +  }
1473 +
1474 +
1475 + #ifdef IS_MPI
1476 +  strcpy( checkPointMsg,
1477 +          "Successfully intialized the mixingRule for Fortran." );
1478 +  MPIcheckPoint();
1479 + #endif // is_mpi
1480 +
1481 + }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines