| 29 | 
  | 
#define FF_EAM  2 | 
| 30 | 
  | 
 | 
| 31 | 
  | 
using namespace std; | 
| 32 | 
+ | 
 | 
| 33 | 
+ | 
/** | 
| 34 | 
+ | 
 * Check whether dividend is divisble by divisor or not | 
| 35 | 
+ | 
 */ | 
| 36 | 
+ | 
bool isDivisible(double dividend, double divisor){ | 
| 37 | 
+ | 
  double tolerance = 0.000001; | 
| 38 | 
+ | 
  double quotient; | 
| 39 | 
+ | 
  double diff; | 
| 40 | 
+ | 
  int intQuotient; | 
| 41 | 
+ | 
   | 
| 42 | 
+ | 
  quotient = dividend / divisor; | 
| 43 | 
+ | 
 | 
| 44 | 
+ | 
  if (quotient < 0) | 
| 45 | 
+ | 
    quotient = -quotient; | 
| 46 | 
+ | 
 | 
| 47 | 
+ | 
  intQuotient = int (quotient + tolerance); | 
| 48 | 
+ | 
 | 
| 49 | 
+ | 
  diff = fabs(fabs(dividend) - intQuotient  * fabs(divisor)); | 
| 50 | 
+ | 
 | 
| 51 | 
+ | 
  if (diff <= tolerance) | 
| 52 | 
+ | 
    return true; | 
| 53 | 
+ | 
  else | 
| 54 | 
+ | 
    return false;   | 
| 55 | 
+ | 
} | 
| 56 | 
  | 
 | 
| 57 | 
  | 
SimSetup::SimSetup(){ | 
| 58 | 
  | 
   | 
| 689 | 
  | 
            " Please give nMol in the components.\n"); | 
| 690 | 
  | 
    painCave.isFatal = 1; | 
| 691 | 
  | 
    simError(); | 
| 692 | 
+ | 
  } | 
| 693 | 
+ | 
 | 
| 694 | 
+ | 
  //check whether sample time, status time, thermal time and reset time are divisble by dt | 
| 695 | 
+ | 
  if (!isDivisible(globals->getSampleTime(), globals->getDt())){ | 
| 696 | 
+ | 
    sprintf(painCave.errMsg, | 
| 697 | 
+ | 
              "Sample time is not divisible by dt \n"); | 
| 698 | 
+ | 
    painCave.isFatal = 0; | 
| 699 | 
+ | 
    simError();     | 
| 700 | 
+ | 
  } | 
| 701 | 
+ | 
 | 
| 702 | 
+ | 
  if (globals->haveStatusTime() && !isDivisible(globals->getSampleTime(), globals->getDt())){ | 
| 703 | 
+ | 
    sprintf(painCave.errMsg, | 
| 704 | 
+ | 
              "Status time is not divisible by dt\n"); | 
| 705 | 
+ | 
    painCave.isFatal = 0; | 
| 706 | 
+ | 
    simError();     | 
| 707 | 
  | 
  } | 
| 708 | 
  | 
 | 
| 709 | 
+ | 
  if (globals->haveThermalTime() && !isDivisible(globals->getThermalTime(), globals->getDt())){ | 
| 710 | 
+ | 
    sprintf(painCave.errMsg, | 
| 711 | 
+ | 
              "Thermal time is not divisible by dt\n"); | 
| 712 | 
+ | 
    painCave.isFatal = 0; | 
| 713 | 
+ | 
    simError();     | 
| 714 | 
+ | 
  }   | 
| 715 | 
+ | 
 | 
| 716 | 
+ | 
  if (globals->haveResetTime() && !isDivisible(globals->getResetTime(), globals->getDt())){ | 
| 717 | 
+ | 
    sprintf(painCave.errMsg, | 
| 718 | 
+ | 
              "Reset time is not divisible by dt\n"); | 
| 719 | 
+ | 
    painCave.isFatal = 0; | 
| 720 | 
+ | 
    simError();     | 
| 721 | 
+ | 
  }  | 
| 722 | 
+ | 
 | 
| 723 | 
  | 
  // set the status, sample, and thermal kick times | 
| 724 | 
  | 
 | 
| 725 | 
  | 
  for (i = 0; i < nInfo; i++){ | 
| 749 | 
  | 
    } | 
| 750 | 
  | 
 | 
| 751 | 
  | 
    // check for the temperature set flag | 
| 752 | 
< | 
 | 
| 752 | 
> | 
     | 
| 753 | 
  | 
    if (globals->haveTempSet()) | 
| 754 | 
  | 
      info[i].setTemp = globals->getTempSet(); | 
| 755 | 
  | 
 | 
| 756 | 
< | 
    // get some of the tricky things that may still be in the globals | 
| 756 | 
> | 
    // check for the extended State init | 
| 757 | 
  | 
 | 
| 758 | 
< | 
    double boxVector[3]; | 
| 759 | 
< | 
    if (globals->haveBox()){ | 
| 760 | 
< | 
      boxVector[0] = globals->getBox(); | 
| 708 | 
< | 
      boxVector[1] = globals->getBox(); | 
| 709 | 
< | 
      boxVector[2] = globals->getBox(); | 
| 710 | 
< | 
 | 
| 711 | 
< | 
      info[i].setBox(boxVector); | 
| 712 | 
< | 
    } | 
| 713 | 
< | 
    else if (globals->haveDensity()){ | 
| 714 | 
< | 
      double vol; | 
| 715 | 
< | 
      vol = (double) tot_nmol / globals->getDensity(); | 
| 716 | 
< | 
      boxVector[0] = pow(vol, (1.0 / 3.0)); | 
| 717 | 
< | 
      boxVector[1] = boxVector[0]; | 
| 718 | 
< | 
      boxVector[2] = boxVector[0]; | 
| 719 | 
< | 
 | 
| 720 | 
< | 
      info[i].setBox(boxVector); | 
| 721 | 
< | 
    } | 
| 722 | 
< | 
    else{ | 
| 723 | 
< | 
      if (!globals->haveBoxX()){ | 
| 724 | 
< | 
        sprintf(painCave.errMsg, | 
| 725 | 
< | 
                "SimSetup error, no periodic BoxX size given.\n"); | 
| 726 | 
< | 
        painCave.isFatal = 1; | 
| 727 | 
< | 
        simError(); | 
| 728 | 
< | 
      } | 
| 729 | 
< | 
      boxVector[0] = globals->getBoxX(); | 
| 730 | 
< | 
 | 
| 731 | 
< | 
      if (!globals->haveBoxY()){ | 
| 732 | 
< | 
        sprintf(painCave.errMsg, | 
| 733 | 
< | 
                "SimSetup error, no periodic BoxY size given.\n"); | 
| 734 | 
< | 
        painCave.isFatal = 1; | 
| 735 | 
< | 
        simError(); | 
| 736 | 
< | 
      } | 
| 737 | 
< | 
      boxVector[1] = globals->getBoxY(); | 
| 738 | 
< | 
 | 
| 739 | 
< | 
      if (!globals->haveBoxZ()){ | 
| 740 | 
< | 
        sprintf(painCave.errMsg, | 
| 741 | 
< | 
                "SimSetup error, no periodic BoxZ size given.\n"); | 
| 742 | 
< | 
        painCave.isFatal = 1; | 
| 743 | 
< | 
        simError(); | 
| 744 | 
< | 
      } | 
| 745 | 
< | 
      boxVector[2] = globals->getBoxZ(); | 
| 746 | 
< | 
 | 
| 747 | 
< | 
      info[i].setBox(boxVector); | 
| 748 | 
< | 
    } | 
| 758 | 
> | 
    info[i].useInitXSstate = globals->getUseInitXSstate(); | 
| 759 | 
> | 
    info[i].orthoTolerance = globals->getOrthoBoxTolerance(); | 
| 760 | 
> | 
     | 
| 761 | 
  | 
  } | 
| 762 | 
< | 
 | 
| 762 | 
> | 
   | 
| 763 | 
  | 
  //setup seed for random number generator | 
| 764 | 
  | 
  int seedValue; | 
| 765 | 
  | 
 | 
| 834 | 
  | 
 | 
| 835 | 
  | 
      if (!globals->haveECR()){ | 
| 836 | 
  | 
        sprintf(painCave.errMsg, | 
| 837 | 
< | 
                "SimSetup Warning: using default value of 1/2 the smallest " | 
| 838 | 
< | 
                "box length for the electrostaticCutoffRadius.\n" | 
| 827 | 
< | 
                "I hope you have a very fast processor!\n"); | 
| 837 | 
> | 
                "SimSetup Warning: using default value of 15.0 angstroms" | 
| 838 | 
> | 
                "box length for the electrostaticCutoffRadius.\n"); | 
| 839 | 
  | 
        painCave.isFatal = 0; | 
| 840 | 
  | 
        simError(); | 
| 841 | 
< | 
        double smallest; | 
| 831 | 
< | 
        smallest = info[i].boxL[0]; | 
| 832 | 
< | 
        if (info[i].boxL[1] <= smallest) | 
| 833 | 
< | 
          smallest = info[i].boxL[1]; | 
| 834 | 
< | 
        if (info[i].boxL[2] <= smallest) | 
| 835 | 
< | 
          smallest = info[i].boxL[2]; | 
| 836 | 
< | 
        theEcr = 0.5 * smallest; | 
| 841 | 
> | 
        theEcr = 15.0; | 
| 842 | 
  | 
      } | 
| 843 | 
  | 
      else{ | 
| 844 | 
  | 
        theEcr = globals->getECR(); | 
| 856 | 
  | 
        theEst = globals->getEST(); | 
| 857 | 
  | 
      } | 
| 858 | 
  | 
 | 
| 859 | 
< | 
      info[i].setEcr(theEcr, theEst); | 
| 859 | 
> | 
      info[i].setDefaultEcr(theEcr, theEst); | 
| 860 | 
  | 
 | 
| 861 | 
  | 
      if (!globals->haveDielectric()){ | 
| 862 | 
  | 
        sprintf(painCave.errMsg, | 
| 871 | 
  | 
      if (usesDipoles){ | 
| 872 | 
  | 
        if (!globals->haveECR()){ | 
| 873 | 
  | 
          sprintf(painCave.errMsg, | 
| 874 | 
< | 
                  "SimSetup Warning: using default value of 1/2 the smallest " | 
| 875 | 
< | 
                  "box length for the electrostaticCutoffRadius.\n" | 
| 876 | 
< | 
                  "I hope you have a very fast processor!\n"); | 
| 877 | 
< | 
          painCave.isFatal = 0; | 
| 878 | 
< | 
          simError(); | 
| 874 | 
< | 
          double smallest; | 
| 875 | 
< | 
          smallest = info[i].boxL[0]; | 
| 876 | 
< | 
          if (info[i].boxL[1] <= smallest) | 
| 877 | 
< | 
            smallest = info[i].boxL[1]; | 
| 878 | 
< | 
          if (info[i].boxL[2] <= smallest) | 
| 879 | 
< | 
            smallest = info[i].boxL[2]; | 
| 880 | 
< | 
          theEcr = 0.5 * smallest; | 
| 874 | 
> | 
                  "SimSetup Warning: using default value of 15.0 angstroms" | 
| 875 | 
> | 
                  "box length for the electrostaticCutoffRadius.\n"); | 
| 876 | 
> | 
          painCave.isFatal = 0; | 
| 877 | 
> | 
          simError(); | 
| 878 | 
> | 
          theEcr = 15.0; | 
| 879 | 
  | 
        } | 
| 880 | 
  | 
        else{ | 
| 881 | 
  | 
          theEcr = globals->getECR(); | 
| 882 | 
  | 
        } | 
| 883 | 
< | 
 | 
| 883 | 
> | 
         | 
| 884 | 
  | 
        if (!globals->haveEST()){ | 
| 885 | 
  | 
          sprintf(painCave.errMsg, | 
| 886 | 
  | 
                  "SimSetup Warning: using default value of 0.05 * the " | 
| 893 | 
  | 
        else{ | 
| 894 | 
  | 
          theEst = globals->getEST(); | 
| 895 | 
  | 
        } | 
| 896 | 
< | 
 | 
| 897 | 
< | 
        info[i].setEcr(theEcr, theEst); | 
| 896 | 
> | 
         | 
| 897 | 
> | 
        info[i].setDefaultEcr(theEcr, theEst); | 
| 898 | 
  | 
      } | 
| 899 | 
  | 
    } | 
| 900 | 
  | 
  } | 
| 903 | 
– | 
 | 
| 901 | 
  | 
#ifdef IS_MPI | 
| 902 | 
  | 
  strcpy(checkPointMsg, "post processing checks out"); | 
| 903 | 
  | 
  MPIcheckPoint(); | 
| 904 | 
  | 
#endif // is_mpi | 
| 905 | 
  | 
} | 
| 906 | 
< | 
 | 
| 906 | 
> | 
   | 
| 907 | 
  | 
void SimSetup::initSystemCoords(void){ | 
| 908 | 
  | 
  int i; | 
| 909 | 
  | 
 | 
| 931 | 
  | 
    delete fileInit; | 
| 932 | 
  | 
  } | 
| 933 | 
  | 
  else{ | 
| 934 | 
< | 
#ifdef IS_MPI  | 
| 938 | 
< | 
 | 
| 934 | 
> | 
     | 
| 935 | 
  | 
    // no init from bass | 
| 936 | 
< | 
 | 
| 936 | 
> | 
     | 
| 937 | 
  | 
    sprintf(painCave.errMsg, | 
| 938 | 
< | 
            "Cannot intialize a parallel simulation without an initial configuration file.\n"); | 
| 938 | 
> | 
            "Cannot intialize a simulation without an initial configuration file.\n"); | 
| 939 | 
  | 
    painCave.isFatal = 1;; | 
| 940 | 
  | 
    simError(); | 
| 941 | 
< | 
 | 
| 946 | 
< | 
#else | 
| 947 | 
< | 
 | 
| 948 | 
< | 
    initFromBass(); | 
| 949 | 
< | 
 | 
| 950 | 
< | 
 | 
| 951 | 
< | 
#endif | 
| 941 | 
> | 
     | 
| 942 | 
  | 
  } | 
| 943 | 
  | 
 | 
| 944 | 
  | 
#ifdef IS_MPI | 
| 1450 | 
  | 
 | 
| 1451 | 
  | 
        if (globals->haveTauThermostat()) | 
| 1452 | 
  | 
          myNPTf->setTauThermostat(globals->getTauThermostat()); | 
| 1453 | 
+ | 
 | 
| 1454 | 
  | 
        else{ | 
| 1455 | 
  | 
          sprintf(painCave.errMsg, | 
| 1456 | 
  | 
                  "SimSetup error: If you use an NPT\n" | 
| 1461 | 
  | 
 | 
| 1462 | 
  | 
        if (globals->haveTauBarostat()) | 
| 1463 | 
  | 
          myNPTf->setTauBarostat(globals->getTauBarostat()); | 
| 1464 | 
+ | 
 | 
| 1465 | 
  | 
        else{ | 
| 1466 | 
  | 
          sprintf(painCave.errMsg, | 
| 1467 | 
  | 
                  "SimSetup error: If you use an NPT\n" |