| 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" |
| 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 |
| 31 |
|
using namespace std; |
| 32 |
|
|
| 33 |
|
SimSetup::SimSetup(){ |
| 34 |
+ |
|
| 35 |
+ |
initSuspend = false; |
| 36 |
|
isInfoArray = 0; |
| 37 |
|
nInfo = 1; |
| 38 |
|
|
| 55 |
|
info = the_info; |
| 56 |
|
nInfo = theNinfo; |
| 57 |
|
isInfoArray = 1; |
| 58 |
+ |
initSuspend = true; |
| 59 |
|
} |
| 60 |
|
|
| 61 |
|
|
| 103 |
|
|
| 104 |
|
sysObjectsCreation(); |
| 105 |
|
|
| 102 |
– |
// check on the post processing info |
| 103 |
– |
|
| 104 |
– |
finalInfoCheck(); |
| 105 |
– |
|
| 106 |
|
// initialize the system coordinates |
| 107 |
|
|
| 108 |
< |
if (!isInfoArray){ |
| 108 |
> |
if ( !initSuspend ){ |
| 109 |
|
initSystemCoords(); |
| 110 |
|
|
| 111 |
|
if( !(globals->getUseInitTime()) ) |
| 112 |
|
info[0].currentTime = 0.0; |
| 113 |
|
} |
| 114 |
|
|
| 115 |
+ |
// check on the post processing info |
| 116 |
+ |
|
| 117 |
+ |
finalInfoCheck(); |
| 118 |
+ |
|
| 119 |
|
// make the output filenames |
| 120 |
|
|
| 121 |
|
makeOutNames(); |
| 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, " |
| 696 |
|
} |
| 697 |
|
|
| 698 |
|
// check for the temperature set flag |
| 699 |
< |
|
| 699 |
> |
|
| 700 |
|
if (globals->haveTempSet()) |
| 701 |
|
info[i].setTemp = globals->getTempSet(); |
| 702 |
< |
|
| 696 |
< |
// get some of the tricky things that may still be in the globals |
| 697 |
< |
|
| 698 |
< |
double boxVector[3]; |
| 699 |
< |
if (globals->haveBox()){ |
| 700 |
< |
boxVector[0] = globals->getBox(); |
| 701 |
< |
boxVector[1] = globals->getBox(); |
| 702 |
< |
boxVector[2] = globals->getBox(); |
| 703 |
< |
|
| 704 |
< |
info[i].setBox(boxVector); |
| 705 |
< |
} |
| 706 |
< |
else if (globals->haveDensity()){ |
| 707 |
< |
double vol; |
| 708 |
< |
vol = (double) tot_nmol / globals->getDensity(); |
| 709 |
< |
boxVector[0] = pow(vol, (1.0 / 3.0)); |
| 710 |
< |
boxVector[1] = boxVector[0]; |
| 711 |
< |
boxVector[2] = boxVector[0]; |
| 712 |
< |
|
| 713 |
< |
info[i].setBox(boxVector); |
| 714 |
< |
} |
| 715 |
< |
else{ |
| 716 |
< |
if (!globals->haveBoxX()){ |
| 717 |
< |
sprintf(painCave.errMsg, |
| 718 |
< |
"SimSetup error, no periodic BoxX size given.\n"); |
| 719 |
< |
painCave.isFatal = 1; |
| 720 |
< |
simError(); |
| 721 |
< |
} |
| 722 |
< |
boxVector[0] = globals->getBoxX(); |
| 723 |
< |
|
| 724 |
< |
if (!globals->haveBoxY()){ |
| 725 |
< |
sprintf(painCave.errMsg, |
| 726 |
< |
"SimSetup error, no periodic BoxY size given.\n"); |
| 727 |
< |
painCave.isFatal = 1; |
| 728 |
< |
simError(); |
| 729 |
< |
} |
| 730 |
< |
boxVector[1] = globals->getBoxY(); |
| 731 |
< |
|
| 732 |
< |
if (!globals->haveBoxZ()){ |
| 733 |
< |
sprintf(painCave.errMsg, |
| 734 |
< |
"SimSetup error, no periodic BoxZ size given.\n"); |
| 735 |
< |
painCave.isFatal = 1; |
| 736 |
< |
simError(); |
| 737 |
< |
} |
| 738 |
< |
boxVector[2] = globals->getBoxZ(); |
| 739 |
< |
|
| 740 |
< |
info[i].setBox(boxVector); |
| 741 |
< |
} |
| 702 |
> |
|
| 703 |
|
} |
| 704 |
< |
|
| 704 |
> |
|
| 705 |
|
//setup seed for random number generator |
| 706 |
|
int seedValue; |
| 707 |
|
|
| 805 |
|
theEst = globals->getEST(); |
| 806 |
|
} |
| 807 |
|
|
| 808 |
< |
info[i].setEcr(theEcr, theEst); |
| 808 |
> |
info[i].setDefaultEcr(theEcr, theEst); |
| 809 |
|
|
| 810 |
|
if (!globals->haveDielectric()){ |
| 811 |
|
sprintf(painCave.errMsg, |
| 850 |
|
theEst = globals->getEST(); |
| 851 |
|
} |
| 852 |
|
|
| 853 |
< |
info[i].setEcr(theEcr, theEst); |
| 853 |
> |
info[i].setDefaultEcr(theEcr, theEst); |
| 854 |
|
} |
| 855 |
|
} |
| 856 |
|
} |
| 1310 |
|
NVT<RealIntegrator>* myNVT = NULL; |
| 1311 |
|
NPTi<NPT<RealIntegrator> >* myNPTi = NULL; |
| 1312 |
|
NPTf<NPT<RealIntegrator> >* myNPTf = NULL; |
| 1313 |
+ |
NPTxyz<NPT<RealIntegrator> >* myNPTxyz = NULL; |
| 1314 |
|
|
| 1315 |
|
for (k = 0; k < nInfo; k++){ |
| 1316 |
|
switch (ensembleCase){ |
| 1414 |
|
|
| 1415 |
|
if (globals->haveTauThermostat()) |
| 1416 |
|
myNPTf->setTauThermostat(globals->getTauThermostat()); |
| 1417 |
+ |
|
| 1418 |
|
else{ |
| 1419 |
|
sprintf(painCave.errMsg, |
| 1420 |
|
"SimSetup error: If you use an NPT\n" |
| 1425 |
|
|
| 1426 |
|
if (globals->haveTauBarostat()) |
| 1427 |
|
myNPTf->setTauBarostat(globals->getTauBarostat()); |
| 1428 |
+ |
|
| 1429 |
|
else{ |
| 1430 |
|
sprintf(painCave.errMsg, |
| 1431 |
|
"SimSetup error: If you use an NPT\n" |
| 1437 |
|
info->the_integrator = myNPTf; |
| 1438 |
|
break; |
| 1439 |
|
|
| 1440 |
+ |
case NPTxyz_ENS: |
| 1441 |
+ |
if (globals->haveZconstraints()){ |
| 1442 |
+ |
setupZConstraint(info[k]); |
| 1443 |
+ |
myNPTxyz = new ZConstraint<NPTxyz<NPT <RealIntegrator> > >(&(info[k]), the_ff); |
| 1444 |
+ |
} |
| 1445 |
+ |
else |
| 1446 |
+ |
myNPTxyz = new NPTxyz<NPT <RealIntegrator> >(&(info[k]), the_ff); |
| 1447 |
+ |
|
| 1448 |
+ |
myNPTxyz->setTargetTemp(globals->getTargetTemp()); |
| 1449 |
+ |
|
| 1450 |
+ |
if (globals->haveTargetPressure()) |
| 1451 |
+ |
myNPTxyz->setTargetPressure(globals->getTargetPressure()); |
| 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"); |
| 1456 |
+ |
painCave.isFatal = 1; |
| 1457 |
+ |
simError(); |
| 1458 |
+ |
} |
| 1459 |
+ |
|
| 1460 |
+ |
if (globals->haveTauThermostat()) |
| 1461 |
+ |
myNPTxyz->setTauThermostat(globals->getTauThermostat()); |
| 1462 |
+ |
else{ |
| 1463 |
+ |
sprintf(painCave.errMsg, |
| 1464 |
+ |
"SimSetup error: If you use an NPT\n" |
| 1465 |
+ |
" ensemble, you must set tauThermostat.\n"); |
| 1466 |
+ |
painCave.isFatal = 1; |
| 1467 |
+ |
simError(); |
| 1468 |
+ |
} |
| 1469 |
+ |
|
| 1470 |
+ |
if (globals->haveTauBarostat()) |
| 1471 |
+ |
myNPTxyz->setTauBarostat(globals->getTauBarostat()); |
| 1472 |
+ |
else{ |
| 1473 |
+ |
sprintf(painCave.errMsg, |
| 1474 |
+ |
"SimSetup error: If you use an NPT\n" |
| 1475 |
+ |
" ensemble, you must set tauBarostat.\n"); |
| 1476 |
+ |
painCave.isFatal = 1; |
| 1477 |
+ |
simError(); |
| 1478 |
+ |
} |
| 1479 |
+ |
|
| 1480 |
+ |
info->the_integrator = myNPTxyz; |
| 1481 |
+ |
break; |
| 1482 |
+ |
|
| 1483 |
|
default: |
| 1484 |
|
sprintf(painCave.errMsg, |
| 1485 |
|
"SimSetup Error. Unrecognized ensemble in case statement.\n"); |