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

Comparing trunk/OOPSE-1.0/libmdtools/SimSetup.cpp (file contents):
Revision 1415 by gezelter, Mon Jul 26 17:50:57 2004 UTC vs.
Revision 1451 by chrisfen, Mon Aug 9 14:50:35 2004 UTC

# Line 57 | Line 57 | bool isDivisible(double dividend, double divisor){
57      return false;  
58   }
59  
60 + string getPrefix(const string& str ){
61 +  string prefix;
62 +  string suffix;
63 +  int pos;
64 +
65 +  pos = str.rfind(".");
66 +
67 +  if (pos >= 0) {
68 +     prefix = str.substr(0, pos);
69 +     suffix = str.substr(pos, str.size());
70 +
71 +     // leave .bass there in case we've reverted to old habits
72 +     if (LowerCase(suffix) == ".md" || LowerCase(suffix) == ".bass")
73 +       return prefix;
74 +     else
75 +       return str;
76 +    
77 +  } else
78 +    return str;
79 + };
80 +
81 +
82   SimSetup::SimSetup(){
83    
84    initSuspend = false;
# Line 74 | Line 96 | SimSetup::~SimSetup(){
96   }
97  
98   SimSetup::~SimSetup(){
99 +  // clean up the forcefield
100 +  the_ff->cleanMe();
101 +
102    delete stamps;
103    delete globals;
104   }
# Line 92 | Line 117 | void SimSetup::parseFile(char* fileName){
117   #endif // is_mpi
118  
119      inFileName = fileName;
120 +
121 +    globals->initalize();
122      set_interface_stamps(stamps, globals);
123  
124   #ifdef IS_MPI
# Line 122 | Line 149 | void SimSetup::createSim(void){
149  
150   void SimSetup::createSim(void){
151  
152 <  // gather all of the information from the Bass file
152 >  // gather all of the information from the meta-data file
153  
154    gatherInfo();
155  
# Line 643 | Line 670 | void SimSetup::makeMolecules(void){
670    sprintf(checkPointMsg, "all molecules initialized succesfully");
671    MPIcheckPoint();
672   #endif // is_mpi
646
647 }
648
649 void SimSetup::initFromBass(void){
650  int i, j, k;
651  int n_cells;
652  double cellx, celly, cellz;
653  double temp1, temp2, temp3;
654  int n_per_extra;
655  int n_extra;
656  int have_extra, done;
657
658  double vel[3];
659  vel[0] = 0.0;
660  vel[1] = 0.0;
661  vel[2] = 0.0;
662
663  temp1 = (double) tot_nmol / 4.0;
664  temp2 = pow(temp1, (1.0 / 3.0));
665  temp3 = ceil(temp2);
666
667  have_extra = 0;
668  if (temp2 < temp3){
669    // we have a non-complete lattice
670    have_extra = 1;
671
672    n_cells = (int) temp3 - 1;
673    cellx = info[0].boxL[0] / temp3;
674    celly = info[0].boxL[1] / temp3;
675    cellz = info[0].boxL[2] / temp3;
676    n_extra = tot_nmol - (4 * n_cells * n_cells * n_cells);
677    temp1 = ((double) n_extra) / (pow(temp3, 3.0) - pow(n_cells, 3.0));
678    n_per_extra = (int) ceil(temp1);
679
680    if (n_per_extra > 4){
681      sprintf(painCave.errMsg,
682              "SimSetup error. There has been an error in constructing"
683              " the non-complete lattice.\n");
684      painCave.isFatal = 1;
685      simError();
686    }
687  }
688  else{
689    n_cells = (int) temp3;
690    cellx = info[0].boxL[0] / temp3;
691    celly = info[0].boxL[1] / temp3;
692    cellz = info[0].boxL[2] / temp3;
693  }
694
695  current_mol = 0;
696  current_comp_mol = 0;
697  current_comp = 0;
698  current_atom_ndx = 0;
699
700  for (i = 0; i < n_cells ; i++){
701    for (j = 0; j < n_cells; j++){
702      for (k = 0; k < n_cells; k++){
703        makeElement(i * cellx, j * celly, k * cellz);
673  
705        makeElement(i * cellx + 0.5 * cellx, j * celly + 0.5 * celly, k * cellz);
706
707        makeElement(i * cellx, j * celly + 0.5 * celly, k * cellz + 0.5 * cellz);
708
709        makeElement(i * cellx + 0.5 * cellx, j * celly, k * cellz + 0.5 * cellz);
710      }
711    }
712  }
713
714  if (have_extra){
715    done = 0;
716
717    int start_ndx;
718    for (i = 0; i < (n_cells + 1) && !done; i++){
719      for (j = 0; j < (n_cells + 1) && !done; j++){
720        if (i < n_cells){
721          if (j < n_cells){
722            start_ndx = n_cells;
723          }
724          else
725            start_ndx = 0;
726        }
727        else
728          start_ndx = 0;
729
730        for (k = start_ndx; k < (n_cells + 1) && !done; k++){
731          makeElement(i * cellx, j * celly, k * cellz);
732          done = (current_mol >= tot_nmol);
733
734          if (!done && n_per_extra > 1){
735            makeElement(i * cellx + 0.5 * cellx, j * celly + 0.5 * celly,
736                        k * cellz);
737            done = (current_mol >= tot_nmol);
738          }
739
740          if (!done && n_per_extra > 2){
741            makeElement(i * cellx, j * celly + 0.5 * celly,
742                        k * cellz + 0.5 * cellz);
743            done = (current_mol >= tot_nmol);
744          }
745
746          if (!done && n_per_extra > 3){
747            makeElement(i * cellx + 0.5 * cellx, j * celly,
748                        k * cellz + 0.5 * cellz);
749            done = (current_mol >= tot_nmol);
750          }
751        }
752      }
753    }
754  }
755
756  for (i = 0; i < info[0].n_atoms; i++){
757    info[0].atoms[i]->setVel(vel);
758  }
674   }
675  
761 void SimSetup::makeElement(double x, double y, double z){
762  int k;
763  AtomStamp* current_atom;
764  DirectionalAtom* dAtom;
765  double rotMat[3][3];
766  double pos[3];
767
768  for (k = 0; k < comp_stamps[current_comp]->getNAtoms(); k++){
769    current_atom = comp_stamps[current_comp]->getAtom(k);
770    if (!current_atom->havePosition()){
771      sprintf(painCave.errMsg,
772              "SimSetup:initFromBass error.\n"
773              "\tComponent %s, atom %s does not have a position specified.\n"
774              "\tThe initialization routine is unable to give a start"
775              " position.\n",
776              comp_stamps[current_comp]->getID(), current_atom->getType());
777      painCave.isFatal = 1;
778      simError();
779    }
780
781    pos[0] = x + current_atom->getPosX();
782    pos[1] = y + current_atom->getPosY();
783    pos[2] = z + current_atom->getPosZ();
784
785    info[0].atoms[current_atom_ndx]->setPos(pos);
786
787    if (info[0].atoms[current_atom_ndx]->isDirectional()){
788      dAtom = (DirectionalAtom *) info[0].atoms[current_atom_ndx];
789
790      rotMat[0][0] = 1.0;
791      rotMat[0][1] = 0.0;
792      rotMat[0][2] = 0.0;
793
794      rotMat[1][0] = 0.0;
795      rotMat[1][1] = 1.0;
796      rotMat[1][2] = 0.0;
797
798      rotMat[2][0] = 0.0;
799      rotMat[2][1] = 0.0;
800      rotMat[2][2] = 1.0;
801
802      dAtom->setA(rotMat);
803    }
804
805    current_atom_ndx++;
806  }
807
808  current_mol++;
809  current_comp_mol++;
810
811  if (current_comp_mol >= components_nmol[current_comp]){
812    current_comp_mol = 0;
813    current_comp++;
814  }
815 }
816
817
676   void SimSetup::gatherInfo(void){
677    int i;
678  
# Line 974 | Line 832 | void SimSetup::gatherInfo(void){
832        }        
833      }
834  
977    
835      for (i=0; i < nInfo; i++) {
836        
837        // check for the temperature set flag
# Line 994 | Line 851 | void SimSetup::gatherInfo(void){
851            info[i].thermIntLambda = globals->getThermIntLambda();
852            info[i].thermIntK = globals->getThermIntK();
853            
854 <          Restraints *myRestraint = new Restraints(tot_nmol, info[i].thermIntLambda, info[i].thermIntK);
854 >          Restraints *myRestraint = new Restraints(info[i].thermIntLambda, info[i].thermIntK);
855            info[i].restraint = myRestraint;
856          }
857          else {
# Line 1003 | Line 860 | void SimSetup::gatherInfo(void){
860                    "\tKeyword useSolidThermInt was set to 'true' but\n"
861                    "\tthermodynamicIntegrationLambda (and/or\n"
862                    "\tthermodynamicIntegrationK) was not specified.\n"
863 <                  "\tPlease provide a lambda value and k value in your .bass file.\n");
863 >                  "\tPlease provide a lambda value and k value in your meta-data file.\n");
864            painCave.isFatal = 1;
865            simError();    
866          }
# Line 1012 | Line 869 | void SimSetup::gatherInfo(void){
869          if (globals->getUseSolidThermInt()) {
870            sprintf( painCave.errMsg,
871                     "SimSetup Warning: It appears that you have both solid and\n"
872 <                   "\tliquid thermodynamic integration activated in your .bass\n"
872 >                   "\tliquid thermodynamic integration activated in your meta-data\n"
873                     "\tfile. To avoid confusion, specify only one technique in\n"
874 <                   "\tyour .bass file. Liquid-state thermodynamic integration\n"
874 >                   "\tyour meta-data file. Liquid-state thermodynamic integration\n"
875                     "\twill be assumed for the current simulation. If this is not\n"
876                     "\twhat you desire, set useSolidThermInt to 'true' and\n"
877 <                   "\tuseLiquidThermInt to 'false' in your .bass file.\n");
877 >                   "\tuseLiquidThermInt to 'false' in your meta-data file.\n");
878            painCave.isFatal = 0;
879            simError();
880          }
# Line 1032 | Line 889 | void SimSetup::gatherInfo(void){
889                    "\tKeyword useLiquidThermInt was set to 'true' but\n"
890                    "\tthermodynamicIntegrationLambda (and/or\n"
891                    "\tthermodynamicIntegrationK) was not specified.\n"
892 <                  "\tPlease provide a lambda value and k value in your .bass file.\n");
892 >                  "\tPlease provide a lambda value and k value in your meta-data file.\n");
893            painCave.isFatal = 1;
894            simError();    
895          }
# Line 1041 | Line 898 | void SimSetup::gatherInfo(void){
898          sprintf(painCave.errMsg,
899                  "SimSetup Warning: If you want to use Thermodynamic\n"
900                  "\tIntegration, set useSolidThermInt or useLiquidThermInt to\n"
901 <                "\t'true' in your .bass file.  These keywords are set to\n"
901 >                "\t'true' in your meta-data file.  These keywords are set to\n"
902                  "\t'false' by default, so your lambda and/or k values are\n"
903                  "\tbeing ignored.\n");
904          painCave.isFatal = 0;
# Line 1136 | Line 993 | void SimSetup::gatherInfo(void){
993    }
994    
995   #ifdef IS_MPI
996 <  strcpy(checkPointMsg, "Successfully gathered all information from Bass\n");
996 >  strcpy(checkPointMsg, "Successfully gathered all information from meta-data file\n");
997    MPIcheckPoint();
998   #endif // is_mpi
999   }
# Line 1273 | Line 1130 | void SimSetup::finalInfoCheck(void){
1130    MPIcheckPoint();
1131   #endif // is_mpi
1132  
1276  // clean up the forcefield
1277  the_ff->cleanMe();
1133   }
1134    
1135   void SimSetup::initSystemCoords(void){
# Line 1305 | Line 1160 | void SimSetup::initSystemCoords(void){
1160    }
1161    else{
1162      
1163 <    // no init from bass
1163 >    // no init from md file
1164      
1165      sprintf(painCave.errMsg,
1166              "Cannot intialize a simulation without an initial configuration file.\n");
# Line 1323 | Line 1178 | void SimSetup::makeOutNames(void){
1178  
1179   void SimSetup::makeOutNames(void){
1180    int k;
1181 +  string prefix;
1182  
1327
1183    for (k = 0; k < nInfo; k++){
1184   #ifdef IS_MPI
1185      if (worldRank == 0){
1186   #endif // is_mpi
1187 <
1188 <      if (globals->haveFinalConfig()){
1189 <        strcpy(info[k].finalName, globals->getFinalConfig());
1190 <      }
1191 <      else{
1337 <        strcpy(info[k].finalName, inFileName);
1338 <        char* endTest;
1339 <        int nameLength = strlen(info[k].finalName);
1340 <        endTest = &(info[k].finalName[nameLength - 5]);
1341 <        if (!strcmp(endTest, ".bass")){
1342 <          strcpy(endTest, ".eor");
1343 <        }
1344 <        else if (!strcmp(endTest, ".BASS")){
1345 <          strcpy(endTest, ".eor");
1346 <        }
1347 <        else{
1348 <          endTest = &(info[k].finalName[nameLength - 4]);
1349 <          if (!strcmp(endTest, ".bss")){
1350 <            strcpy(endTest, ".eor");
1351 <          }
1352 <          else if (!strcmp(endTest, ".mdl")){
1353 <            strcpy(endTest, ".eor");
1354 <          }
1355 <          else{
1356 <            strcat(info[k].finalName, ".eor");
1357 <          }
1358 <        }
1359 <      }
1187 >      
1188 >      if(globals->haveFinalConfig())
1189 >        prefix = getPrefix(globals->getFinalConfig());  
1190 >      else
1191 >        prefix = getPrefix(inFileName);
1192  
1193 <      // make the sample and status out names
1194 <
1195 <      strcpy(info[k].sampleName, inFileName);
1364 <      char* endTest;
1365 <      int nameLength = strlen(info[k].sampleName);
1366 <      endTest = &(info[k].sampleName[nameLength - 5]);
1367 <      if (!strcmp(endTest, ".bass")){
1368 <        strcpy(endTest, ".dump");
1369 <      }
1370 <      else if (!strcmp(endTest, ".BASS")){
1371 <        strcpy(endTest, ".dump");
1372 <      }
1373 <      else{
1374 <        endTest = &(info[k].sampleName[nameLength - 4]);
1375 <        if (!strcmp(endTest, ".bss")){
1376 <          strcpy(endTest, ".dump");
1377 <        }
1378 <        else if (!strcmp(endTest, ".mdl")){
1379 <          strcpy(endTest, ".dump");
1380 <        }
1381 <        else{
1382 <          strcat(info[k].sampleName, ".dump");
1383 <        }
1384 <      }
1385 <
1386 <      strcpy(info[k].statusName, inFileName);
1387 <      nameLength = strlen(info[k].statusName);
1388 <      endTest = &(info[k].statusName[nameLength - 5]);
1389 <      if (!strcmp(endTest, ".bass")){
1390 <        strcpy(endTest, ".stat");
1391 <      }
1392 <      else if (!strcmp(endTest, ".BASS")){
1393 <        strcpy(endTest, ".stat");
1394 <      }
1395 <      else{
1396 <        endTest = &(info[k].statusName[nameLength - 4]);
1397 <        if (!strcmp(endTest, ".bss")){
1398 <          strcpy(endTest, ".stat");
1399 <        }
1400 <        else if (!strcmp(endTest, ".mdl")){
1401 <          strcpy(endTest, ".stat");
1402 <        }
1403 <        else{
1404 <          strcat(info[k].statusName, ".stat");
1405 <        }
1406 <      }
1407 <
1408 <      strcpy(info[k].rawPotName, inFileName);
1409 <      nameLength = strlen(info[k].rawPotName);
1410 <      endTest = &(info[k].rawPotName[nameLength - 5]);
1411 <      if (!strcmp(endTest, ".bass")){
1412 <        strcpy(endTest, ".raw");
1413 <      }
1414 <      else if (!strcmp(endTest, ".BASS")){
1415 <        strcpy(endTest, ".raw");
1416 <      }
1417 <      else{
1418 <        endTest = &(info[k].rawPotName[nameLength - 4]);
1419 <        if (!strcmp(endTest, ".bss")){
1420 <          strcpy(endTest, ".raw");
1421 <        }
1422 <        else if (!strcmp(endTest, ".mdl")){
1423 <          strcpy(endTest, ".raw");
1424 <        }
1425 <        else{
1426 <          strcat(info[k].rawPotName, ".raw");
1427 <        }
1428 <      }
1193 >      info[k].finalName = prefix + ".eor";      
1194 >      info[k].sampleName = prefix + ".dump";
1195 >      info[k].statusName = prefix + ".stat";
1196  
1197   #ifdef IS_MPI
1198  
# Line 1853 | Line 1620 | void SimSetup::makeIntegrator(void){
1620          else{
1621            sprintf(painCave.errMsg,
1622                    "SimSetup error: If you use a constant pressure\n"
1623 <                  "\tensemble, you must set targetPressure in the BASS file.\n");
1623 >                  "\tensemble, you must set targetPressure in the meta-data file.\n");
1624            painCave.isFatal = 1;
1625            simError();
1626          }
# Line 1904 | Line 1671 | void SimSetup::makeIntegrator(void){
1671          else{
1672            sprintf(painCave.errMsg,
1673                    "SimSetup error: If you use a constant pressure\n"
1674 <                  "\tensemble, you must set targetPressure in the BASS file.\n");
1674 >                  "\tensemble, you must set targetPressure in the meta-data file.\n");
1675            painCave.isFatal = 1;
1676            simError();
1677          }    
# Line 1957 | Line 1724 | void SimSetup::makeIntegrator(void){
1724          else{
1725            sprintf(painCave.errMsg,
1726                    "SimSetup error: If you use a constant pressure\n"
1727 <                  "\tensemble, you must set targetPressure in the BASS file.\n");
1727 >                  "\tensemble, you must set targetPressure in the meta-data file.\n");
1728            painCave.isFatal = 1;
1729            simError();
1730          }    

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines