| 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 NPTim_ENS 4 |
| 26 |
< |
#define NPTfm_ENS 5 |
| 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 |
|
|
| 94 |
|
#endif // is_mpi |
| 95 |
|
|
| 96 |
|
void SimSetup::createSim(void){ |
| 95 |
– |
int i, j, k, globalAtomIndex; |
| 97 |
|
|
| 98 |
|
// gather all of the information from the Bass file |
| 99 |
|
|
| 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 |
| 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; |
| 557 |
|
|
| 558 |
|
|
| 559 |
|
void SimSetup::gatherInfo(void){ |
| 560 |
< |
int i, j, k; |
| 560 |
> |
int i; |
| 561 |
|
|
| 562 |
|
ensembleCase = -1; |
| 563 |
|
ffCase = -1; |
| 608 |
|
else if (!strcasecmp(ensemble, "NPTf")){ |
| 609 |
|
ensembleCase = NPTf_ENS; |
| 610 |
|
} |
| 611 |
< |
else if (!strcasecmp(ensemble, "NPTim")){ |
| 612 |
< |
ensembleCase = NPTim_ENS; |
| 609 |
< |
} |
| 610 |
< |
else if (!strcasecmp(ensemble, "NPTfm")){ |
| 611 |
< |
ensembleCase = NPTfm_ENS; |
| 611 |
> |
else if (!strcasecmp(ensemble, "NPTxyz")){ |
| 612 |
> |
ensembleCase = NPTxyz_ENS; |
| 613 |
|
} |
| 614 |
|
else{ |
| 615 |
|
sprintf(painCave.errMsg, |
| 687 |
|
|
| 688 |
|
if (globals->haveThermalTime()){ |
| 689 |
|
info[i].thermalTime = globals->getThermalTime(); |
| 690 |
+ |
} |
| 691 |
+ |
|
| 692 |
+ |
info[i].resetIntegrator = 0; |
| 693 |
+ |
if( globals->haveResetTime() ){ |
| 694 |
+ |
info[i].resetTime = globals->getResetTime(); |
| 695 |
+ |
info[i].resetIntegrator = 1; |
| 696 |
|
} |
| 697 |
|
|
| 698 |
|
// check for the temperature set flag |
| 923 |
|
if (worldRank == 0){ |
| 924 |
|
#endif //is_mpi |
| 925 |
|
inName = globals->getInitialConfig(); |
| 919 |
– |
double* tempDouble = new double[1000000]; |
| 926 |
|
fileInit = new InitializeFromFile(inName); |
| 927 |
|
#ifdef IS_MPI |
| 928 |
|
} |
| 940 |
|
|
| 941 |
|
sprintf(painCave.errMsg, |
| 942 |
|
"Cannot intialize a parallel simulation without an initial configuration file.\n"); |
| 943 |
< |
painCave.isFatal; |
| 943 |
> |
painCave.isFatal = 1;; |
| 944 |
|
simError(); |
| 945 |
|
|
| 946 |
|
#else |
| 1166 |
|
} |
| 1167 |
|
|
| 1168 |
|
void SimSetup::calcSysValues(void){ |
| 1169 |
< |
int i, j, k; |
| 1169 |
> |
int i; |
| 1170 |
|
|
| 1171 |
|
int* molMembershipArray; |
| 1172 |
|
|
| 1265 |
|
|
| 1266 |
|
|
| 1267 |
|
void SimSetup::makeSysArrays(void){ |
| 1268 |
< |
int i, j, k, l; |
| 1268 |
> |
|
| 1269 |
> |
#ifndef IS_MPI |
| 1270 |
> |
int k, j; |
| 1271 |
> |
#endif // is_mpi |
| 1272 |
> |
int i, l; |
| 1273 |
|
|
| 1274 |
|
Atom** the_atoms; |
| 1275 |
|
Molecule* the_molecules; |
| 1352 |
|
void SimSetup::makeIntegrator(void){ |
| 1353 |
|
int k; |
| 1354 |
|
|
| 1355 |
+ |
NVE<RealIntegrator>* myNVE = NULL; |
| 1356 |
|
NVT<RealIntegrator>* myNVT = NULL; |
| 1357 |
< |
NPTi<RealIntegrator>* myNPTi = NULL; |
| 1358 |
< |
NPTf<RealIntegrator>* myNPTf = NULL; |
| 1359 |
< |
NPTim<RealIntegrator>* myNPTim = NULL; |
| 1349 |
< |
NPTfm<RealIntegrator>* myNPTfm = NULL; |
| 1357 |
> |
NPTi<NPT<RealIntegrator> >* myNPTi = NULL; |
| 1358 |
> |
NPTf<NPT<RealIntegrator> >* myNPTf = NULL; |
| 1359 |
> |
NPTxyz<NPT<RealIntegrator> >* myNPTxyz = NULL; |
| 1360 |
|
|
| 1361 |
|
for (k = 0; k < nInfo; k++){ |
| 1362 |
|
switch (ensembleCase){ |
| 1363 |
|
case NVE_ENS: |
| 1364 |
|
if (globals->haveZconstraints()){ |
| 1365 |
|
setupZConstraint(info[k]); |
| 1366 |
< |
new ZConstraint<NVE<RealIntegrator> >(&(info[k]), the_ff); |
| 1366 |
> |
myNVE = new ZConstraint<NVE<RealIntegrator> >(&(info[k]), the_ff); |
| 1367 |
|
} |
| 1368 |
< |
else |
| 1369 |
< |
new NVE<RealIntegrator>(&(info[k]), the_ff); |
| 1368 |
> |
else{ |
| 1369 |
> |
myNVE = new NVE<RealIntegrator>(&(info[k]), the_ff); |
| 1370 |
> |
} |
| 1371 |
> |
|
| 1372 |
> |
info->the_integrator = myNVE; |
| 1373 |
|
break; |
| 1374 |
|
|
| 1375 |
|
case NVT_ENS: |
| 1391 |
|
painCave.isFatal = 1; |
| 1392 |
|
simError(); |
| 1393 |
|
} |
| 1394 |
+ |
|
| 1395 |
+ |
info->the_integrator = myNVT; |
| 1396 |
|
break; |
| 1397 |
|
|
| 1398 |
|
case NPTi_ENS: |
| 1399 |
|
if (globals->haveZconstraints()){ |
| 1400 |
|
setupZConstraint(info[k]); |
| 1401 |
< |
myNPTi = new ZConstraint<NPTi<RealIntegrator> >(&(info[k]), the_ff); |
| 1401 |
> |
myNPTi = new ZConstraint<NPTi<NPT <RealIntegrator> > >(&(info[k]), the_ff); |
| 1402 |
|
} |
| 1403 |
|
else |
| 1404 |
< |
myNPTi = new NPTi<RealIntegrator>(&(info[k]), the_ff); |
| 1404 |
> |
myNPTi = new NPTi<NPT<RealIntegrator> >(&(info[k]), the_ff); |
| 1405 |
|
|
| 1406 |
|
myNPTi->setTargetTemp(globals->getTargetTemp()); |
| 1407 |
|
|
| 1434 |
|
painCave.isFatal = 1; |
| 1435 |
|
simError(); |
| 1436 |
|
} |
| 1437 |
+ |
|
| 1438 |
+ |
info->the_integrator = myNPTi; |
| 1439 |
|
break; |
| 1440 |
|
|
| 1441 |
|
case NPTf_ENS: |
| 1442 |
|
if (globals->haveZconstraints()){ |
| 1443 |
|
setupZConstraint(info[k]); |
| 1444 |
< |
myNPTf = new ZConstraint<NPTf<RealIntegrator> >(&(info[k]), the_ff); |
| 1444 |
> |
myNPTf = new ZConstraint<NPTf<NPT <RealIntegrator> > >(&(info[k]), the_ff); |
| 1445 |
|
} |
| 1446 |
|
else |
| 1447 |
< |
myNPTf = new NPTf<RealIntegrator>(&(info[k]), the_ff); |
| 1447 |
> |
myNPTf = new NPTf<NPT <RealIntegrator> >(&(info[k]), the_ff); |
| 1448 |
|
|
| 1449 |
|
myNPTf->setTargetTemp(globals->getTargetTemp()); |
| 1450 |
|
|
| 1477 |
|
painCave.isFatal = 1; |
| 1478 |
|
simError(); |
| 1479 |
|
} |
| 1463 |
– |
break; |
| 1480 |
|
|
| 1481 |
< |
case NPTim_ENS: |
| 1466 |
< |
if (globals->haveZconstraints()){ |
| 1467 |
< |
setupZConstraint(info[k]); |
| 1468 |
< |
myNPTim = new ZConstraint<NPTim<RealIntegrator> >(&(info[k]), the_ff); |
| 1469 |
< |
} |
| 1470 |
< |
else |
| 1471 |
< |
myNPTim = new NPTim<RealIntegrator>(&(info[k]), the_ff); |
| 1472 |
< |
|
| 1473 |
< |
myNPTim->setTargetTemp(globals->getTargetTemp()); |
| 1474 |
< |
|
| 1475 |
< |
if (globals->haveTargetPressure()) |
| 1476 |
< |
myNPTim->setTargetPressure(globals->getTargetPressure()); |
| 1477 |
< |
else{ |
| 1478 |
< |
sprintf(painCave.errMsg, |
| 1479 |
< |
"SimSetup error: If you use a constant pressure\n" |
| 1480 |
< |
" ensemble, you must set targetPressure in the BASS file.\n"); |
| 1481 |
< |
painCave.isFatal = 1; |
| 1482 |
< |
simError(); |
| 1483 |
< |
} |
| 1484 |
< |
|
| 1485 |
< |
if (globals->haveTauThermostat()) |
| 1486 |
< |
myNPTim->setTauThermostat(globals->getTauThermostat()); |
| 1487 |
< |
else{ |
| 1488 |
< |
sprintf(painCave.errMsg, |
| 1489 |
< |
"SimSetup error: If you use an NPT\n" |
| 1490 |
< |
" ensemble, you must set tauThermostat.\n"); |
| 1491 |
< |
painCave.isFatal = 1; |
| 1492 |
< |
simError(); |
| 1493 |
< |
} |
| 1494 |
< |
|
| 1495 |
< |
if (globals->haveTauBarostat()) |
| 1496 |
< |
myNPTim->setTauBarostat(globals->getTauBarostat()); |
| 1497 |
< |
else{ |
| 1498 |
< |
sprintf(painCave.errMsg, |
| 1499 |
< |
"SimSetup error: If you use an NPT\n" |
| 1500 |
< |
" ensemble, you must set tauBarostat.\n"); |
| 1501 |
< |
painCave.isFatal = 1; |
| 1502 |
< |
simError(); |
| 1503 |
< |
} |
| 1481 |
> |
info->the_integrator = myNPTf; |
| 1482 |
|
break; |
| 1483 |
|
|
| 1484 |
< |
case NPTfm_ENS: |
| 1484 |
> |
case NPTxyz_ENS: |
| 1485 |
|
if (globals->haveZconstraints()){ |
| 1486 |
|
setupZConstraint(info[k]); |
| 1487 |
< |
myNPTfm = new ZConstraint<NPTfm<RealIntegrator> >(&(info[k]), the_ff); |
| 1487 |
> |
myNPTxyz = new ZConstraint<NPTxyz<NPT <RealIntegrator> > >(&(info[k]), the_ff); |
| 1488 |
|
} |
| 1489 |
|
else |
| 1490 |
< |
myNPTfm = new NPTfm<RealIntegrator>(&(info[k]), the_ff); |
| 1490 |
> |
myNPTxyz = new NPTxyz<NPT <RealIntegrator> >(&(info[k]), the_ff); |
| 1491 |
|
|
| 1492 |
< |
myNPTfm->setTargetTemp(globals->getTargetTemp()); |
| 1492 |
> |
myNPTxyz->setTargetTemp(globals->getTargetTemp()); |
| 1493 |
|
|
| 1494 |
|
if (globals->haveTargetPressure()) |
| 1495 |
< |
myNPTfm->setTargetPressure(globals->getTargetPressure()); |
| 1495 |
> |
myNPTxyz->setTargetPressure(globals->getTargetPressure()); |
| 1496 |
|
else{ |
| 1497 |
|
sprintf(painCave.errMsg, |
| 1498 |
|
"SimSetup error: If you use a constant pressure\n" |
| 1499 |
|
" ensemble, you must set targetPressure in the BASS file.\n"); |
| 1500 |
|
painCave.isFatal = 1; |
| 1501 |
|
simError(); |
| 1502 |
< |
} |
| 1502 |
> |
} |
| 1503 |
|
|
| 1504 |
|
if (globals->haveTauThermostat()) |
| 1505 |
< |
myNPTfm->setTauThermostat(globals->getTauThermostat()); |
| 1505 |
> |
myNPTxyz->setTauThermostat(globals->getTauThermostat()); |
| 1506 |
|
else{ |
| 1507 |
|
sprintf(painCave.errMsg, |
| 1508 |
|
"SimSetup error: If you use an NPT\n" |
| 1512 |
|
} |
| 1513 |
|
|
| 1514 |
|
if (globals->haveTauBarostat()) |
| 1515 |
< |
myNPTfm->setTauBarostat(globals->getTauBarostat()); |
| 1515 |
> |
myNPTxyz->setTauBarostat(globals->getTauBarostat()); |
| 1516 |
|
else{ |
| 1517 |
|
sprintf(painCave.errMsg, |
| 1518 |
|
"SimSetup error: If you use an NPT\n" |
| 1520 |
|
painCave.isFatal = 1; |
| 1521 |
|
simError(); |
| 1522 |
|
} |
| 1523 |
+ |
|
| 1524 |
+ |
info->the_integrator = myNPTxyz; |
| 1525 |
|
break; |
| 1526 |
|
|
| 1527 |
|
default: |
| 1595 |
|
} |
| 1596 |
|
theInfo.addProperty(zconsTol); |
| 1597 |
|
|
| 1598 |
< |
//set Force Substraction Policy |
| 1598 |
> |
//set Force Subtraction Policy |
| 1599 |
|
StringData* zconsForcePolicy = new StringData(); |
| 1600 |
|
zconsForcePolicy->setID(ZCONSFORCEPOLICY_ID); |
| 1601 |
|
|
| 1604 |
|
} |
| 1605 |
|
else{ |
| 1606 |
|
sprintf(painCave.errMsg, |
| 1607 |
< |
"ZConstraint Warning: User does not set force substraction policy, " |
| 1607 |
> |
"ZConstraint Warning: User does not set force Subtraction policy, " |
| 1608 |
|
"PolicyByMass is used\n"); |
| 1609 |
|
painCave.isFatal = 0; |
| 1610 |
|
simError(); |