# | Line 53 | Line 53 | |
---|---|---|
53 | #include "math/Vector3.hpp" | |
54 | #include "primitives/Molecule.hpp" | |
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" | |
# | Line 81 | Line 83 | namespace oopse { | |
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; | |
# | Line 106 | Line 108 | namespace oopse { | |
108 | } | |
109 | ||
110 | nGroups += nCutoffGroupsInStamp * nMolWithSameStamp; | |
111 | + | |
112 | nCutoffAtoms += nAtomsInGroups * nMolWithSameStamp; | |
113 | ||
114 | //calculate atoms in rigid bodies | |
# | Line 122 | Line 125 | namespace oopse { | |
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; | |
130 | – | |
131 | – | //every free atom (atom does not belong to rigid bodies) is an integrable object |
132 | – | //therefore the total number of integrable objects in the system is equal to |
133 | – | //the total number of atoms minus number of atoms belong to rigid body defined in meta-data |
134 | – | //file plus the number of rigid bodies defined in meta-data file |
135 | – | nGlobalIntegrableObjects_ = nGlobalAtoms_ - nRigidAtoms + nGlobalRigidBodies_; |
134 | ||
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 | |
# | Line 463 | Line 469 | namespace oopse { | |
469 | //setup fortran force field | |
470 | /** @deprecate */ | |
471 | int isError = 0; | |
472 | < | initFortranFF( &fInfo_.SIM_uses_RF, &fInfo_.SIM_uses_UW, |
473 | < | &fInfo_.SIM_uses_DW, &isError ); |
472 | > | |
473 | > | setupElectrostaticSummationMethod( isError ); |
474 | > | |
475 | if(isError){ | |
476 | sprintf( painCave.errMsg, | |
477 | "ForceField error: There was an error initializing the forceField in fortran.\n" ); | |
# | Line 519 | Line 526 | namespace oopse { | |
526 | int useDirectionalAtom = 0; | |
527 | int useElectrostatics = 0; | |
528 | //usePBC and useRF are from simParams | |
529 | < | int usePBC = simParams_->getPBC(); |
530 | < | int useRF = simParams_->getUseRF(); |
531 | < | int useUW = simParams_->getUseUndampedWolf(); |
532 | < | int useDW = simParams_->getUseDampedWolf(); |
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(); | |
# | Line 588 | Line 606 | namespace oopse { | |
606 | temp = useRF; | |
607 | MPI_Allreduce(&temp, &useRF, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); | |
608 | ||
591 | – | temp = useUW; |
592 | – | MPI_Allreduce(&temp, &useUW, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); |
593 | – | |
594 | – | temp = useDW; |
595 | – | MPI_Allreduce(&temp, &useDW, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); |
596 | – | |
609 | #endif | |
610 | ||
611 | fInfo_.SIM_uses_PBC = usePBC; | |
# | Line 609 | Line 621 | namespace oopse { | |
621 | fInfo_.SIM_uses_Shapes = useShape; | |
622 | fInfo_.SIM_uses_FLARB = useFLARB; | |
623 | fInfo_.SIM_uses_RF = useRF; | |
612 | – | fInfo_.SIM_uses_UW = useUW; |
613 | – | fInfo_.SIM_uses_DW = useDW; |
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 { | |
# | Line 623 | Line 633 | namespace oopse { | |
633 | "\tsetting a dielectric constant!\n"); | |
634 | painCave.isFatal = 1; | |
635 | simError(); | |
636 | < | } |
627 | < | |
628 | < | } else { |
629 | < | fInfo_.dielect = 0.0; |
636 | > | } |
637 | } | |
631 | – | |
638 | } | |
639 | ||
640 | void SimInfo::setupFortranSim() { | |
# | Line 662 | Line 668 | namespace oopse { | |
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 | } | |
# | Line 798 | Line 808 | namespace oopse { | |
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" | |
# | Line 807 | Line 817 | namespace oopse { | |
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 | } | |
# | Line 851 | Line 861 | namespace oopse { | |
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 { | |
# | Line 869 | Line 880 | namespace oopse { | |
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 ) { |
895 | + | |
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 { |
911 | + | if (myMethod == "UNDAMPED_WOLF") { |
912 | + | esm = UNDAMPED_WOLF; |
913 | + | } else { |
914 | + | if (myMethod == "DAMPED_WOLF") { |
915 | + | esm = DAMPED_WOLF; |
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.", alphaVal); |
920 | + | painCave.isFatal = 0; |
921 | + | simError(); |
922 | + | } |
923 | + | } else { |
924 | + | if (myMethod == "REACTION_FIELD") { |
925 | + | esm = REACTION_FIELD; |
926 | + | } else { |
927 | + | // throw error |
928 | + | sprintf( painCave.errMsg, |
929 | + | "SimInfo error: Unknown electrostaticSummationMethod. (Input file specified %s .)\n\telectrostaticSummationMethod must be one of: \"none\", \"undamped_wolf\", \"damped_wolf\", or \"reaction_field\".", myMethod.c_str() ); |
930 | + | painCave.isFatal = 1; |
931 | + | simError(); |
932 | + | } |
933 | + | } |
934 | + | } |
935 | + | } |
936 | + | } |
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) { |
– | Removed lines |
+ | Added lines |
< | Changed lines |
> | Changed lines |