| 53 | 
  | 
#include "brains/SimInfo.hpp" | 
| 54 | 
  | 
#include "math/Vector3.hpp" | 
| 55 | 
  | 
#include "primitives/Molecule.hpp" | 
| 56 | 
+ | 
#include "primitives/StuntDouble.hpp" | 
| 57 | 
  | 
#include "UseTheForce/fCutoffPolicy.h" | 
| 58 | 
  | 
#include "UseTheForce/DarkSide/fElectrostaticSummationMethod.h" | 
| 59 | 
  | 
#include "UseTheForce/DarkSide/fElectrostaticScreeningMethod.h" | 
| 90 | 
  | 
    nGlobalIntegrableObjects_(0), nGlobalRigidBodies_(0), | 
| 91 | 
  | 
    nAtoms_(0), nBonds_(0),  nBends_(0), nTorsions_(0), nRigidBodies_(0), | 
| 92 | 
  | 
    nIntegrableObjects_(0),  nCutoffGroups_(0), nConstraints_(0), | 
| 93 | 
< | 
    sman_(NULL), fortranInitialized_(false) { | 
| 93 | 
> | 
    sman_(NULL), fortranInitialized_(false), calcBoxDipole_(false) { | 
| 94 | 
  | 
 | 
| 95 | 
  | 
      MoleculeStamp* molStamp; | 
| 96 | 
  | 
      int nMolWithSameStamp; | 
| 601 | 
  | 
    /** @deprecate */     | 
| 602 | 
  | 
    int isError = 0; | 
| 603 | 
  | 
     | 
| 604 | 
+ | 
    setupCutoff(); | 
| 605 | 
+ | 
     | 
| 606 | 
  | 
    setupElectrostaticSummationMethod( isError ); | 
| 607 | 
  | 
    setupSwitchingFunction(); | 
| 608 | 
+ | 
    setupAccumulateBoxDipole(); | 
| 609 | 
  | 
 | 
| 610 | 
  | 
    if(isError){ | 
| 611 | 
  | 
      sprintf( painCave.errMsg, | 
| 613 | 
  | 
      painCave.isFatal = 1; | 
| 614 | 
  | 
      simError(); | 
| 615 | 
  | 
    } | 
| 612 | 
– | 
   | 
| 613 | 
– | 
     | 
| 614 | 
– | 
    setupCutoff(); | 
| 616 | 
  | 
 | 
| 617 | 
  | 
    calcNdf(); | 
| 618 | 
  | 
    calcNdfRaw(); | 
| 662 | 
  | 
    int usePBC = simParams_->getUsePeriodicBoundaryConditions(); | 
| 663 | 
  | 
    int useRF; | 
| 664 | 
  | 
    int useSF; | 
| 665 | 
+ | 
    int useSP; | 
| 666 | 
+ | 
    int useBoxDipole; | 
| 667 | 
  | 
    std::string myMethod; | 
| 668 | 
  | 
 | 
| 669 | 
  | 
    // set the useRF logical | 
| 670 | 
  | 
    useRF = 0; | 
| 671 | 
  | 
    useSF = 0; | 
| 672 | 
+ | 
    useSP = 0; | 
| 673 | 
  | 
 | 
| 674 | 
  | 
 | 
| 675 | 
  | 
    if (simParams_->haveElectrostaticSummationMethod()) { | 
| 676 | 
  | 
      std::string myMethod = simParams_->getElectrostaticSummationMethod(); | 
| 677 | 
  | 
      toUpper(myMethod); | 
| 678 | 
< | 
      if (myMethod == "REACTION_FIELD") { | 
| 679 | 
< | 
        useRF=1; | 
| 680 | 
< | 
      } else { | 
| 681 | 
< | 
        if (myMethod == "SHIFTED_FORCE") { | 
| 682 | 
< | 
          useSF = 1; | 
| 683 | 
< | 
        } | 
| 678 | 
> | 
      if (myMethod == "REACTION_FIELD"){ | 
| 679 | 
> | 
        useRF = 1; | 
| 680 | 
> | 
      } else if (myMethod == "SHIFTED_FORCE"){ | 
| 681 | 
> | 
        useSF = 1; | 
| 682 | 
> | 
      } else if (myMethod == "SHIFTED_POTENTIAL"){ | 
| 683 | 
> | 
        useSP = 1; | 
| 684 | 
  | 
      } | 
| 685 | 
  | 
    } | 
| 686 | 
+ | 
     | 
| 687 | 
+ | 
    if (simParams_->haveAccumulateBoxDipole())  | 
| 688 | 
+ | 
      if (simParams_->getAccumulateBoxDipole()) | 
| 689 | 
+ | 
        useBoxDipole = 1; | 
| 690 | 
  | 
 | 
| 691 | 
  | 
    //loop over all of the atom types | 
| 692 | 
  | 
    for (i = atomTypes.begin(); i != atomTypes.end(); ++i) { | 
| 757 | 
  | 
    MPI_Allreduce(&temp, &useRF, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);     | 
| 758 | 
  | 
 | 
| 759 | 
  | 
    temp = useSF; | 
| 760 | 
< | 
    MPI_Allreduce(&temp, &useSF, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);     | 
| 760 | 
> | 
    MPI_Allreduce(&temp, &useSF, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    | 
| 761 | 
  | 
 | 
| 762 | 
+ | 
    temp = useSP; | 
| 763 | 
+ | 
    MPI_Allreduce(&temp, &useSP, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); | 
| 764 | 
+ | 
 | 
| 765 | 
+ | 
    temp = useBoxDipole; | 
| 766 | 
+ | 
    MPI_Allreduce(&temp, &useBoxDipole, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);  | 
| 767 | 
+ | 
 | 
| 768 | 
  | 
#endif | 
| 769 | 
  | 
 | 
| 770 | 
  | 
    fInfo_.SIM_uses_PBC = usePBC;     | 
| 782 | 
  | 
    fInfo_.SIM_uses_FLARB = useFLARB; | 
| 783 | 
  | 
    fInfo_.SIM_uses_RF = useRF; | 
| 784 | 
  | 
    fInfo_.SIM_uses_SF = useSF; | 
| 785 | 
< | 
 | 
| 786 | 
< | 
    if( myMethod == "REACTION_FIELD") { | 
| 773 | 
< | 
       | 
| 774 | 
< | 
      if (simParams_->haveDielectric()) { | 
| 775 | 
< | 
        fInfo_.dielect = simParams_->getDielectric(); | 
| 776 | 
< | 
      } else { | 
| 777 | 
< | 
        sprintf(painCave.errMsg, | 
| 778 | 
< | 
                "SimSetup Error: No Dielectric constant was set.\n" | 
| 779 | 
< | 
                "\tYou are trying to use Reaction Field without" | 
| 780 | 
< | 
                "\tsetting a dielectric constant!\n"); | 
| 781 | 
< | 
        painCave.isFatal = 1; | 
| 782 | 
< | 
        simError(); | 
| 783 | 
< | 
      }       | 
| 784 | 
< | 
    } | 
| 785 | 
< | 
 | 
| 785 | 
> | 
    fInfo_.SIM_uses_SP = useSP; | 
| 786 | 
> | 
    fInfo_.SIM_uses_BoxDipole = useBoxDipole; | 
| 787 | 
  | 
  } | 
| 788 | 
  | 
 | 
| 789 | 
  | 
  void SimInfo::setupFortranSim() { | 
| 1059 | 
  | 
    int sm = UNDAMPED; | 
| 1060 | 
  | 
    RealType alphaVal; | 
| 1061 | 
  | 
    RealType dielectric; | 
| 1062 | 
< | 
 | 
| 1062 | 
> | 
     | 
| 1063 | 
  | 
    errorOut = isError; | 
| 1063 | 
– | 
    alphaVal = simParams_->getDampingAlpha(); | 
| 1064 | 
– | 
    dielectric = simParams_->getDielectric(); | 
| 1064 | 
  | 
 | 
| 1065 | 
  | 
    if (simParams_->haveElectrostaticSummationMethod()) { | 
| 1066 | 
  | 
      std::string myMethod = simParams_->getElectrostaticSummationMethod(); | 
| 1077 | 
  | 
            if (myMethod == "SHIFTED_FORCE") {             | 
| 1078 | 
  | 
              esm = SHIFTED_FORCE; | 
| 1079 | 
  | 
            } else { | 
| 1080 | 
< | 
              if (myMethod == "REACTION_FIELD") {              | 
| 1080 | 
> | 
              if (myMethod == "REACTION_FIELD") { | 
| 1081 | 
  | 
                esm = REACTION_FIELD; | 
| 1082 | 
+ | 
                dielectric = simParams_->getDielectric(); | 
| 1083 | 
+ | 
                if (!simParams_->haveDielectric()) { | 
| 1084 | 
+ | 
                  // throw warning | 
| 1085 | 
+ | 
                  sprintf( painCave.errMsg, | 
| 1086 | 
+ | 
                           "SimInfo warning: dielectric was not specified in the input file\n\tfor the reaction field correction method.\n" | 
| 1087 | 
+ | 
                           "\tA default value of %f will be used for the dielectric.\n", dielectric); | 
| 1088 | 
+ | 
                  painCave.isFatal = 0; | 
| 1089 | 
+ | 
                  simError(); | 
| 1090 | 
+ | 
                } | 
| 1091 | 
  | 
              } else { | 
| 1092 | 
  | 
                // throw error         | 
| 1093 | 
  | 
                sprintf( painCave.errMsg, | 
| 1114 | 
  | 
        if (myScreen == "DAMPED") { | 
| 1115 | 
  | 
          sm = DAMPED; | 
| 1116 | 
  | 
          if (!simParams_->haveDampingAlpha()) { | 
| 1117 | 
< | 
            //throw error | 
| 1117 | 
> | 
            // first set a cutoff dependent alpha value | 
| 1118 | 
> | 
            // we assume alpha depends linearly with rcut from 0 to 20.5 ang | 
| 1119 | 
> | 
            alphaVal = 0.5125 - rcut_* 0.025; | 
| 1120 | 
> | 
            // for values rcut > 20.5, alpha is zero | 
| 1121 | 
> | 
            if (alphaVal < 0) alphaVal = 0; | 
| 1122 | 
> | 
 | 
| 1123 | 
> | 
            // throw warning | 
| 1124 | 
  | 
            sprintf( painCave.errMsg, | 
| 1125 | 
  | 
                     "SimInfo warning: dampingAlpha was not specified in the input file.\n" | 
| 1126 | 
< | 
                     "\tA default value of %f (1/ang) will be used.\n", alphaVal); | 
| 1126 | 
> | 
                     "\tA default value of %f (1/ang) will be used for the cutoff of\n\t%f (ang).\n", alphaVal, rcut_); | 
| 1127 | 
  | 
            painCave.isFatal = 0; | 
| 1128 | 
  | 
            simError(); | 
| 1129 | 
  | 
          } | 
| 1175 | 
  | 
 | 
| 1176 | 
  | 
  } | 
| 1177 | 
  | 
 | 
| 1178 | 
+ | 
  void SimInfo::setupAccumulateBoxDipole() {     | 
| 1179 | 
+ | 
 | 
| 1180 | 
+ | 
    // we only call setAccumulateBoxDipole if the accumulateBoxDipole parameter is true | 
| 1181 | 
+ | 
    if ( simParams_->haveAccumulateBoxDipole() )  | 
| 1182 | 
+ | 
      if ( simParams_->getAccumulateBoxDipole() ) { | 
| 1183 | 
+ | 
        setAccumulateBoxDipole(); | 
| 1184 | 
+ | 
        calcBoxDipole_ = true; | 
| 1185 | 
+ | 
      } | 
| 1186 | 
+ | 
 | 
| 1187 | 
+ | 
  } | 
| 1188 | 
+ | 
 | 
| 1189 | 
  | 
  void SimInfo::addProperty(GenericData* genData) { | 
| 1190 | 
  | 
    properties_.addProperty(genData);   | 
| 1191 | 
  | 
  } | 
| 1435 | 
  | 
      return angularMomentum; | 
| 1436 | 
  | 
   } | 
| 1437 | 
  | 
    | 
| 1438 | 
< | 
    | 
| 1438 | 
> | 
  StuntDouble* SimInfo::getIOIndexToIntegrableObject(int index) { | 
| 1439 | 
> | 
    return IOIndexToIntegrableObject.at(index); | 
| 1440 | 
> | 
  } | 
| 1441 | 
> | 
   | 
| 1442 | 
> | 
  void SimInfo::setIOIndexToIntegrableObject(const std::vector<StuntDouble*>& v) { | 
| 1443 | 
> | 
    IOIndexToIntegrableObject= v; | 
| 1444 | 
> | 
  } | 
| 1445 | 
> | 
 | 
| 1446 | 
> | 
/* | 
| 1447 | 
> | 
   void SimInfo::setStuntDoubleFromGlobalIndex(std::vector<StuntDouble*> v) { | 
| 1448 | 
> | 
      assert( v.size() == nAtoms_ + nRigidBodies_); | 
| 1449 | 
> | 
      sdByGlobalIndex_ = v; | 
| 1450 | 
> | 
    } | 
| 1451 | 
> | 
 | 
| 1452 | 
> | 
    StuntDouble* SimInfo::getStuntDoubleFromGlobalIndex(int index) { | 
| 1453 | 
> | 
      //assert(index < nAtoms_ + nRigidBodies_); | 
| 1454 | 
> | 
      return sdByGlobalIndex_.at(index); | 
| 1455 | 
> | 
    }    | 
| 1456 | 
> | 
*/    | 
| 1457 | 
  | 
}//end namespace oopse | 
| 1458 | 
  | 
 |