# | 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 | – | #define FF_DUFF 0 |
27 | – | #define FF_LJ 1 |
28 | – | #define FF_EAM 2 |
26 | ||
27 | + | #define FF_DUFF 0 |
28 | + | #define FF_LJ 1 |
29 | + | #define FF_EAM 2 |
30 | + | #define FF_H2O 3 |
31 | + | |
32 | using namespace std; | |
33 | ||
34 | + | /** |
35 | + | * Check whether dividend is divisble by divisor or not |
36 | + | */ |
37 | + | bool isDivisible(double dividend, double divisor){ |
38 | + | double tolerance = 0.000001; |
39 | + | double quotient; |
40 | + | double diff; |
41 | + | int intQuotient; |
42 | + | |
43 | + | quotient = dividend / divisor; |
44 | + | |
45 | + | if (quotient < 0) |
46 | + | quotient = -quotient; |
47 | + | |
48 | + | intQuotient = int (quotient + tolerance); |
49 | + | |
50 | + | diff = fabs(fabs(dividend) - intQuotient * fabs(divisor)); |
51 | + | |
52 | + | if (diff <= tolerance) |
53 | + | return true; |
54 | + | else |
55 | + | return false; |
56 | + | } |
57 | + | |
58 | SimSetup::SimSetup(){ | |
59 | + | |
60 | + | initSuspend = false; |
61 | isInfoArray = 0; | |
62 | nInfo = 1; | |
63 | ||
# | Line 52 | Line 80 | void SimSetup::setSimInfo(SimInfo* the_info, int theNi | |
80 | info = the_info; | |
81 | nInfo = theNinfo; | |
82 | isInfoArray = 1; | |
83 | + | initSuspend = true; |
84 | } | |
85 | ||
86 | ||
# | Line 90 | Line 119 | void SimSetup::createSim(void){ | |
119 | #endif // is_mpi | |
120 | ||
121 | void SimSetup::createSim(void){ | |
93 | – | int i, j, k, globalAtomIndex; |
122 | ||
123 | // gather all of the information from the Bass file | |
124 | ||
# | Line 106 | Line 134 | void SimSetup::createSim(void){ | |
134 | ||
135 | // initialize the system coordinates | |
136 | ||
137 | < | if (!isInfoArray){ |
137 | > | if ( !initSuspend ){ |
138 | initSystemCoords(); | |
139 | + | |
140 | + | if( !(globals->getUseInitTime()) ) |
141 | + | info[0].currentTime = 0.0; |
142 | } | |
143 | ||
144 | // make the output filenames | |
# | Line 129 | Line 160 | void SimSetup::makeMolecules(void){ | |
160 | ||
161 | ||
162 | void SimSetup::makeMolecules(void){ | |
163 | < | int k, l; |
163 | > | int k; |
164 | int i, j, exI, exJ, tempEx, stampID, atomOffset, excludeOffset; | |
165 | molInit molInfo; | |
166 | DirectionalAtom* dAtom; | |
# | Line 144 | Line 175 | void SimSetup::makeMolecules(void){ | |
175 | bend_set* theBends; | |
176 | torsion_set* theTorsions; | |
177 | ||
147 | – | |
178 | //init the forceField paramters | |
179 | ||
180 | the_ff->readParams(); | |
# | Line 152 | Line 182 | void SimSetup::makeMolecules(void){ | |
182 | ||
183 | // init the atoms | |
184 | ||
185 | + | double phi, theta, psi; |
186 | + | double sux, suy, suz; |
187 | + | double Axx, Axy, Axz, Ayx, Ayy, Ayz, Azx, Azy, Azz; |
188 | double ux, uy, uz, u, uSqr; | |
189 | ||
190 | for (k = 0; k < nInfo; k++){ | |
# | Line 188 | Line 221 | void SimSetup::makeMolecules(void){ | |
221 | info[k].n_oriented++; | |
222 | molInfo.myAtoms[j] = dAtom; | |
223 | ||
224 | < | ux = currentAtom->getOrntX(); |
225 | < | uy = currentAtom->getOrntY(); |
226 | < | uz = currentAtom->getOrntZ(); |
224 | > | // Directional Atoms have standard unit vectors which are oriented |
225 | > | // in space using the three Euler angles. We assume the standard |
226 | > | // unit vector was originally along the z axis below. |
227 | ||
228 | + | phi = currentAtom->getEulerPhi() * M_PI / 180.0; |
229 | + | theta = currentAtom->getEulerTheta() * M_PI / 180.0; |
230 | + | psi = currentAtom->getEulerPsi()* M_PI / 180.0; |
231 | + | |
232 | + | Axx = (cos(phi) * cos(psi)) - (sin(phi) * cos(theta) * sin(psi)); |
233 | + | Axy = (sin(phi) * cos(psi)) + (cos(phi) * cos(theta) * sin(psi)); |
234 | + | Axz = sin(theta) * sin(psi); |
235 | + | |
236 | + | Ayx = -(cos(phi) * sin(psi)) - (sin(phi) * cos(theta) * cos(psi)); |
237 | + | Ayy = -(sin(phi) * sin(psi)) + (cos(phi) * cos(theta) * cos(psi)); |
238 | + | Ayz = sin(theta) * cos(psi); |
239 | + | |
240 | + | Azx = sin(phi) * sin(theta); |
241 | + | Azy = -cos(phi) * sin(theta); |
242 | + | Azz = cos(theta); |
243 | + | |
244 | + | sux = 0.0; |
245 | + | suy = 0.0; |
246 | + | suz = 1.0; |
247 | + | |
248 | + | ux = (Axx * sux) + (Ayx * suy) + (Azx * suz); |
249 | + | uy = (Axy * sux) + (Ayy * suy) + (Azy * suz); |
250 | + | uz = (Axz * sux) + (Ayz * suy) + (Azz * suz); |
251 | + | |
252 | uSqr = (ux * ux) + (uy * uy) + (uz * uz); | |
253 | ||
254 | u = sqrt(uSqr); | |
# | Line 551 | Line 608 | void SimSetup::gatherInfo(void){ | |
608 | ||
609 | ||
610 | void SimSetup::gatherInfo(void){ | |
611 | < | int i, j, k; |
611 | > | int i; |
612 | ||
613 | ensembleCase = -1; | |
614 | ffCase = -1; | |
# | Line 579 | Line 636 | void SimSetup::gatherInfo(void){ | |
636 | else if (!strcasecmp(force_field, "EAM")){ | |
637 | ffCase = FF_EAM; | |
638 | } | |
639 | + | else if (!strcasecmp(force_field, "WATER")){ |
640 | + | ffCase = FF_H2O; |
641 | + | } |
642 | else{ | |
643 | sprintf(painCave.errMsg, "SimSetup Error. Unrecognized force field -> %s\n", | |
644 | force_field); | |
# | Line 602 | Line 662 | void SimSetup::gatherInfo(void){ | |
662 | else if (!strcasecmp(ensemble, "NPTf")){ | |
663 | ensembleCase = NPTf_ENS; | |
664 | } | |
665 | + | else if (!strcasecmp(ensemble, "NPTxyz")){ |
666 | + | ensembleCase = NPTxyz_ENS; |
667 | + | } |
668 | else{ | |
669 | sprintf(painCave.errMsg, | |
670 | < | "SimSetup Warning. Unrecognized Ensemble -> %s, " |
671 | < | "reverting to NVE for this simulation.\n", |
670 | > | "SimSetup Warning. Unrecognized Ensemble -> %s \n" |
671 | > | "\treverting to NVE for this simulation.\n", |
672 | ensemble); | |
673 | painCave.isFatal = 0; | |
674 | simError(); | |
# | Line 637 | Line 700 | void SimSetup::gatherInfo(void){ | |
700 | if (!the_components[i]->haveNMol()){ | |
701 | // we have a problem | |
702 | sprintf(painCave.errMsg, | |
703 | < | "SimSetup Error. No global NMol or component NMol" |
704 | < | " given. Cannot calculate the number of atoms.\n"); |
703 | > | "SimSetup Error. No global NMol or component NMol given.\n" |
704 | > | "\tCannot calculate the number of atoms.\n"); |
705 | painCave.isFatal = 1; | |
706 | simError(); | |
707 | } | |
# | Line 658 | Line 721 | void SimSetup::gatherInfo(void){ | |
721 | simError(); | |
722 | } | |
723 | ||
724 | + | //check whether sample time, status time, thermal time and reset time are divisble by dt |
725 | + | if (!isDivisible(globals->getSampleTime(), globals->getDt())){ |
726 | + | sprintf(painCave.errMsg, |
727 | + | "Sample time is not divisible by dt.\n" |
728 | + | "\tThis will result in samples that are not uniformly\n" |
729 | + | "\tdistributed in time. If this is a problem, change\n" |
730 | + | "\tyour sampleTime variable.\n"); |
731 | + | painCave.isFatal = 0; |
732 | + | simError(); |
733 | + | } |
734 | + | |
735 | + | if (globals->haveStatusTime() && !isDivisible(globals->getSampleTime(), globals->getDt())){ |
736 | + | sprintf(painCave.errMsg, |
737 | + | "Status time is not divisible by dt.\n" |
738 | + | "\tThis will result in status reports that are not uniformly\n" |
739 | + | "\tdistributed in time. If this is a problem, change \n" |
740 | + | "\tyour statusTime variable.\n"); |
741 | + | painCave.isFatal = 0; |
742 | + | simError(); |
743 | + | } |
744 | + | |
745 | + | if (globals->haveThermalTime() && !isDivisible(globals->getThermalTime(), globals->getDt())){ |
746 | + | sprintf(painCave.errMsg, |
747 | + | "Thermal time is not divisible by dt.\n" |
748 | + | "\tThis will result in thermalizations that are not uniformly\n" |
749 | + | "\tdistributed in time. If this is a problem, change \n" |
750 | + | "\tyour thermalTime variable.\n"); |
751 | + | painCave.isFatal = 0; |
752 | + | simError(); |
753 | + | } |
754 | + | |
755 | + | if (globals->haveResetTime() && !isDivisible(globals->getResetTime(), globals->getDt())){ |
756 | + | sprintf(painCave.errMsg, |
757 | + | "Reset time is not divisible by dt.\n" |
758 | + | "\tThis will result in integrator resets that are not uniformly\n" |
759 | + | "\tdistributed in time. If this is a problem, change\n" |
760 | + | "\tyour resetTime variable.\n"); |
761 | + | painCave.isFatal = 0; |
762 | + | simError(); |
763 | + | } |
764 | + | |
765 | // set the status, sample, and thermal kick times | |
766 | ||
767 | for (i = 0; i < nInfo; i++){ | |
# | Line 687 | Line 791 | void SimSetup::gatherInfo(void){ | |
791 | } | |
792 | ||
793 | // check for the temperature set flag | |
794 | < | |
794 | > | |
795 | if (globals->haveTempSet()) | |
796 | info[i].setTemp = globals->getTempSet(); | |
797 | ||
798 | < | // get some of the tricky things that may still be in the globals |
798 | > | // check for the extended State init |
799 | ||
800 | < | double boxVector[3]; |
801 | < | if (globals->haveBox()){ |
802 | < | 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 | < | } |
800 | > | info[i].useInitXSstate = globals->getUseInitXSstate(); |
801 | > | info[i].orthoTolerance = globals->getOrthoBoxTolerance(); |
802 | > | |
803 | } | |
804 | < | |
804 | > | |
805 | //setup seed for random number generator | |
806 | int seedValue; | |
807 | ||
# | Line 780 | Line 843 | void SimSetup::gatherInfo(void){ | |
843 | } | |
844 | ||
845 | #ifdef IS_MPI | |
846 | < | strcpy(checkPointMsg, "Succesfully gathered all information from Bass\n"); |
846 | > | strcpy(checkPointMsg, "Successfully gathered all information from Bass\n"); |
847 | MPIcheckPoint(); | |
848 | #endif // is_mpi | |
849 | } | |
# | Line 813 | Line 876 | void SimSetup::finalInfoCheck(void){ | |
876 | ||
877 | if (!globals->haveECR()){ | |
878 | sprintf(painCave.errMsg, | |
879 | < | "SimSetup Warning: using default value of 1/2 the smallest " |
880 | < | "box length for the electrostaticCutoffRadius.\n" |
881 | < | "I hope you have a very fast processor!\n"); |
879 | > | "SimSetup Warning: No value was set for electrostaticCutoffRadius.\n" |
880 | > | "\tOOPSE will use a default value of 15.0 angstroms" |
881 | > | "\tfor the electrostaticCutoffRadius.\n"); |
882 | painCave.isFatal = 0; | |
883 | simError(); | |
884 | < | 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; |
884 | > | theEcr = 15.0; |
885 | } | |
886 | else{ | |
887 | theEcr = globals->getECR(); | |
# | Line 832 | Line 889 | void SimSetup::finalInfoCheck(void){ | |
889 | ||
890 | if (!globals->haveEST()){ | |
891 | sprintf(painCave.errMsg, | |
892 | < | "SimSetup Warning: using default value of 0.05 * the " |
893 | < | "electrostaticCutoffRadius for the electrostaticSkinThickness\n"); |
892 | > | "SimSetup Warning: No value was set for electrostaticSkinThickness.\n" |
893 | > | "\tOOPSE will use a default value of\n" |
894 | > | "\t0.05 * electrostaticCutoffRadius\n" |
895 | > | "\tfor the electrostaticSkinThickness\n"); |
896 | painCave.isFatal = 0; | |
897 | simError(); | |
898 | theEst = 0.05 * theEcr; | |
# | Line 842 | Line 901 | void SimSetup::finalInfoCheck(void){ | |
901 | theEst = globals->getEST(); | |
902 | } | |
903 | ||
904 | < | info[i].setEcr(theEcr, theEst); |
904 | > | info[i].setDefaultEcr(theEcr, theEst); |
905 | ||
906 | if (!globals->haveDielectric()){ | |
907 | sprintf(painCave.errMsg, | |
908 | < | "SimSetup Error: You are trying to use Reaction Field without" |
909 | < | "setting a dielectric constant!\n"); |
908 | > | "SimSetup Error: No Dielectric constant was set.\n" |
909 | > | "\tYou are trying to use Reaction Field without" |
910 | > | "\tsetting a dielectric constant!\n"); |
911 | painCave.isFatal = 1; | |
912 | simError(); | |
913 | } | |
# | Line 857 | Line 917 | void SimSetup::finalInfoCheck(void){ | |
917 | if (usesDipoles){ | |
918 | if (!globals->haveECR()){ | |
919 | sprintf(painCave.errMsg, | |
920 | < | "SimSetup Warning: using default value of 1/2 the smallest " |
921 | < | "box length for the electrostaticCutoffRadius.\n" |
922 | < | "I hope you have a very fast processor!\n"); |
923 | < | painCave.isFatal = 0; |
924 | < | simError(); |
925 | < | 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; |
920 | > | "SimSetup Warning: No value was set for electrostaticCutoffRadius.\n" |
921 | > | "\tOOPSE will use a default value of 15.0 angstroms" |
922 | > | "\tfor the electrostaticCutoffRadius.\n"); |
923 | > | painCave.isFatal = 0; |
924 | > | simError(); |
925 | > | theEcr = 15.0; |
926 | } | |
927 | else{ | |
928 | theEcr = globals->getECR(); | |
929 | } | |
930 | < | |
930 | > | |
931 | if (!globals->haveEST()){ | |
932 | sprintf(painCave.errMsg, | |
933 | < | "SimSetup Warning: using default value of 0.05 * the " |
934 | < | "electrostaticCutoffRadius for the " |
935 | < | "electrostaticSkinThickness\n"); |
933 | > | "SimSetup Warning: No value was set for electrostaticSkinThickness.\n" |
934 | > | "\tOOPSE will use a default value of\n" |
935 | > | "\t0.05 * electrostaticCutoffRadius\n" |
936 | > | "\tfor the electrostaticSkinThickness\n"); |
937 | painCave.isFatal = 0; | |
938 | simError(); | |
939 | theEst = 0.05 * theEcr; | |
# | Line 886 | Line 941 | void SimSetup::finalInfoCheck(void){ | |
941 | else{ | |
942 | theEst = globals->getEST(); | |
943 | } | |
944 | < | |
945 | < | info[i].setEcr(theEcr, theEst); |
944 | > | |
945 | > | info[i].setDefaultEcr(theEcr, theEst); |
946 | } | |
947 | } | |
948 | } | |
894 | – | |
949 | #ifdef IS_MPI | |
950 | strcpy(checkPointMsg, "post processing checks out"); | |
951 | MPIcheckPoint(); | |
952 | #endif // is_mpi | |
953 | } | |
954 | < | |
954 | > | |
955 | void SimSetup::initSystemCoords(void){ | |
956 | int i; | |
957 | ||
# | Line 914 | Line 968 | void SimSetup::initSystemCoords(void){ | |
968 | if (worldRank == 0){ | |
969 | #endif //is_mpi | |
970 | inName = globals->getInitialConfig(); | |
917 | – | double* tempDouble = new double[1000000]; |
971 | fileInit = new InitializeFromFile(inName); | |
972 | #ifdef IS_MPI | |
973 | } | |
# | Line 926 | Line 979 | void SimSetup::initSystemCoords(void){ | |
979 | delete fileInit; | |
980 | } | |
981 | else{ | |
982 | < | #ifdef IS_MPI |
930 | < | |
982 | > | |
983 | // no init from bass | |
984 | < | |
984 | > | |
985 | sprintf(painCave.errMsg, | |
986 | < | "Cannot intialize a parallel simulation without an initial configuration file.\n"); |
987 | < | painCave.isFatal; |
986 | > | "Cannot intialize a simulation without an initial configuration file.\n"); |
987 | > | painCave.isFatal = 1;; |
988 | simError(); | |
989 | < | |
938 | < | #else |
939 | < | |
940 | < | initFromBass(); |
941 | < | |
942 | < | |
943 | < | #endif |
989 | > | |
990 | } | |
991 | ||
992 | #ifdef IS_MPI | |
# | Line 1094 | Line 1140 | void SimSetup::createFF(void){ | |
1140 | the_ff = new EAM_FF(); | |
1141 | break; | |
1142 | ||
1143 | + | case FF_H2O: |
1144 | + | the_ff = new WATER(); |
1145 | + | break; |
1146 | + | |
1147 | default: | |
1148 | sprintf(painCave.errMsg, | |
1149 | "SimSetup Error. Unrecognized force field in case statement.\n"); | |
# | Line 1158 | Line 1208 | void SimSetup::calcSysValues(void){ | |
1208 | } | |
1209 | ||
1210 | void SimSetup::calcSysValues(void){ | |
1211 | < | int i, j, k; |
1211 | > | int i; |
1212 | ||
1213 | int* molMembershipArray; | |
1214 | ||
# | Line 1236 | Line 1286 | void SimSetup::mpiMolDivide(void){ | |
1286 | ||
1287 | if (local_atoms != info[0].n_atoms){ | |
1288 | sprintf(painCave.errMsg, | |
1289 | < | "SimSetup error: mpiSim's localAtom (%d) and SimSetup's" |
1290 | < | " localAtom (%d) are not equal.\n", |
1289 | > | "SimSetup error: mpiSim's localAtom (%d) and SimSetup's\n" |
1290 | > | "\tlocalAtom (%d) are not equal.\n", |
1291 | info[0].n_atoms, local_atoms); | |
1292 | painCave.isFatal = 1; | |
1293 | simError(); | |
# | Line 1257 | Line 1307 | void SimSetup::makeSysArrays(void){ | |
1307 | ||
1308 | ||
1309 | void SimSetup::makeSysArrays(void){ | |
1310 | < | int i, j, k, l; |
1310 | > | |
1311 | > | #ifndef IS_MPI |
1312 | > | int k, j; |
1313 | > | #endif // is_mpi |
1314 | > | int i, l; |
1315 | ||
1316 | Atom** the_atoms; | |
1317 | Molecule* the_molecules; | |
# | Line 1344 | Line 1398 | void SimSetup::makeIntegrator(void){ | |
1398 | NVT<RealIntegrator>* myNVT = NULL; | |
1399 | NPTi<NPT<RealIntegrator> >* myNPTi = NULL; | |
1400 | NPTf<NPT<RealIntegrator> >* myNPTf = NULL; | |
1401 | + | NPTxyz<NPT<RealIntegrator> >* myNPTxyz = NULL; |
1402 | ||
1403 | for (k = 0; k < nInfo; k++){ | |
1404 | switch (ensembleCase){ | |
# | Line 1374 | Line 1429 | void SimSetup::makeIntegrator(void){ | |
1429 | else{ | |
1430 | sprintf(painCave.errMsg, | |
1431 | "SimSetup error: If you use the NVT\n" | |
1432 | < | " ensemble, you must set tauThermostat.\n"); |
1432 | > | "\tensemble, you must set tauThermostat.\n"); |
1433 | painCave.isFatal = 1; | |
1434 | simError(); | |
1435 | } | |
# | Line 1397 | Line 1452 | void SimSetup::makeIntegrator(void){ | |
1452 | else{ | |
1453 | sprintf(painCave.errMsg, | |
1454 | "SimSetup error: If you use a constant pressure\n" | |
1455 | < | " ensemble, you must set targetPressure in the BASS file.\n"); |
1455 | > | "\tensemble, you must set targetPressure in the BASS file.\n"); |
1456 | painCave.isFatal = 1; | |
1457 | simError(); | |
1458 | } | |
# | Line 1407 | Line 1462 | void SimSetup::makeIntegrator(void){ | |
1462 | else{ | |
1463 | sprintf(painCave.errMsg, | |
1464 | "SimSetup error: If you use an NPT\n" | |
1465 | < | " ensemble, you must set tauThermostat.\n"); |
1465 | > | "\tensemble, you must set tauThermostat.\n"); |
1466 | painCave.isFatal = 1; | |
1467 | simError(); | |
1468 | } | |
# | Line 1417 | Line 1472 | void SimSetup::makeIntegrator(void){ | |
1472 | else{ | |
1473 | sprintf(painCave.errMsg, | |
1474 | "SimSetup error: If you use an NPT\n" | |
1475 | < | " ensemble, you must set tauBarostat.\n"); |
1475 | > | "\tensemble, you must set tauBarostat.\n"); |
1476 | painCave.isFatal = 1; | |
1477 | simError(); | |
1478 | } | |
# | Line 1440 | Line 1495 | void SimSetup::makeIntegrator(void){ | |
1495 | else{ | |
1496 | sprintf(painCave.errMsg, | |
1497 | "SimSetup error: If you use a constant pressure\n" | |
1498 | < | " ensemble, you must set targetPressure in the BASS file.\n"); |
1498 | > | "\tensemble, you must set targetPressure in the BASS file.\n"); |
1499 | painCave.isFatal = 1; | |
1500 | simError(); | |
1501 | } | |
1502 | ||
1503 | if (globals->haveTauThermostat()) | |
1504 | myNPTf->setTauThermostat(globals->getTauThermostat()); | |
1505 | + | |
1506 | else{ | |
1507 | sprintf(painCave.errMsg, | |
1508 | "SimSetup error: If you use an NPT\n" | |
1509 | < | " ensemble, you must set tauThermostat.\n"); |
1509 | > | "\tensemble, you must set tauThermostat.\n"); |
1510 | painCave.isFatal = 1; | |
1511 | simError(); | |
1512 | } | |
1513 | ||
1514 | if (globals->haveTauBarostat()) | |
1515 | myNPTf->setTauBarostat(globals->getTauBarostat()); | |
1516 | + | |
1517 | else{ | |
1518 | sprintf(painCave.errMsg, | |
1519 | "SimSetup error: If you use an NPT\n" | |
1520 | < | " ensemble, you must set tauBarostat.\n"); |
1520 | > | "\tensemble, you must set tauBarostat.\n"); |
1521 | painCave.isFatal = 1; | |
1522 | simError(); | |
1523 | } | |
1524 | ||
1525 | info->the_integrator = myNPTf; | |
1526 | + | break; |
1527 | + | |
1528 | + | case NPTxyz_ENS: |
1529 | + | if (globals->haveZconstraints()){ |
1530 | + | setupZConstraint(info[k]); |
1531 | + | myNPTxyz = new ZConstraint<NPTxyz<NPT <RealIntegrator> > >(&(info[k]), the_ff); |
1532 | + | } |
1533 | + | else |
1534 | + | myNPTxyz = new NPTxyz<NPT <RealIntegrator> >(&(info[k]), the_ff); |
1535 | + | |
1536 | + | myNPTxyz->setTargetTemp(globals->getTargetTemp()); |
1537 | + | |
1538 | + | if (globals->haveTargetPressure()) |
1539 | + | myNPTxyz->setTargetPressure(globals->getTargetPressure()); |
1540 | + | else{ |
1541 | + | sprintf(painCave.errMsg, |
1542 | + | "SimSetup error: If you use a constant pressure\n" |
1543 | + | "\tensemble, you must set targetPressure in the BASS file.\n"); |
1544 | + | painCave.isFatal = 1; |
1545 | + | simError(); |
1546 | + | } |
1547 | + | |
1548 | + | if (globals->haveTauThermostat()) |
1549 | + | myNPTxyz->setTauThermostat(globals->getTauThermostat()); |
1550 | + | else{ |
1551 | + | sprintf(painCave.errMsg, |
1552 | + | "SimSetup error: If you use an NPT\n" |
1553 | + | "\tensemble, you must set tauThermostat.\n"); |
1554 | + | painCave.isFatal = 1; |
1555 | + | simError(); |
1556 | + | } |
1557 | + | |
1558 | + | if (globals->haveTauBarostat()) |
1559 | + | myNPTxyz->setTauBarostat(globals->getTauBarostat()); |
1560 | + | else{ |
1561 | + | sprintf(painCave.errMsg, |
1562 | + | "SimSetup error: If you use an NPT\n" |
1563 | + | "\tensemble, you must set tauBarostat.\n"); |
1564 | + | painCave.isFatal = 1; |
1565 | + | simError(); |
1566 | + | } |
1567 | + | |
1568 | + | info->the_integrator = myNPTxyz; |
1569 | break; | |
1570 | ||
1571 | default: | |
# | Line 1513 | Line 1613 | void SimSetup::setupZConstraint(SimInfo& theInfo){ | |
1613 | } | |
1614 | else{ | |
1615 | sprintf(painCave.errMsg, | |
1616 | < | "ZConstraint error: If you use an ZConstraint\n" |
1617 | < | " , you must set sample time.\n"); |
1616 | > | "ZConstraint error: If you use a ZConstraint,\n" |
1617 | > | "\tyou must set zconsTime.\n"); |
1618 | painCave.isFatal = 1; | |
1619 | simError(); | |
1620 | } | |
# | Line 1529 | Line 1629 | void SimSetup::setupZConstraint(SimInfo& theInfo){ | |
1629 | else{ | |
1630 | double defaultZConsTol = 0.01; | |
1631 | sprintf(painCave.errMsg, | |
1632 | < | "ZConstraint Waring: Tolerance for z-constraint methodl is not specified\n" |
1633 | < | " , default value %f is used.\n", |
1632 | > | "ZConstraint Warning: Tolerance for z-constraint method is not specified.\n" |
1633 | > | "\tOOPSE will use a default value of %f.\n" |
1634 | > | "\tTo set the tolerance, use the zconsTol variable.\n", |
1635 | defaultZConsTol); | |
1636 | painCave.isFatal = 0; | |
1637 | simError(); | |
# | Line 1548 | Line 1649 | void SimSetup::setupZConstraint(SimInfo& theInfo){ | |
1649 | } | |
1650 | else{ | |
1651 | sprintf(painCave.errMsg, | |
1652 | < | "ZConstraint Warning: User does not set force Subtraction policy, " |
1653 | < | "PolicyByMass is used\n"); |
1652 | > | "ZConstraint Warning: No force subtraction policy was set.\n" |
1653 | > | "\tOOPSE will use PolicyByMass.\n" |
1654 | > | "\tTo set the policy, use the zconsForcePolicy variable.\n"); |
1655 | painCave.isFatal = 0; | |
1656 | simError(); | |
1657 | zconsForcePolicy->setData("BYMASS"); | |
# | Line 1593 | Line 1695 | void SimSetup::setupZConstraint(SimInfo& theInfo){ | |
1695 | //check the uniqueness of index | |
1696 | if(!zconsParaData->isIndexUnique()){ | |
1697 | sprintf(painCave.errMsg, | |
1698 | < | "ZConstraint Error: molIndex is not unique\n"); |
1698 | > | "ZConstraint Error: molIndex is not unique!\n"); |
1699 | painCave.isFatal = 1; | |
1700 | simError(); | |
1701 | } |
– | Removed lines |
+ | Added lines |
< | Changed lines |
> | Changed lines |