# | Line 1 | Line 1 | |
---|---|---|
1 | #include <algorithm> | |
2 | < | #include <cstdlib> |
2 | > | #include <stdlib.h> |
3 | #include <iostream> | |
4 | < | #include <cmath> |
4 | > | #include <math.h> |
5 | #include <string> | |
6 | #include <sprng.h> | |
7 | – | |
7 | #include "SimSetup.hpp" | |
8 | #include "ReadWrite.hpp" | |
9 | #include "parse_me.h" | |
# | Line 22 | Line 21 | |
21 | #define NVT_ENS 1 | |
22 | #define NPTi_ENS 2 | |
23 | #define NPTf_ENS 3 | |
24 | + | #define NPTxyz_ENS 4 |
25 | ||
26 | + | |
27 | #define FF_DUFF 0 | |
28 | #define FF_LJ 1 | |
29 | #define FF_EAM 2 | |
# | Line 30 | Line 31 | SimSetup::SimSetup(){ | |
31 | using namespace std; | |
32 | ||
33 | SimSetup::SimSetup(){ | |
34 | + | |
35 | + | initSuspend = false; |
36 | isInfoArray = 0; | |
37 | nInfo = 1; | |
38 | ||
# | Line 52 | Line 55 | void SimSetup::setSimInfo(SimInfo* the_info, int theNi | |
55 | info = the_info; | |
56 | nInfo = theNinfo; | |
57 | isInfoArray = 1; | |
58 | + | initSuspend = true; |
59 | } | |
60 | ||
61 | ||
# | Line 90 | Line 94 | void SimSetup::createSim(void){ | |
94 | #endif // is_mpi | |
95 | ||
96 | void SimSetup::createSim(void){ | |
93 | – | int i, j, k, globalAtomIndex; |
97 | ||
98 | // gather all of the information from the Bass file | |
99 | ||
# | Line 106 | Line 109 | void SimSetup::createSim(void){ | |
109 | ||
110 | // initialize the system coordinates | |
111 | ||
112 | < | if (!isInfoArray){ |
112 | > | if ( !initSuspend ){ |
113 | initSystemCoords(); | |
114 | + | |
115 | + | if( !(globals->getUseInitTime()) ) |
116 | + | info[0].currentTime = 0.0; |
117 | } | |
118 | ||
119 | // make the output filenames | |
# | Line 129 | Line 135 | void SimSetup::makeMolecules(void){ | |
135 | ||
136 | ||
137 | void SimSetup::makeMolecules(void){ | |
138 | < | int k, l; |
138 | > | int k; |
139 | int i, j, exI, exJ, tempEx, stampID, atomOffset, excludeOffset; | |
140 | molInit molInfo; | |
141 | DirectionalAtom* dAtom; | |
# | Line 551 | Line 557 | void SimSetup::gatherInfo(void){ | |
557 | ||
558 | ||
559 | void SimSetup::gatherInfo(void){ | |
560 | < | int i, j, k; |
560 | > | int i; |
561 | ||
562 | ensembleCase = -1; | |
563 | ffCase = -1; | |
# | Line 602 | Line 608 | void SimSetup::gatherInfo(void){ | |
608 | else if (!strcasecmp(ensemble, "NPTf")){ | |
609 | ensembleCase = NPTf_ENS; | |
610 | } | |
611 | + | else if (!strcasecmp(ensemble, "NPTxyz")){ |
612 | + | ensembleCase = NPTxyz_ENS; |
613 | + | } |
614 | else{ | |
615 | sprintf(painCave.errMsg, | |
616 | "SimSetup Warning. Unrecognized Ensemble -> %s, " | |
# | Line 687 | Line 696 | void SimSetup::gatherInfo(void){ | |
696 | } | |
697 | ||
698 | // check for the temperature set flag | |
699 | < | |
699 | > | |
700 | if (globals->haveTempSet()) | |
701 | info[i].setTemp = globals->getTempSet(); | |
702 | ||
703 | < | // get some of the tricky things that may still be in the globals |
703 | > | // check for the extended State init |
704 | ||
705 | < | double boxVector[3]; |
706 | < | if (globals->haveBox()){ |
707 | < | boxVector[0] = globals->getBox(); |
699 | < | boxVector[1] = globals->getBox(); |
700 | < | boxVector[2] = globals->getBox(); |
701 | < | |
702 | < | info[i].setBox(boxVector); |
703 | < | } |
704 | < | else if (globals->haveDensity()){ |
705 | < | double vol; |
706 | < | vol = (double) tot_nmol / globals->getDensity(); |
707 | < | boxVector[0] = pow(vol, (1.0 / 3.0)); |
708 | < | boxVector[1] = boxVector[0]; |
709 | < | boxVector[2] = boxVector[0]; |
710 | < | |
711 | < | info[i].setBox(boxVector); |
712 | < | } |
713 | < | else{ |
714 | < | if (!globals->haveBoxX()){ |
715 | < | sprintf(painCave.errMsg, |
716 | < | "SimSetup error, no periodic BoxX size given.\n"); |
717 | < | painCave.isFatal = 1; |
718 | < | simError(); |
719 | < | } |
720 | < | boxVector[0] = globals->getBoxX(); |
721 | < | |
722 | < | if (!globals->haveBoxY()){ |
723 | < | sprintf(painCave.errMsg, |
724 | < | "SimSetup error, no periodic BoxY size given.\n"); |
725 | < | painCave.isFatal = 1; |
726 | < | simError(); |
727 | < | } |
728 | < | boxVector[1] = globals->getBoxY(); |
729 | < | |
730 | < | if (!globals->haveBoxZ()){ |
731 | < | sprintf(painCave.errMsg, |
732 | < | "SimSetup error, no periodic BoxZ size given.\n"); |
733 | < | painCave.isFatal = 1; |
734 | < | simError(); |
735 | < | } |
736 | < | boxVector[2] = globals->getBoxZ(); |
737 | < | |
738 | < | info[i].setBox(boxVector); |
739 | < | } |
705 | > | info[i].useInitXSstate = globals->getUseInitXSstate(); |
706 | > | info[i].orthoTolerance = globals->getOrthoBoxTolerance(); |
707 | > | |
708 | } | |
709 | < | |
709 | > | |
710 | //setup seed for random number generator | |
711 | int seedValue; | |
712 | ||
# | Line 813 | Line 781 | void SimSetup::finalInfoCheck(void){ | |
781 | ||
782 | if (!globals->haveECR()){ | |
783 | sprintf(painCave.errMsg, | |
784 | < | "SimSetup Warning: using default value of 1/2 the smallest " |
785 | < | "box length for the electrostaticCutoffRadius.\n" |
818 | < | "I hope you have a very fast processor!\n"); |
784 | > | "SimSetup Warning: using default value of 15.0 angstroms" |
785 | > | "box length for the electrostaticCutoffRadius.\n"); |
786 | painCave.isFatal = 0; | |
787 | simError(); | |
788 | < | double smallest; |
822 | < | smallest = info[i].boxL[0]; |
823 | < | if (info[i].boxL[1] <= smallest) |
824 | < | smallest = info[i].boxL[1]; |
825 | < | if (info[i].boxL[2] <= smallest) |
826 | < | smallest = info[i].boxL[2]; |
827 | < | theEcr = 0.5 * smallest; |
788 | > | theEcr = 15.0; |
789 | } | |
790 | else{ | |
791 | theEcr = globals->getECR(); | |
# | Line 842 | Line 803 | void SimSetup::finalInfoCheck(void){ | |
803 | theEst = globals->getEST(); | |
804 | } | |
805 | ||
806 | < | info[i].setEcr(theEcr, theEst); |
806 | > | info[i].setDefaultEcr(theEcr, theEst); |
807 | ||
808 | if (!globals->haveDielectric()){ | |
809 | sprintf(painCave.errMsg, | |
# | Line 857 | Line 818 | void SimSetup::finalInfoCheck(void){ | |
818 | if (usesDipoles){ | |
819 | if (!globals->haveECR()){ | |
820 | sprintf(painCave.errMsg, | |
821 | < | "SimSetup Warning: using default value of 1/2 the smallest " |
822 | < | "box length for the electrostaticCutoffRadius.\n" |
823 | < | "I hope you have a very fast processor!\n"); |
824 | < | painCave.isFatal = 0; |
825 | < | simError(); |
865 | < | double smallest; |
866 | < | smallest = info[i].boxL[0]; |
867 | < | if (info[i].boxL[1] <= smallest) |
868 | < | smallest = info[i].boxL[1]; |
869 | < | if (info[i].boxL[2] <= smallest) |
870 | < | smallest = info[i].boxL[2]; |
871 | < | theEcr = 0.5 * smallest; |
821 | > | "SimSetup Warning: using default value of 15.0 angstroms" |
822 | > | "box length for the electrostaticCutoffRadius.\n"); |
823 | > | painCave.isFatal = 0; |
824 | > | simError(); |
825 | > | theEcr = 15.0; |
826 | } | |
827 | else{ | |
828 | theEcr = globals->getECR(); | |
829 | } | |
830 | < | |
830 | > | |
831 | if (!globals->haveEST()){ | |
832 | sprintf(painCave.errMsg, | |
833 | "SimSetup Warning: using default value of 0.05 * the " | |
# | Line 886 | Line 840 | void SimSetup::finalInfoCheck(void){ | |
840 | else{ | |
841 | theEst = globals->getEST(); | |
842 | } | |
843 | < | |
844 | < | info[i].setEcr(theEcr, theEst); |
843 | > | |
844 | > | info[i].setDefaultEcr(theEcr, theEst); |
845 | } | |
846 | } | |
847 | } | |
894 | – | |
848 | #ifdef IS_MPI | |
849 | strcpy(checkPointMsg, "post processing checks out"); | |
850 | MPIcheckPoint(); | |
851 | #endif // is_mpi | |
852 | } | |
853 | < | |
853 | > | |
854 | void SimSetup::initSystemCoords(void){ | |
855 | int i; | |
856 | ||
# | Line 914 | Line 867 | void SimSetup::initSystemCoords(void){ | |
867 | if (worldRank == 0){ | |
868 | #endif //is_mpi | |
869 | inName = globals->getInitialConfig(); | |
917 | – | double* tempDouble = new double[1000000]; |
870 | fileInit = new InitializeFromFile(inName); | |
871 | #ifdef IS_MPI | |
872 | } | |
# | Line 926 | Line 878 | void SimSetup::initSystemCoords(void){ | |
878 | delete fileInit; | |
879 | } | |
880 | else{ | |
881 | < | #ifdef IS_MPI |
930 | < | |
881 | > | |
882 | // no init from bass | |
883 | < | |
883 | > | |
884 | sprintf(painCave.errMsg, | |
885 | < | "Cannot intialize a parallel simulation without an initial configuration file.\n"); |
886 | < | painCave.isFatal; |
885 | > | "Cannot intialize a simulation without an initial configuration file.\n"); |
886 | > | painCave.isFatal = 1;; |
887 | simError(); | |
888 | < | |
938 | < | #else |
939 | < | |
940 | < | initFromBass(); |
941 | < | |
942 | < | |
943 | < | #endif |
888 | > | |
889 | } | |
890 | ||
891 | #ifdef IS_MPI | |
# | Line 1158 | Line 1103 | void SimSetup::calcSysValues(void){ | |
1103 | } | |
1104 | ||
1105 | void SimSetup::calcSysValues(void){ | |
1106 | < | int i, j, k; |
1106 | > | int i; |
1107 | ||
1108 | int* molMembershipArray; | |
1109 | ||
# | Line 1257 | Line 1202 | void SimSetup::makeSysArrays(void){ | |
1202 | ||
1203 | ||
1204 | void SimSetup::makeSysArrays(void){ | |
1205 | < | int i, j, k, l; |
1205 | > | |
1206 | > | #ifndef IS_MPI |
1207 | > | int k, j; |
1208 | > | #endif // is_mpi |
1209 | > | int i, l; |
1210 | ||
1211 | Atom** the_atoms; | |
1212 | Molecule* the_molecules; | |
# | Line 1344 | Line 1293 | void SimSetup::makeIntegrator(void){ | |
1293 | NVT<RealIntegrator>* myNVT = NULL; | |
1294 | NPTi<NPT<RealIntegrator> >* myNPTi = NULL; | |
1295 | NPTf<NPT<RealIntegrator> >* myNPTf = NULL; | |
1296 | + | NPTxyz<NPT<RealIntegrator> >* myNPTxyz = NULL; |
1297 | ||
1298 | for (k = 0; k < nInfo; k++){ | |
1299 | switch (ensembleCase){ | |
# | Line 1447 | Line 1397 | void SimSetup::makeIntegrator(void){ | |
1397 | ||
1398 | if (globals->haveTauThermostat()) | |
1399 | myNPTf->setTauThermostat(globals->getTauThermostat()); | |
1400 | + | |
1401 | else{ | |
1402 | sprintf(painCave.errMsg, | |
1403 | "SimSetup error: If you use an NPT\n" | |
# | Line 1457 | Line 1408 | void SimSetup::makeIntegrator(void){ | |
1408 | ||
1409 | if (globals->haveTauBarostat()) | |
1410 | myNPTf->setTauBarostat(globals->getTauBarostat()); | |
1411 | + | |
1412 | else{ | |
1413 | sprintf(painCave.errMsg, | |
1414 | "SimSetup error: If you use an NPT\n" | |
# | Line 1468 | Line 1420 | void SimSetup::makeIntegrator(void){ | |
1420 | info->the_integrator = myNPTf; | |
1421 | break; | |
1422 | ||
1423 | + | case NPTxyz_ENS: |
1424 | + | if (globals->haveZconstraints()){ |
1425 | + | setupZConstraint(info[k]); |
1426 | + | myNPTxyz = new ZConstraint<NPTxyz<NPT <RealIntegrator> > >(&(info[k]), the_ff); |
1427 | + | } |
1428 | + | else |
1429 | + | myNPTxyz = new NPTxyz<NPT <RealIntegrator> >(&(info[k]), the_ff); |
1430 | + | |
1431 | + | myNPTxyz->setTargetTemp(globals->getTargetTemp()); |
1432 | + | |
1433 | + | if (globals->haveTargetPressure()) |
1434 | + | myNPTxyz->setTargetPressure(globals->getTargetPressure()); |
1435 | + | else{ |
1436 | + | sprintf(painCave.errMsg, |
1437 | + | "SimSetup error: If you use a constant pressure\n" |
1438 | + | " ensemble, you must set targetPressure in the BASS file.\n"); |
1439 | + | painCave.isFatal = 1; |
1440 | + | simError(); |
1441 | + | } |
1442 | + | |
1443 | + | if (globals->haveTauThermostat()) |
1444 | + | myNPTxyz->setTauThermostat(globals->getTauThermostat()); |
1445 | + | else{ |
1446 | + | sprintf(painCave.errMsg, |
1447 | + | "SimSetup error: If you use an NPT\n" |
1448 | + | " ensemble, you must set tauThermostat.\n"); |
1449 | + | painCave.isFatal = 1; |
1450 | + | simError(); |
1451 | + | } |
1452 | + | |
1453 | + | if (globals->haveTauBarostat()) |
1454 | + | myNPTxyz->setTauBarostat(globals->getTauBarostat()); |
1455 | + | else{ |
1456 | + | sprintf(painCave.errMsg, |
1457 | + | "SimSetup error: If you use an NPT\n" |
1458 | + | " ensemble, you must set tauBarostat.\n"); |
1459 | + | painCave.isFatal = 1; |
1460 | + | simError(); |
1461 | + | } |
1462 | + | |
1463 | + | info->the_integrator = myNPTxyz; |
1464 | + | break; |
1465 | + | |
1466 | default: | |
1467 | sprintf(painCave.errMsg, | |
1468 | "SimSetup Error. Unrecognized ensemble in case statement.\n"); |
– | Removed lines |
+ | Added lines |
< | Changed lines |
> | Changed lines |