| 54 | 
  | 
#include "math/Vector3.hpp" | 
| 55 | 
  | 
#include "primitives/Molecule.hpp" | 
| 56 | 
  | 
#include "primitives/StuntDouble.hpp" | 
| 57 | 
– | 
#include "UseTheForce/DarkSide/neighborLists_interface.h" | 
| 57 | 
  | 
#include "utils/MemoryUtils.hpp" | 
| 58 | 
  | 
#include "utils/simError.h" | 
| 59 | 
  | 
#include "selection/SelectionManager.hpp" | 
| 61 | 
  | 
#include "UseTheForce/ForceField.hpp" | 
| 62 | 
  | 
#include "nonbonded/SwitchingFunction.hpp" | 
| 63 | 
  | 
 | 
| 65 | 
– | 
#ifdef IS_MPI | 
| 66 | 
– | 
#include "UseTheForce/mpiComponentPlan.h" | 
| 67 | 
– | 
#include "UseTheForce/DarkSide/simParallel_interface.h" | 
| 68 | 
– | 
#endif  | 
| 69 | 
– | 
 | 
| 64 | 
  | 
using namespace std; | 
| 65 | 
  | 
namespace OpenMD { | 
| 66 | 
  | 
   | 
| 71 | 
  | 
    nGlobalIntegrableObjects_(0), nGlobalRigidBodies_(0), | 
| 72 | 
  | 
    nAtoms_(0), nBonds_(0),  nBends_(0), nTorsions_(0), nInversions_(0),  | 
| 73 | 
  | 
    nRigidBodies_(0), nIntegrableObjects_(0), nCutoffGroups_(0),  | 
| 74 | 
< | 
    nConstraints_(0), sman_(NULL), fortranInitialized_(false),  | 
| 74 | 
> | 
    nConstraints_(0), sman_(NULL), topologyDone_(false),  | 
| 75 | 
  | 
    calcBoxDipole_(false), useAtomicVirial_(true) {     | 
| 76 | 
  | 
     | 
| 77 | 
  | 
    MoleculeStamp* molStamp; | 
| 125 | 
  | 
    //equal to the total number of atoms minus number of atoms belong to  | 
| 126 | 
  | 
    //cutoff group defined in meta-data file plus the number of cutoff  | 
| 127 | 
  | 
    //groups defined in meta-data file | 
| 128 | 
+ | 
 | 
| 129 | 
  | 
    nGlobalCutoffGroups_ = nGlobalAtoms_ - nCutoffAtoms + nGroups; | 
| 130 | 
  | 
     | 
| 131 | 
  | 
    //every free atom (atom does not belong to rigid bodies) is an  | 
| 269 | 
  | 
#endif | 
| 270 | 
  | 
    return fdf_; | 
| 271 | 
  | 
  } | 
| 272 | 
+ | 
   | 
| 273 | 
+ | 
  unsigned int SimInfo::getNLocalCutoffGroups(){ | 
| 274 | 
+ | 
    int nLocalCutoffAtoms = 0; | 
| 275 | 
+ | 
    Molecule* mol; | 
| 276 | 
+ | 
    MoleculeIterator mi; | 
| 277 | 
+ | 
    CutoffGroup* cg; | 
| 278 | 
+ | 
    Molecule::CutoffGroupIterator ci; | 
| 279 | 
+ | 
     | 
| 280 | 
+ | 
    for (mol = beginMolecule(mi); mol != NULL; mol  = nextMolecule(mi)) { | 
| 281 | 
+ | 
       | 
| 282 | 
+ | 
      for (cg = mol->beginCutoffGroup(ci); cg != NULL;  | 
| 283 | 
+ | 
           cg = mol->nextCutoffGroup(ci)) { | 
| 284 | 
+ | 
        nLocalCutoffAtoms += cg->getNumAtom(); | 
| 285 | 
+ | 
         | 
| 286 | 
+ | 
      }         | 
| 287 | 
+ | 
    } | 
| 288 | 
+ | 
     | 
| 289 | 
+ | 
    return nAtoms_ - nLocalCutoffAtoms + nCutoffGroups_; | 
| 290 | 
+ | 
  } | 
| 291 | 
  | 
     | 
| 292 | 
  | 
  void SimInfo::calcNdfRaw() { | 
| 293 | 
  | 
    int ndfRaw_local; | 
| 784 | 
  | 
    temp = usesElectrostatic; | 
| 785 | 
  | 
    MPI_Allreduce(&temp, &usesElectrostaticAtoms_, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);  | 
| 786 | 
  | 
#endif | 
| 773 | 
– | 
    fInfo_.SIM_uses_PBC = usesPeriodicBoundaries_;     | 
| 774 | 
– | 
    fInfo_.SIM_uses_DirectionalAtoms = usesDirectionalAtoms_; | 
| 775 | 
– | 
    fInfo_.SIM_uses_MetallicAtoms = usesMetallicAtoms_; | 
| 776 | 
– | 
    fInfo_.SIM_requires_SkipCorrection = usesElectrostaticAtoms_; | 
| 777 | 
– | 
    fInfo_.SIM_requires_SelfCorrection = usesElectrostaticAtoms_; | 
| 778 | 
– | 
    fInfo_.SIM_uses_AtomicVirial = usesAtomicVirial_; | 
| 787 | 
  | 
  } | 
| 788 | 
  | 
 | 
| 781 | 
– | 
  void SimInfo::setupFortran() { | 
| 782 | 
– | 
    int isError; | 
| 783 | 
– | 
    int nExclude, nOneTwo, nOneThree, nOneFour; | 
| 784 | 
– | 
    vector<int> fortranGlobalGroupMembership; | 
| 785 | 
– | 
     | 
| 786 | 
– | 
    isError = 0; | 
| 789 | 
  | 
 | 
| 790 | 
< | 
    //globalGroupMembership_ is filled by SimCreator     | 
| 791 | 
< | 
    for (int i = 0; i < nGlobalAtoms_; i++) { | 
| 792 | 
< | 
      fortranGlobalGroupMembership.push_back(globalGroupMembership_[i] + 1); | 
| 790 | 
> | 
  vector<int> SimInfo::getGlobalAtomIndices() { | 
| 791 | 
> | 
    SimInfo::MoleculeIterator mi; | 
| 792 | 
> | 
    Molecule* mol; | 
| 793 | 
> | 
    Molecule::AtomIterator ai; | 
| 794 | 
> | 
    Atom* atom; | 
| 795 | 
> | 
 | 
| 796 | 
> | 
    vector<int> GlobalAtomIndices(getNAtoms(), 0); | 
| 797 | 
> | 
     | 
| 798 | 
> | 
    for (mol = beginMolecule(mi); mol != NULL; mol  = nextMolecule(mi)) { | 
| 799 | 
> | 
       | 
| 800 | 
> | 
      for (atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) { | 
| 801 | 
> | 
        GlobalAtomIndices[atom->getLocalIndex()] = atom->getGlobalIndex(); | 
| 802 | 
> | 
      } | 
| 803 | 
  | 
    } | 
| 804 | 
+ | 
    return GlobalAtomIndices; | 
| 805 | 
+ | 
  } | 
| 806 | 
  | 
 | 
| 807 | 
+ | 
 | 
| 808 | 
+ | 
  vector<int> SimInfo::getGlobalGroupIndices() { | 
| 809 | 
+ | 
    SimInfo::MoleculeIterator mi; | 
| 810 | 
+ | 
    Molecule* mol; | 
| 811 | 
+ | 
    Molecule::CutoffGroupIterator ci; | 
| 812 | 
+ | 
    CutoffGroup* cg; | 
| 813 | 
+ | 
 | 
| 814 | 
+ | 
    vector<int> GlobalGroupIndices; | 
| 815 | 
+ | 
     | 
| 816 | 
+ | 
    for (mol = beginMolecule(mi); mol != NULL; mol  = nextMolecule(mi)) { | 
| 817 | 
+ | 
       | 
| 818 | 
+ | 
      //local index of cutoff group is trivial, it only depends on the | 
| 819 | 
+ | 
      //order of travesing | 
| 820 | 
+ | 
      for (cg = mol->beginCutoffGroup(ci); cg != NULL;  | 
| 821 | 
+ | 
           cg = mol->nextCutoffGroup(ci)) { | 
| 822 | 
+ | 
        GlobalGroupIndices.push_back(cg->getGlobalIndex()); | 
| 823 | 
+ | 
      }         | 
| 824 | 
+ | 
    } | 
| 825 | 
+ | 
    return GlobalGroupIndices; | 
| 826 | 
+ | 
  } | 
| 827 | 
+ | 
 | 
| 828 | 
+ | 
 | 
| 829 | 
+ | 
  void SimInfo::prepareTopology() { | 
| 830 | 
+ | 
    int nExclude, nOneTwo, nOneThree, nOneFour; | 
| 831 | 
+ | 
 | 
| 832 | 
  | 
    //calculate mass ratio of cutoff group | 
| 794 | 
– | 
    vector<RealType> mfact; | 
| 833 | 
  | 
    SimInfo::MoleculeIterator mi; | 
| 834 | 
  | 
    Molecule* mol; | 
| 835 | 
  | 
    Molecule::CutoffGroupIterator ci; | 
| 838 | 
  | 
    Atom* atom; | 
| 839 | 
  | 
    RealType totalMass; | 
| 840 | 
  | 
 | 
| 841 | 
< | 
    //to avoid memory reallocation, reserve enough space for mfact | 
| 842 | 
< | 
    mfact.reserve(getNCutoffGroups()); | 
| 841 | 
> | 
    //to avoid memory reallocation, reserve enough space for massFactors_ | 
| 842 | 
> | 
    massFactors_.clear(); | 
| 843 | 
> | 
    massFactors_.reserve(getNCutoffGroups()); | 
| 844 | 
  | 
     | 
| 845 | 
  | 
    for(mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) {         | 
| 846 | 
< | 
      for (cg = mol->beginCutoffGroup(ci); cg != NULL; cg = mol->nextCutoffGroup(ci)) { | 
| 846 | 
> | 
      for (cg = mol->beginCutoffGroup(ci); cg != NULL;  | 
| 847 | 
> | 
           cg = mol->nextCutoffGroup(ci)) { | 
| 848 | 
  | 
 | 
| 849 | 
  | 
        totalMass = cg->getMass(); | 
| 850 | 
  | 
        for(atom = cg->beginAtom(ai); atom != NULL; atom = cg->nextAtom(ai)) { | 
| 851 | 
  | 
          // Check for massless groups - set mfact to 1 if true | 
| 852 | 
  | 
          if (totalMass != 0) | 
| 853 | 
< | 
            mfact.push_back(atom->getMass()/totalMass); | 
| 853 | 
> | 
            massFactors_.push_back(atom->getMass()/totalMass); | 
| 854 | 
  | 
          else | 
| 855 | 
< | 
            mfact.push_back( 1.0 ); | 
| 855 | 
> | 
            massFactors_.push_back( 1.0 ); | 
| 856 | 
  | 
        } | 
| 857 | 
  | 
      }        | 
| 858 | 
  | 
    } | 
| 859 | 
  | 
 | 
| 860 | 
< | 
    //fill ident array of local atoms (it is actually ident of | 
| 821 | 
< | 
    //AtomType, it is so confusing !!!) | 
| 822 | 
< | 
    vector<int> identArray; | 
| 860 | 
> | 
    // Build the identArray_ | 
| 861 | 
  | 
 | 
| 862 | 
< | 
    //to avoid memory reallocation, reserve enough space identArray | 
| 863 | 
< | 
    identArray.reserve(getNAtoms()); | 
| 826 | 
< | 
     | 
| 862 | 
> | 
    identArray_.clear(); | 
| 863 | 
> | 
    identArray_.reserve(getNAtoms());     | 
| 864 | 
  | 
    for(mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) {         | 
| 865 | 
  | 
      for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) { | 
| 866 | 
< | 
        identArray.push_back(atom->getIdent()); | 
| 866 | 
> | 
        identArray_.push_back(atom->getIdent()); | 
| 867 | 
  | 
      } | 
| 868 | 
  | 
    }     | 
| 832 | 
– | 
 | 
| 833 | 
– | 
    //fill molMembershipArray | 
| 834 | 
– | 
    //molMembershipArray is filled by SimCreator     | 
| 835 | 
– | 
    vector<int> molMembershipArray(nGlobalAtoms_); | 
| 836 | 
– | 
    for (int i = 0; i < nGlobalAtoms_; i++) { | 
| 837 | 
– | 
      molMembershipArray[i] = globalMolMembership_[i] + 1; | 
| 838 | 
– | 
    } | 
| 869 | 
  | 
     | 
| 870 | 
< | 
    //setup fortran simulation | 
| 870 | 
> | 
    //scan topology  | 
| 871 | 
  | 
 | 
| 872 | 
  | 
    nExclude = excludedInteractions_.getSize(); | 
| 873 | 
  | 
    nOneTwo = oneTwoInteractions_.getSize(); | 
| 879 | 
  | 
    int* oneThreeList = oneThreeInteractions_.getPairList(); | 
| 880 | 
  | 
    int* oneFourList = oneFourInteractions_.getPairList(); | 
| 881 | 
  | 
 | 
| 882 | 
< | 
    setFortranSim( &fInfo_, &nGlobalAtoms_, &nAtoms_, &identArray[0],  | 
| 883 | 
< | 
                   &nExclude, excludeList,  | 
| 884 | 
< | 
                   &nOneTwo, oneTwoList, | 
| 885 | 
< | 
                   &nOneThree, oneThreeList, | 
| 886 | 
< | 
                   &nOneFour, oneFourList, | 
| 887 | 
< | 
                   &molMembershipArray[0], &mfact[0], &nCutoffGroups_,  | 
| 888 | 
< | 
                   &fortranGlobalGroupMembership[0], &isError);  | 
| 882 | 
> | 
    //setFortranSim( &fInfo_, &nGlobalAtoms_, &nAtoms_, &identArray_[0],  | 
| 883 | 
> | 
    //               &nExclude, excludeList,  | 
| 884 | 
> | 
    //               &nOneTwo, oneTwoList, | 
| 885 | 
> | 
    //               &nOneThree, oneThreeList, | 
| 886 | 
> | 
    //               &nOneFour, oneFourList, | 
| 887 | 
> | 
    //               &molMembershipArray[0], &mfact[0], &nCutoffGroups_,  | 
| 888 | 
> | 
    //               &fortranGlobalGroupMembership[0], &isError);  | 
| 889 | 
  | 
     | 
| 890 | 
< | 
    if( isError ){ | 
| 861 | 
< | 
       | 
| 862 | 
< | 
      sprintf( painCave.errMsg, | 
| 863 | 
< | 
               "There was an error setting the simulation information in fortran.\n" ); | 
| 864 | 
< | 
      painCave.isFatal = 1; | 
| 865 | 
< | 
      painCave.severity = OPENMD_ERROR; | 
| 866 | 
< | 
      simError(); | 
| 867 | 
< | 
    } | 
| 868 | 
< | 
     | 
| 869 | 
< | 
     | 
| 870 | 
< | 
    sprintf( checkPointMsg, | 
| 871 | 
< | 
             "succesfully sent the simulation information to fortran.\n"); | 
| 872 | 
< | 
     | 
| 873 | 
< | 
    errorCheckPoint(); | 
| 874 | 
< | 
     | 
| 875 | 
< | 
    // Setup number of neighbors in neighbor list if present | 
| 876 | 
< | 
    if (simParams_->haveNeighborListNeighbors()) { | 
| 877 | 
< | 
      int nlistNeighbors = simParams_->getNeighborListNeighbors(); | 
| 878 | 
< | 
      setNeighbors(&nlistNeighbors); | 
| 879 | 
< | 
    } | 
| 880 | 
< | 
    | 
| 881 | 
< | 
#ifdef IS_MPI     | 
| 882 | 
< | 
    //SimInfo is responsible for creating localToGlobalAtomIndex and | 
| 883 | 
< | 
    //localToGlobalGroupIndex | 
| 884 | 
< | 
    vector<int> localToGlobalAtomIndex(getNAtoms(), 0); | 
| 885 | 
< | 
    vector<int> localToGlobalCutoffGroupIndex; | 
| 886 | 
< | 
    mpiSimData parallelData; | 
| 887 | 
< | 
 | 
| 888 | 
< | 
    for (mol = beginMolecule(mi); mol != NULL; mol  = nextMolecule(mi)) { | 
| 889 | 
< | 
 | 
| 890 | 
< | 
      //local index(index in DataStorge) of atom is important | 
| 891 | 
< | 
      for (atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) { | 
| 892 | 
< | 
        localToGlobalAtomIndex[atom->getLocalIndex()] = atom->getGlobalIndex() + 1; | 
| 893 | 
< | 
      } | 
| 894 | 
< | 
 | 
| 895 | 
< | 
      //local index of cutoff group is trivial, it only depends on the order of travesing | 
| 896 | 
< | 
      for (cg = mol->beginCutoffGroup(ci); cg != NULL; cg = mol->nextCutoffGroup(ci)) { | 
| 897 | 
< | 
        localToGlobalCutoffGroupIndex.push_back(cg->getGlobalIndex() + 1); | 
| 898 | 
< | 
      }         | 
| 899 | 
< | 
         | 
| 900 | 
< | 
    } | 
| 901 | 
< | 
 | 
| 902 | 
< | 
    //fill up mpiSimData struct | 
| 903 | 
< | 
    parallelData.nMolGlobal = getNGlobalMolecules(); | 
| 904 | 
< | 
    parallelData.nMolLocal = getNMolecules(); | 
| 905 | 
< | 
    parallelData.nAtomsGlobal = getNGlobalAtoms(); | 
| 906 | 
< | 
    parallelData.nAtomsLocal = getNAtoms(); | 
| 907 | 
< | 
    parallelData.nGroupsGlobal = getNGlobalCutoffGroups(); | 
| 908 | 
< | 
    parallelData.nGroupsLocal = getNCutoffGroups(); | 
| 909 | 
< | 
    parallelData.myNode = worldRank; | 
| 910 | 
< | 
    MPI_Comm_size(MPI_COMM_WORLD, &(parallelData.nProcessors)); | 
| 911 | 
< | 
 | 
| 912 | 
< | 
    //pass mpiSimData struct and index arrays to fortran | 
| 913 | 
< | 
    setFsimParallel(¶llelData, &(parallelData.nAtomsLocal), | 
| 914 | 
< | 
                    &localToGlobalAtomIndex[0],  &(parallelData.nGroupsLocal), | 
| 915 | 
< | 
                    &localToGlobalCutoffGroupIndex[0], &isError); | 
| 916 | 
< | 
 | 
| 917 | 
< | 
    if (isError) { | 
| 918 | 
< | 
      sprintf(painCave.errMsg, | 
| 919 | 
< | 
              "mpiRefresh errror: fortran didn't like something we gave it.\n"); | 
| 920 | 
< | 
      painCave.isFatal = 1; | 
| 921 | 
< | 
      simError(); | 
| 922 | 
< | 
    } | 
| 923 | 
< | 
 | 
| 924 | 
< | 
    sprintf(checkPointMsg, " mpiRefresh successful.\n"); | 
| 925 | 
< | 
    errorCheckPoint(); | 
| 926 | 
< | 
#endif | 
| 927 | 
< | 
    fortranInitialized_ = true; | 
| 890 | 
> | 
    topologyDone_ = true; | 
| 891 | 
  | 
  } | 
| 892 | 
  | 
 | 
| 893 | 
  | 
  void SimInfo::addProperty(GenericData* genData) { | 
| 924 | 
  | 
    Molecule* mol; | 
| 925 | 
  | 
    RigidBody* rb; | 
| 926 | 
  | 
    Atom* atom; | 
| 927 | 
+ | 
    CutoffGroup* cg; | 
| 928 | 
  | 
    SimInfo::MoleculeIterator mi; | 
| 929 | 
  | 
    Molecule::RigidBodyIterator rbIter; | 
| 930 | 
< | 
    Molecule::AtomIterator atomIter;; | 
| 930 | 
> | 
    Molecule::AtomIterator atomIter; | 
| 931 | 
> | 
    Molecule::CutoffGroupIterator cgIter; | 
| 932 | 
  | 
  | 
| 933 | 
  | 
    for (mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) { | 
| 934 | 
  | 
         | 
| 939 | 
  | 
      for (rb = mol->beginRigidBody(rbIter); rb != NULL; rb = mol->nextRigidBody(rbIter)) { | 
| 940 | 
  | 
        rb->setSnapshotManager(sman_); | 
| 941 | 
  | 
      } | 
| 942 | 
+ | 
 | 
| 943 | 
+ | 
      for (cg = mol->beginCutoffGroup(cgIter); cg != NULL; cg = mol->nextCutoffGroup(cgIter)) { | 
| 944 | 
+ | 
        cg->setSnapshotManager(sman_); | 
| 945 | 
+ | 
      } | 
| 946 | 
  | 
    }     | 
| 947 | 
  | 
     | 
| 948 | 
  | 
  } |