ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/branches/devel_omp/src/parallel/ForceMatrixDecomposition.cpp
(Generate patch)

Comparing branches/devel_omp/src/parallel/ForceMatrixDecomposition.cpp (file contents):
Revision 1599 by mciznick, Fri Jul 29 19:03:36 2011 UTC vs.
Revision 1614 by mciznick, Tue Aug 23 20:55:51 2011 UTC

# Line 1001 | Line 1001 | void ForceMatrixDecomposition::fillInteractionData(Int
1001  
1002   // filling interaction blocks with pointers
1003   void ForceMatrixDecomposition::fillInteractionData(InteractionData &idat, int atom1, int atom2) {
1004 +
1005 +        idat.excluded = excludeAtomPair(atom1, atom2);
1006 +
1007 + #ifdef IS_MPI
1008 +        idat.atypes = make_pair( atypesRow[atom1], atypesCol[atom2]);
1009 +        //idat.atypes = make_pair( ff_->getAtomType(identsRow[atom1]),
1010 +        //                         ff_->getAtomType(identsCol[atom2]) );
1011 +
1012 +        if (storageLayout_ & DataStorage::dslAmat)
1013 +        {
1014 +                idat.A1 = &(atomRowData.aMat[atom1]);
1015 +                idat.A2 = &(atomColData.aMat[atom2]);
1016 +        }
1017 +
1018 +        if (storageLayout_ & DataStorage::dslElectroFrame)
1019 +        {
1020 +                idat.eFrame1 = &(atomRowData.electroFrame[atom1]);
1021 +                idat.eFrame2 = &(atomColData.electroFrame[atom2]);
1022 +        }
1023 +
1024 +        if (storageLayout_ & DataStorage::dslTorque)
1025 +        {
1026 +                idat.t1 = &(atomRowData.torque[atom1]);
1027 +                idat.t2 = &(atomColData.torque[atom2]);
1028 +        }
1029 +
1030 +        if (storageLayout_ & DataStorage::dslDensity)
1031 +        {
1032 +                idat.rho1 = &(atomRowData.density[atom1]);
1033 +                idat.rho2 = &(atomColData.density[atom2]);
1034 +        }
1035 +
1036 +        if (storageLayout_ & DataStorage::dslFunctional)
1037 +        {
1038 +                idat.frho1 = &(atomRowData.functional[atom1]);
1039 +                idat.frho2 = &(atomColData.functional[atom2]);
1040 +        }
1041 +
1042 +        if (storageLayout_ & DataStorage::dslFunctionalDerivative)
1043 +        {
1044 +                idat.dfrho1 = &(atomRowData.functionalDerivative[atom1]);
1045 +                idat.dfrho2 = &(atomColData.functionalDerivative[atom2]);
1046 +        }
1047 +
1048 +        if (storageLayout_ & DataStorage::dslParticlePot)
1049 +        {
1050 +                idat.particlePot1 = &(atomRowData.particlePot[atom1]);
1051 +                idat.particlePot2 = &(atomColData.particlePot[atom2]);
1052 +        }
1053 +
1054 +        if (storageLayout_ & DataStorage::dslSkippedCharge)
1055 +        {
1056 +                idat.skippedCharge1 = &(atomRowData.skippedCharge[atom1]);
1057 +                idat.skippedCharge2 = &(atomColData.skippedCharge[atom2]);
1058 +        }
1059 +
1060 + #else
1061 +
1062 +        idat.atypes = make_pair(atypesLocal[atom1], atypesLocal[atom2]);
1063 +        //idat.atypes = make_pair( ff_->getAtomType(idents[atom1]),
1064 +        //                         ff_->getAtomType(idents[atom2]) );
1065 +
1066 +        if (storageLayout_ & DataStorage::dslAmat)
1067 +        {
1068 +                idat.A1 = &(snap_->atomData.aMat[atom1]);
1069 +                idat.A2 = &(snap_->atomData.aMat[atom2]);
1070 +        }
1071 +
1072 +        if (storageLayout_ & DataStorage::dslElectroFrame)
1073 +        {
1074 +                idat.eFrame1 = &(snap_->atomData.electroFrame[atom1]);
1075 +                idat.eFrame2 = &(snap_->atomData.electroFrame[atom2]);
1076 +        }
1077 +
1078 +        if (storageLayout_ & DataStorage::dslTorque)
1079 +        {
1080 +                idat.t1 = &(snap_->atomData.torque[atom1]);
1081 +                idat.t2 = &(snap_->atomData.torque[atom2]);
1082 +        }
1083  
1084 +        if (storageLayout_ & DataStorage::dslDensity)
1085 +        {
1086 +                idat.rho1 = &(snap_->atomData.density[atom1]);
1087 +                idat.rho2 = &(snap_->atomData.density[atom2]);
1088 +        }
1089 +
1090 +        if (storageLayout_ & DataStorage::dslFunctional)
1091 +        {
1092 +                idat.frho1 = &(snap_->atomData.functional[atom1]);
1093 +                idat.frho2 = &(snap_->atomData.functional[atom2]);
1094 +        }
1095 +
1096 +        if (storageLayout_ & DataStorage::dslFunctionalDerivative)
1097 +        {
1098 +                idat.dfrho1 = &(snap_->atomData.functionalDerivative[atom1]);
1099 +                idat.dfrho2 = &(snap_->atomData.functionalDerivative[atom2]);
1100 +        }
1101 +
1102 +        if (storageLayout_ & DataStorage::dslParticlePot)
1103 +        {
1104 +                idat.particlePot1 = &(snap_->atomData.particlePot[atom1]);
1105 +                idat.particlePot2 = &(snap_->atomData.particlePot[atom2]);
1106 +        }
1107 +
1108 +        if (storageLayout_ & DataStorage::dslSkippedCharge)
1109 +        {
1110 +                idat.skippedCharge1 = &(snap_->atomData.skippedCharge[atom1]);
1111 +                idat.skippedCharge2 = &(snap_->atomData.skippedCharge[atom2]);
1112 +        }
1113 + #endif
1114 + }
1115 +
1116 + // filling interaction blocks with pointers
1117 + void ForceMatrixDecomposition::fillInteractionDataOMP(InteractionDataPrv &idat, int atom1, int atom2) {
1118 +
1119          idat.excluded = excludeAtomPair(atom1, atom2);
1120  
1121   #ifdef IS_MPI
# Line 1126 | Line 1240 | void ForceMatrixDecomposition::unpackInteractionData(I
1240          snap_->atomData.force[atom1] += *(idat.f1);
1241          snap_->atomData.force[atom2] -= *(idat.f1);
1242   #endif
1243 +
1244 + }
1245 +
1246 + void ForceMatrixDecomposition::unpackInteractionDataOMP(InteractionDataPrv &idat, int atom1, int atom2) {
1247 +                pairwisePot += idat.pot;
1248  
1249 +                snap_->atomData.force[atom1] += idat.f1;
1250 +                snap_->atomData.force[atom2] -= idat.f1;
1251   }
1252  
1253   void ForceMatrixDecomposition::reorderGroupCutoffs(vector<int> &order) {
# Line 1172 | Line 1293 | void ForceMatrixDecomposition::reorderPosition(vector<
1293                          cg->setLocalIndex(mapPos[cg->getLocalIndex()]);
1294                  }
1295          }
1175
1176        /*      if (info_->getNCutoffGroups() > 0)
1177         {
1178         for (mol = info_->beginMolecule(mi); mol != NULL; mol = info_->nextMolecule(mi))
1179         {
1180         for (cg = mol->beginCutoffGroup(ci); cg != NULL; cg = mol->nextCutoffGroup(ci))
1181         {
1182         printf("gbI:%d locI:%d x:%f y:%f z:%f\n", cg->getGlobalIndex(), cg->getLocalIndex(),
1183         cgConfig->position[cg->getLocalIndex()].x(), cgConfig->position[cg->getLocalIndex()].y(),
1184         cgConfig->position[cg->getLocalIndex()].z());
1185         }
1186         }
1187         } else
1188         {
1189         // center of mass of the group is the same as position of the atom
1190         // if cutoff group does not exist
1191         printf("ERROR!!!!!!!!!!!!!!!!!!!!!!!!!!!\n");
1192         //                     cgConfig->position = config->position;
1193         }*/
1296   }
1297  
1298   void ForceMatrixDecomposition::reorderGroupList(vector<int> &order) {
# Line 1209 | Line 1311 | void ForceMatrixDecomposition::reorderMemory(vector<ve
1311  
1312   void ForceMatrixDecomposition::reorderMemory(vector<vector<CutoffGroup *> > &H_c_l) {
1313          int n = 0;
1212        //      printf("Reorder memory time:%f!!!!!!!!!!!!!!!!!!!!!!!!!!!\n",
1213        //                      info_->getSnapshotManager()->getCurrentSnapshot()->getTime());
1314  
1315          /* record the reordered atom indices */
1316          vector<int> k = vector<int> (nGroups_);
# Line 1231 | Line 1331 | vector<vector<CutoffGroup *> > ForceMatrixDecompositio
1331   }
1332  
1333   vector<vector<CutoffGroup *> > ForceMatrixDecomposition::buildLayerBasedNeighborList() {
1234        //      printf("buildLayerBasedNeighborList; nGroups:%d\n", nGroups_);
1334          // Na = nGroups_
1335          /* cell occupancy counter */
1336          //      vector<int> k_c;
# Line 1323 | Line 1422 | vector<vector<CutoffGroup *> > ForceMatrixDecompositio
1422          {
1423                  if (reorderFreqCounter == neighborListReorderFreq)
1424                  {
1326                        //printf("neighborListReorderFreq:%d\n", neighborListReorderFreq);
1425                          reorderMemory(H_c_l);
1426                          reorderFreqCounter = 1;
1427                  } else
# Line 1333 | Line 1431 | vector<vector<CutoffGroup *> > ForceMatrixDecompositio
1431          }
1432  
1433          int m;
1434 <        /* the neighbor matrix */
1434 >        /* The neighbor matrix */
1435          vector<vector<CutoffGroup *> > neighborMatW = vector<vector<CutoffGroup *> > (nGroups_);
1436  
1437          groupCutoffs cuts;
1438          CutoffGroup *cg1;
1439  
1440 <        /* loops over objects(atoms, rigidBodies, cutoffGroups, etc.) */
1440 >        /* Loops over objects(atoms, rigidBodies, cutoffGroups, etc.) */
1441          for (mol = info_->beginMolecule(mi); mol != NULL; mol = info_->nextMolecule(mi))
1442          {
1443                  for (cg1 = mol->beginCutoffGroup(ci); cg1 != NULL; cg1 = mol->nextCutoffGroup(ci))
# Line 1378 | Line 1476 | vector<vector<CutoffGroup *> > ForceMatrixDecompositio
1476                                  }
1477  
1478                                  int c2 = Vlinear(c2v, nCells_);
1479 <                                /* loops over layers l to access the neighbor atoms */
1479 >                                /* Loops over layers l to access the neighbor atoms */
1480                                  for (vector<CutoffGroup *>::iterator cg2 = H_c_l[c2].begin(); cg2 != H_c_l[c2].end(); ++cg2)
1481                                  {
1482 <                                        //                              if i'' = 0 then break // doesn't apply to vector implementation of matrix
1482 >                                        //                              if i'' = 0 then break // doesn't apply to vector implementation of the matrix
1483                                          //                              if(i != *j)
1484                                          if (c2 != c1 || (*cg2)->getGlobalIndex() < cg1->getGlobalIndex())
1485                                          {
# Line 1390 | Line 1488 | vector<vector<CutoffGroup *> > ForceMatrixDecompositio
1488                                                  cuts = getGroupCutoffs(cg1->getGlobalIndex(), (*cg2)->getGlobalIndex());
1489                                                  if (dr.lengthSquare() < cuts.third)
1490                                                  {
1491 <                                                        /* transposed version of Rapaport W mat, to occupy successive memory locations on CPU */
1491 >                                                        /* Transposed version of Rapaport W mat, to occupy successive memory locations on CPU */
1492                                                          neighborMatW[cg1->getGlobalIndex()].push_back((*cg2));
1493                                                  }
1494                                          }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines