| 55 | 
  | 
#include "UseTheForce/fCutoffPolicy.h" | 
| 56 | 
  | 
#include "UseTheForce/DarkSide/fElectrostaticSummationMethod.h" | 
| 57 | 
  | 
#include "UseTheForce/doForces_interface.h" | 
| 58 | 
+ | 
#include "UseTheForce/DarkSide/electrostatic_interface.h" | 
| 59 | 
  | 
#include "UseTheForce/notifyCutoffs_interface.h" | 
| 60 | 
  | 
#include "utils/MemoryUtils.hpp" | 
| 61 | 
  | 
#include "utils/simError.h" | 
| 83 | 
  | 
      MoleculeStamp* molStamp; | 
| 84 | 
  | 
      int nMolWithSameStamp; | 
| 85 | 
  | 
      int nCutoffAtoms = 0; // number of atoms belong to cutoff groups | 
| 86 | 
< | 
      int nGroups = 0;          //total cutoff groups defined in meta-data file | 
| 86 | 
> | 
      int nGroups = 0;      //total cutoff groups defined in meta-data file | 
| 87 | 
  | 
      CutoffGroupStamp* cgStamp;     | 
| 88 | 
  | 
      RigidBodyStamp* rbStamp; | 
| 89 | 
  | 
      int nRigidAtoms = 0; | 
| 108 | 
  | 
        } | 
| 109 | 
  | 
 | 
| 110 | 
  | 
        nGroups += nCutoffGroupsInStamp * nMolWithSameStamp; | 
| 111 | 
+ | 
 | 
| 112 | 
  | 
        nCutoffAtoms += nAtomsInGroups * nMolWithSameStamp;             | 
| 113 | 
  | 
 | 
| 114 | 
  | 
        //calculate atoms in rigid bodies | 
| 125 | 
  | 
         | 
| 126 | 
  | 
      } | 
| 127 | 
  | 
 | 
| 128 | 
< | 
      //every free atom (atom does not belong to cutoff groups) is a cutoff group | 
| 129 | 
< | 
      //therefore the total number of cutoff groups in the system is equal to  | 
| 130 | 
< | 
      //the total number of atoms minus number of atoms belong to cutoff group defined in meta-data | 
| 131 | 
< | 
      //file plus the number of cutoff groups defined in meta-data file | 
| 128 | 
> | 
      //every free atom (atom does not belong to cutoff groups) is a cutoff  | 
| 129 | 
> | 
      //group therefore the total number of cutoff groups in the system is  | 
| 130 | 
> | 
      //equal to the total number of atoms minus number of atoms belong to  | 
| 131 | 
> | 
      //cutoff group defined in meta-data file plus the number of cutoff  | 
| 132 | 
> | 
      //groups defined in meta-data file | 
| 133 | 
  | 
      nGlobalCutoffGroups_ = nGlobalAtoms_ - nCutoffAtoms + nGroups; | 
| 134 | 
  | 
 | 
| 135 | 
< | 
      //every free atom (atom does not belong to rigid bodies) is an integrable object | 
| 136 | 
< | 
      //therefore the total number of  integrable objects in the system is equal to  | 
| 137 | 
< | 
      //the total number of atoms minus number of atoms belong to  rigid body defined in meta-data | 
| 138 | 
< | 
      //file plus the number of  rigid bodies defined in meta-data file | 
| 139 | 
< | 
      nGlobalIntegrableObjects_ = nGlobalAtoms_ - nRigidAtoms + nGlobalRigidBodies_; | 
| 140 | 
< | 
 | 
| 135 | 
> | 
      //every free atom (atom does not belong to rigid bodies) is an  | 
| 136 | 
> | 
      //integrable object therefore the total number of integrable objects  | 
| 137 | 
> | 
      //in the system is equal to the total number of atoms minus number of  | 
| 138 | 
> | 
      //atoms belong to rigid body defined in meta-data file plus the number  | 
| 139 | 
> | 
      //of rigid bodies defined in meta-data file | 
| 140 | 
> | 
      nGlobalIntegrableObjects_ = nGlobalAtoms_ - nRigidAtoms  | 
| 141 | 
> | 
                                                + nGlobalRigidBodies_; | 
| 142 | 
> | 
   | 
| 143 | 
  | 
      nGlobalMols_ = molStampIds_.size(); | 
| 144 | 
  | 
 | 
| 145 | 
  | 
#ifdef IS_MPI     | 
| 526 | 
  | 
    int useDirectionalAtom = 0;     | 
| 527 | 
  | 
    int useElectrostatics = 0; | 
| 528 | 
  | 
    //usePBC and useRF are from simParams | 
| 529 | 
< | 
    int usePBC = simParams_->getPBC(); | 
| 529 | 
> | 
    int usePBC = simParams_->getUsePeriodicBoundaryConditions(); | 
| 530 | 
> | 
    int useRF; | 
| 531 | 
> | 
    std::string myMethod; | 
| 532 | 
  | 
 | 
| 533 | 
+ | 
    // set the useRF logical | 
| 534 | 
+ | 
    useRF = 0; | 
| 535 | 
+ | 
 | 
| 536 | 
+ | 
 | 
| 537 | 
+ | 
    if (simParams_->haveElectrostaticSummationMethod()) { | 
| 538 | 
+ | 
      std::string myMethod = simParams_->getElectrostaticSummationMethod(); | 
| 539 | 
+ | 
      toUpper(myMethod); | 
| 540 | 
+ | 
      if (myMethod == "REACTION_FIELD") { | 
| 541 | 
+ | 
        useRF=1; | 
| 542 | 
+ | 
      } | 
| 543 | 
+ | 
    } | 
| 544 | 
+ | 
 | 
| 545 | 
  | 
    //loop over all of the atom types | 
| 546 | 
  | 
    for (i = atomTypes.begin(); i != atomTypes.end(); ++i) { | 
| 547 | 
  | 
      useLennardJones |= (*i)->isLennardJones(); | 
| 603 | 
  | 
    temp = useFLARB; | 
| 604 | 
  | 
    MPI_Allreduce(&temp, &useFLARB, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);     | 
| 605 | 
  | 
 | 
| 606 | 
+ | 
    temp = useRF; | 
| 607 | 
+ | 
    MPI_Allreduce(&temp, &useRF, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);     | 
| 608 | 
+ | 
 | 
| 609 | 
  | 
#endif | 
| 610 | 
  | 
 | 
| 611 | 
  | 
    fInfo_.SIM_uses_PBC = usePBC;     | 
| 620 | 
  | 
    fInfo_.SIM_uses_EAM = useEAM; | 
| 621 | 
  | 
    fInfo_.SIM_uses_Shapes = useShape; | 
| 622 | 
  | 
    fInfo_.SIM_uses_FLARB = useFLARB; | 
| 623 | 
+ | 
    fInfo_.SIM_uses_RF = useRF; | 
| 624 | 
  | 
 | 
| 625 | 
< | 
    if( fInfo_.SIM_uses_Dipoles && fInfo_.SIM_uses_RF) { | 
| 626 | 
< | 
 | 
| 625 | 
> | 
    if( myMethod == "REACTION_FIELD") { | 
| 626 | 
> | 
       | 
| 627 | 
  | 
      if (simParams_->haveDielectric()) { | 
| 628 | 
  | 
        fInfo_.dielect = simParams_->getDielectric(); | 
| 629 | 
  | 
      } else { | 
| 633 | 
  | 
                "\tsetting a dielectric constant!\n"); | 
| 634 | 
  | 
        painCave.isFatal = 1; | 
| 635 | 
  | 
        simError(); | 
| 636 | 
< | 
      } | 
| 614 | 
< | 
         | 
| 615 | 
< | 
    } else { | 
| 616 | 
< | 
      fInfo_.dielect = 0.0; | 
| 636 | 
> | 
      }       | 
| 637 | 
  | 
    } | 
| 618 | 
– | 
 | 
| 638 | 
  | 
  } | 
| 639 | 
  | 
 | 
| 640 | 
  | 
  void SimInfo::setupFortranSim() { | 
| 668 | 
  | 
 | 
| 669 | 
  | 
        totalMass = cg->getMass(); | 
| 670 | 
  | 
        for(atom = cg->beginAtom(ai); atom != NULL; atom = cg->nextAtom(ai)) { | 
| 671 | 
< | 
          mfact.push_back(atom->getMass()/totalMass); | 
| 671 | 
> | 
          // Check for massless groups - set mfact to 1 if true | 
| 672 | 
> | 
          if (totalMass != 0) | 
| 673 | 
> | 
            mfact.push_back(atom->getMass()/totalMass); | 
| 674 | 
> | 
          else | 
| 675 | 
> | 
            mfact.push_back( 1.0 ); | 
| 676 | 
  | 
        } | 
| 677 | 
  | 
 | 
| 678 | 
  | 
      }        | 
| 808 | 
  | 
     | 
| 809 | 
  | 
    if (fInfo_.SIM_uses_Charges | fInfo_.SIM_uses_Dipoles | fInfo_.SIM_uses_RF) { | 
| 810 | 
  | 
         | 
| 811 | 
< | 
      if (!simParams_->haveRcut()){ | 
| 811 | 
> | 
      if (!simParams_->haveCutoffRadius()){ | 
| 812 | 
  | 
        sprintf(painCave.errMsg, | 
| 813 | 
  | 
                "SimCreator Warning: No value was set for the cutoffRadius.\n" | 
| 814 | 
  | 
                "\tOOPSE will use a default value of 15.0 angstroms" | 
| 817 | 
  | 
        simError(); | 
| 818 | 
  | 
        rcut = 15.0; | 
| 819 | 
  | 
      } else{ | 
| 820 | 
< | 
        rcut = simParams_->getRcut(); | 
| 820 | 
> | 
        rcut = simParams_->getCutoffRadius(); | 
| 821 | 
  | 
      } | 
| 822 | 
  | 
 | 
| 823 | 
< | 
      if (!simParams_->haveRsw()){ | 
| 823 | 
> | 
      if (!simParams_->haveSwitchingRadius()){ | 
| 824 | 
  | 
        sprintf(painCave.errMsg, | 
| 825 | 
  | 
                "SimCreator Warning: No value was set for switchingRadius.\n" | 
| 826 | 
  | 
                "\tOOPSE will use a default value of\n" | 
| 827 | 
< | 
                "\t0.95 * cutoffRadius for the switchingRadius\n"); | 
| 827 | 
> | 
                "\t0.85 * cutoffRadius for the switchingRadius\n"); | 
| 828 | 
  | 
        painCave.isFatal = 0; | 
| 829 | 
  | 
        simError(); | 
| 830 | 
< | 
        rsw = 0.95 * rcut; | 
| 830 | 
> | 
        rsw = 0.85 * rcut; | 
| 831 | 
  | 
      } else{ | 
| 832 | 
< | 
        rsw = simParams_->getRsw(); | 
| 832 | 
> | 
        rsw = simParams_->getSwitchingRadius(); | 
| 833 | 
  | 
      } | 
| 834 | 
  | 
 | 
| 835 | 
  | 
    } else { | 
| 836 | 
  | 
      // if charge, dipole or reaction field is not used and the cutofff radius is not specified in | 
| 837 | 
  | 
      //meta-data file, the maximum cutoff radius calculated from forcefiled will be used | 
| 838 | 
  | 
         | 
| 839 | 
< | 
      if (simParams_->haveRcut()) { | 
| 840 | 
< | 
        rcut = simParams_->getRcut(); | 
| 839 | 
> | 
      if (simParams_->haveCutoffRadius()) { | 
| 840 | 
> | 
        rcut = simParams_->getCutoffRadius(); | 
| 841 | 
  | 
      } else { | 
| 842 | 
  | 
        //set cutoff radius to the maximum cutoff radius based on atom types in the whole system | 
| 843 | 
  | 
        rcut = calcMaxCutoffRadius(); | 
| 844 | 
  | 
      } | 
| 845 | 
  | 
 | 
| 846 | 
< | 
      if (simParams_->haveRsw()) { | 
| 847 | 
< | 
        rsw  = simParams_->getRsw(); | 
| 846 | 
> | 
      if (simParams_->haveSwitchingRadius()) { | 
| 847 | 
> | 
        rsw  = simParams_->getSwitchingRadius(); | 
| 848 | 
  | 
      } else { | 
| 849 | 
  | 
        rsw = rcut; | 
| 850 | 
  | 
      } | 
| 861 | 
  | 
    int cp =  TRADITIONAL_CUTOFF_POLICY; | 
| 862 | 
  | 
    if (simParams_->haveCutoffPolicy()) { | 
| 863 | 
  | 
      std::string myPolicy = simParams_->getCutoffPolicy(); | 
| 864 | 
+ | 
      toUpper(myPolicy); | 
| 865 | 
  | 
      if (myPolicy == "MIX") { | 
| 866 | 
  | 
        cp = MIX_CUTOFF_POLICY; | 
| 867 | 
  | 
      } else { | 
| 880 | 
  | 
        }            | 
| 881 | 
  | 
      } | 
| 882 | 
  | 
    } | 
| 883 | 
+ | 
 | 
| 884 | 
+ | 
 | 
| 885 | 
+ | 
    if (simParams_->haveSkinThickness()) { | 
| 886 | 
+ | 
      double skinThickness = simParams_->getSkinThickness(); | 
| 887 | 
+ | 
    } | 
| 888 | 
+ | 
 | 
| 889 | 
  | 
    notifyFortranCutoffs(&rcut_, &rsw_, &rnblist, &cp); | 
| 890 | 
+ | 
    // also send cutoff notification to electrostatics | 
| 891 | 
+ | 
    setElectrostaticCutoffRadius(&rcut_, &rsw_); | 
| 892 | 
  | 
  } | 
| 893 | 
  | 
 | 
| 894 | 
  | 
  void SimInfo::setupElectrostaticSummationMethod( int isError ) {     | 
| 896 | 
  | 
    int errorOut; | 
| 897 | 
  | 
    int esm =  NONE; | 
| 898 | 
  | 
    double alphaVal; | 
| 899 | 
+ | 
    double dielectric; | 
| 900 | 
  | 
 | 
| 901 | 
  | 
    errorOut = isError; | 
| 902 | 
+ | 
    alphaVal = simParams_->getDampingAlpha(); | 
| 903 | 
+ | 
    dielectric = simParams_->getDielectric(); | 
| 904 | 
  | 
 | 
| 905 | 
  | 
    if (simParams_->haveElectrostaticSummationMethod()) { | 
| 906 | 
  | 
      std::string myMethod = simParams_->getElectrostaticSummationMethod(); | 
| 907 | 
+ | 
      toUpper(myMethod); | 
| 908 | 
  | 
      if (myMethod == "NONE") { | 
| 909 | 
  | 
        esm = NONE; | 
| 910 | 
  | 
      } else { | 
| 916 | 
  | 
            if (!simParams_->haveDampingAlpha()) { | 
| 917 | 
  | 
              //throw error | 
| 918 | 
  | 
              sprintf( painCave.errMsg, | 
| 919 | 
< | 
                       "SimInfo warning: dampingAlpha was not specified in the input file. A default value of %f (1/ang) will be used for the Damped Wolf Method.", simParams_->getDampingAlpha()); | 
| 919 | 
> | 
                       "SimInfo warning: dampingAlpha was not specified in the input file. A default value of %f (1/ang) will be used for the Damped Wolf Method.", alphaVal); | 
| 920 | 
  | 
              painCave.isFatal = 0; | 
| 921 | 
  | 
              simError(); | 
| 922 | 
  | 
            } | 
| 887 | 
– | 
            alphaVal = simParams_->getDampingAlpha(); | 
| 923 | 
  | 
          } else { | 
| 924 | 
< | 
            if (myMethod == "REACTION_FIELD") { | 
| 924 | 
> | 
            if (myMethod == "REACTION_FIELD") {        | 
| 925 | 
  | 
              esm = REACTION_FIELD; | 
| 926 | 
  | 
            } else { | 
| 927 | 
  | 
              // throw error         | 
| 934 | 
  | 
        } | 
| 935 | 
  | 
      } | 
| 936 | 
  | 
    } | 
| 937 | 
< | 
    initFortranFF( &esm, &alphaVal, &errorOut ); | 
| 937 | 
> | 
    // let's pass some summation method variables to fortran | 
| 938 | 
> | 
    setElectrostaticSummationMethod( &esm ); | 
| 939 | 
> | 
    setDampedWolfAlpha( &alphaVal ); | 
| 940 | 
> | 
    setReactionFieldDielectric( &dielectric ); | 
| 941 | 
> | 
    initFortranFF( &esm, &errorOut ); | 
| 942 | 
  | 
  } | 
| 943 | 
  | 
 | 
| 944 | 
  | 
  void SimInfo::addProperty(GenericData* genData) { |