| 112 |  | AtomCommIntRow->gather(idents, identsRow); | 
| 113 |  | AtomCommIntColumn->gather(idents, identsCol); | 
| 114 |  |  | 
| 115 | + | // allocate memory for the parallel objects | 
| 116 | + | atypesRow.resize(nAtomsInRow_); | 
| 117 | + | atypesCol.resize(nAtomsInCol_); | 
| 118 | + |  | 
| 119 | + | for (int i = 0; i < nAtomsInRow_; i++) | 
| 120 | + | atypesRow[i] = ff_->getAtomType(identsRow[i]); | 
| 121 | + | for (int i = 0; i < nAtomsInCol_; i++) | 
| 122 | + | atypesCol[i] = ff_->getAtomType(identsCol[i]); | 
| 123 | + |  | 
| 124 | + | pot_row.resize(nAtomsInRow_); | 
| 125 | + | pot_col.resize(nAtomsInCol_); | 
| 126 | + |  | 
| 127 | + | AtomRowToGlobal.resize(nAtomsInRow_); | 
| 128 | + | AtomColToGlobal.resize(nAtomsInCol_); | 
| 129 |  | AtomCommIntRow->gather(AtomLocalToGlobal, AtomRowToGlobal); | 
| 130 |  | AtomCommIntColumn->gather(AtomLocalToGlobal, AtomColToGlobal); | 
| 131 |  |  | 
| 132 | + | cgRowToGlobal.resize(nGroupsInRow_); | 
| 133 | + | cgColToGlobal.resize(nGroupsInCol_); | 
| 134 |  | cgCommIntRow->gather(cgLocalToGlobal, cgRowToGlobal); | 
| 135 |  | cgCommIntColumn->gather(cgLocalToGlobal, cgColToGlobal); | 
| 136 |  |  | 
| 137 | + | massFactorsRow.resize(nAtomsInRow_); | 
| 138 | + | massFactorsCol.resize(nAtomsInCol_); | 
| 139 |  | AtomCommRealRow->gather(massFactors, massFactorsRow); | 
| 140 |  | AtomCommRealColumn->gather(massFactors, massFactorsCol); | 
| 141 |  |  | 
| 194 |  | } | 
| 195 |  |  | 
| 196 |  | #endif | 
| 197 | + |  | 
| 198 | + | // allocate memory for the parallel objects | 
| 199 | + | atypesLocal.resize(nLocal_); | 
| 200 | + |  | 
| 201 | + | for (int i = 0; i < nLocal_; i++) | 
| 202 | + | atypesLocal[i] = ff_->getAtomType(idents[i]); | 
| 203 |  |  | 
| 204 |  | groupList_.clear(); | 
| 205 |  | groupList_.resize(nGroups_); | 
| 253 |  | void ForceMatrixDecomposition::createGtypeCutoffMap() { | 
| 254 |  |  | 
| 255 |  | RealType tol = 1e-6; | 
| 256 | + | largestRcut_ = 0.0; | 
| 257 |  | RealType rc; | 
| 258 |  | int atid; | 
| 259 |  | set<AtomType*> atypes = info_->getSimulatedAtomTypes(); | 
| 260 | + |  | 
| 261 |  | map<int, RealType> atypeCutoff; | 
| 262 |  |  | 
| 263 |  | for (set<AtomType*>::iterator at = atypes.begin(); | 
| 265 |  | atid = (*at)->getIdent(); | 
| 266 |  | if (userChoseCutoff_) | 
| 267 |  | atypeCutoff[atid] = userCutoff_; | 
| 268 | < | else | 
| 268 | > | else | 
| 269 |  | atypeCutoff[atid] = interactionMan_->getSuggestedCutoffRadius(*at); | 
| 270 |  | } | 
| 271 | < |  | 
| 271 | > |  | 
| 272 |  | vector<RealType> gTypeCutoffs; | 
| 273 |  | // first we do a single loop over the cutoff groups to find the | 
| 274 |  | // largest cutoff for any atypes present in this group. | 
| 328 |  | vector<RealType> groupCutoff(nGroups_, 0.0); | 
| 329 |  | groupToGtype.resize(nGroups_); | 
| 330 |  | for (int cg1 = 0; cg1 < nGroups_; cg1++) { | 
| 305 | – |  | 
| 331 |  | groupCutoff[cg1] = 0.0; | 
| 332 |  | vector<int> atomList = getAtomsInGroupRow(cg1); | 
| 308 | – |  | 
| 333 |  | for (vector<int>::iterator ia = atomList.begin(); | 
| 334 |  | ia != atomList.end(); ++ia) { | 
| 335 |  | int atom1 = (*ia); | 
| 336 |  | atid = idents[atom1]; | 
| 337 | < | if (atypeCutoff[atid] > groupCutoff[cg1]) { | 
| 337 | > | if (atypeCutoff[atid] > groupCutoff[cg1]) | 
| 338 |  | groupCutoff[cg1] = atypeCutoff[atid]; | 
| 315 | – | } | 
| 339 |  | } | 
| 340 | < |  | 
| 340 | > |  | 
| 341 |  | bool gTypeFound = false; | 
| 342 |  | for (int gt = 0; gt < gTypeCutoffs.size(); gt++) { | 
| 343 |  | if (abs(groupCutoff[cg1] - gTypeCutoffs[gt]) < tol) { | 
| 345 |  | gTypeFound = true; | 
| 346 |  | } | 
| 347 |  | } | 
| 348 | < | if (!gTypeFound) { | 
| 348 | > | if (!gTypeFound) { | 
| 349 |  | gTypeCutoffs.push_back( groupCutoff[cg1] ); | 
| 350 |  | groupToGtype[cg1] = gTypeCutoffs.size() - 1; | 
| 351 |  | } | 
| 354 |  |  | 
| 355 |  | // Now we find the maximum group cutoff value present in the simulation | 
| 356 |  |  | 
| 357 | < | RealType groupMax = *max_element(gTypeCutoffs.begin(), gTypeCutoffs.end()); | 
| 357 | > | RealType groupMax = *max_element(gTypeCutoffs.begin(), | 
| 358 | > | gTypeCutoffs.end()); | 
| 359 |  |  | 
| 360 |  | #ifdef IS_MPI | 
| 361 | < | MPI::COMM_WORLD.Allreduce(&groupMax, &groupMax, 1, MPI::REALTYPE, MPI::MAX); | 
| 361 | > | MPI::COMM_WORLD.Allreduce(&groupMax, &groupMax, 1, MPI::REALTYPE, | 
| 362 | > | MPI::MAX); | 
| 363 |  | #endif | 
| 364 |  |  | 
| 365 |  | RealType tradRcut = groupMax; | 
| 389 |  |  | 
| 390 |  | pair<int,int> key = make_pair(i,j); | 
| 391 |  | gTypeCutoffMap[key].first = thisRcut; | 
| 367 | – |  | 
| 392 |  | if (thisRcut > largestRcut_) largestRcut_ = thisRcut; | 
| 369 | – |  | 
| 393 |  | gTypeCutoffMap[key].second = thisRcut*thisRcut; | 
| 371 | – |  | 
| 394 |  | gTypeCutoffMap[key].third = pow(thisRcut + skinThickness_, 2); | 
| 373 | – |  | 
| 395 |  | // sanity check | 
| 396 |  |  | 
| 397 |  | if (userChoseCutoff_) { | 
| 451 |  | Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); | 
| 452 |  |  | 
| 453 |  | if (storageLayout_ & DataStorage::dslParticlePot) { | 
| 454 | < | fill(atomRowData.particlePot.begin(), atomRowData.particlePot.end(), 0.0); | 
| 455 | < | fill(atomColData.particlePot.begin(), atomColData.particlePot.end(), 0.0); | 
| 454 | > | fill(atomRowData.particlePot.begin(), atomRowData.particlePot.end(), | 
| 455 | > | 0.0); | 
| 456 | > | fill(atomColData.particlePot.begin(), atomColData.particlePot.end(), | 
| 457 | > | 0.0); | 
| 458 |  | } | 
| 459 |  |  | 
| 460 |  | if (storageLayout_ & DataStorage::dslDensity) { | 
| 463 |  | } | 
| 464 |  |  | 
| 465 |  | if (storageLayout_ & DataStorage::dslFunctional) { | 
| 466 | < | fill(atomRowData.functional.begin(), atomRowData.functional.end(), 0.0); | 
| 467 | < | fill(atomColData.functional.begin(), atomColData.functional.end(), 0.0); | 
| 466 | > | fill(atomRowData.functional.begin(), atomRowData.functional.end(), | 
| 467 | > | 0.0); | 
| 468 | > | fill(atomColData.functional.begin(), atomColData.functional.end(), | 
| 469 | > | 0.0); | 
| 470 |  | } | 
| 471 |  |  | 
| 472 |  | if (storageLayout_ & DataStorage::dslFunctionalDerivative) { | 
| 483 |  | atomColData.skippedCharge.end(), 0.0); | 
| 484 |  | } | 
| 485 |  |  | 
| 486 | < | #else | 
| 487 | < |  | 
| 486 | > | #endif | 
| 487 | > | // even in parallel, we need to zero out the local arrays: | 
| 488 | > |  | 
| 489 |  | if (storageLayout_ & DataStorage::dslParticlePot) { | 
| 490 |  | fill(snap_->atomData.particlePot.begin(), | 
| 491 |  | snap_->atomData.particlePot.end(), 0.0); | 
| 507 |  | fill(snap_->atomData.skippedCharge.begin(), | 
| 508 |  | snap_->atomData.skippedCharge.end(), 0.0); | 
| 509 |  | } | 
| 484 | – | #endif | 
| 510 |  |  | 
| 511 |  | } | 
| 512 |  |  | 
| 543 |  | AtomCommMatrixColumn->gather(snap_->atomData.electroFrame, | 
| 544 |  | atomColData.electroFrame); | 
| 545 |  | } | 
| 546 | + |  | 
| 547 |  | #endif | 
| 548 |  | } | 
| 549 |  |  | 
| 610 |  | AtomCommVectorColumn->scatter(atomColData.force, frc_tmp); | 
| 611 |  | for (int i = 0; i < n; i++) | 
| 612 |  | snap_->atomData.force[i] += frc_tmp[i]; | 
| 613 | < |  | 
| 588 | < |  | 
| 613 | > |  | 
| 614 |  | if (storageLayout_ & DataStorage::dslTorque) { | 
| 615 |  |  | 
| 616 |  | int nt = snap_->atomData.torque.size(); | 
| 634 |  |  | 
| 635 |  | AtomCommRealRow->scatter(atomRowData.skippedCharge, skch_tmp); | 
| 636 |  | for (int i = 0; i < ns; i++) { | 
| 637 | < | snap_->atomData.skippedCharge[i] = skch_tmp[i]; | 
| 637 | > | snap_->atomData.skippedCharge[i] += skch_tmp[i]; | 
| 638 |  | skch_tmp[i] = 0.0; | 
| 639 |  | } | 
| 640 |  |  | 
| 845 |  | idat.excluded = excludeAtomPair(atom1, atom2); | 
| 846 |  |  | 
| 847 |  | #ifdef IS_MPI | 
| 848 | < |  | 
| 849 | < | idat.atypes = make_pair( ff_->getAtomType(identsRow[atom1]), | 
| 850 | < | ff_->getAtomType(identsCol[atom2]) ); | 
| 848 | > | idat.atypes = make_pair( atypesRow[atom1], atypesCol[atom2]); | 
| 849 | > | //idat.atypes = make_pair( ff_->getAtomType(identsRow[atom1]), | 
| 850 | > | //                         ff_->getAtomType(identsCol[atom2]) ); | 
| 851 |  |  | 
| 852 |  | if (storageLayout_ & DataStorage::dslAmat) { | 
| 853 |  | idat.A1 = &(atomRowData.aMat[atom1]); | 
| 891 |  |  | 
| 892 |  | #else | 
| 893 |  |  | 
| 894 | < | idat.atypes = make_pair( ff_->getAtomType(idents[atom1]), | 
| 895 | < | ff_->getAtomType(idents[atom2]) ); | 
| 894 | > | idat.atypes = make_pair( atypesLocal[atom1], atypesLocal[atom2]); | 
| 895 | > | //idat.atypes = make_pair( ff_->getAtomType(idents[atom1]), | 
| 896 | > | //                         ff_->getAtomType(idents[atom2]) ); | 
| 897 |  |  | 
| 898 |  | if (storageLayout_ & DataStorage::dslAmat) { | 
| 899 |  | idat.A1 = &(snap_->atomData.aMat[atom1]); |