| 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 |
| 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) { |
| 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) { |
| 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_); |
| 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; |
| 1422 |
|
{ |
| 1423 |
|
if (reorderFreqCounter == neighborListReorderFreq) |
| 1424 |
|
{ |
| 1326 |
– |
//printf("neighborListReorderFreq:%d\n", neighborListReorderFreq); |
| 1425 |
|
reorderMemory(H_c_l); |
| 1426 |
|
reorderFreqCounter = 1; |
| 1427 |
|
} else |
| 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)) |
| 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 |
|
{ |
| 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 |
|
} |