# | Line 95 | Line 95 | namespace OpenMD { | |
---|---|---|
95 | storageLayout_ = sman_->getStorageLayout(); | |
96 | ff_ = info_->getForceField(); | |
97 | nLocal_ = snap_->getNumberOfAtoms(); | |
98 | < | |
98 | > | |
99 | nGroups_ = info_->getNLocalCutoffGroups(); | |
100 | // gather the information for atomtype IDs (atids): | |
101 | idents = info_->getIdentArray(); | |
# | Line 109 | Line 109 | namespace OpenMD { | |
109 | PairList* oneTwo = info_->getOneTwoInteractions(); | |
110 | PairList* oneThree = info_->getOneThreeInteractions(); | |
111 | PairList* oneFour = info_->getOneFourInteractions(); | |
112 | < | |
112 | > | |
113 | > | if (needVelocities_) |
114 | > | snap_->cgData.setStorageLayout(DataStorage::dslPosition | |
115 | > | DataStorage::dslVelocity); |
116 | > | else |
117 | > | snap_->cgData.setStorageLayout(DataStorage::dslPosition); |
118 | > | |
119 | #ifdef IS_MPI | |
120 | ||
121 | MPI::Intracomm row = rowComm.getComm(); | |
# | Line 145 | Line 151 | namespace OpenMD { | |
151 | cgRowData.resize(nGroupsInRow_); | |
152 | cgRowData.setStorageLayout(DataStorage::dslPosition); | |
153 | cgColData.resize(nGroupsInCol_); | |
154 | < | cgColData.setStorageLayout(DataStorage::dslPosition); |
155 | < | |
154 | > | if (needVelocities_) |
155 | > | // we only need column velocities if we need them. |
156 | > | cgColData.setStorageLayout(DataStorage::dslPosition | |
157 | > | DataStorage::dslVelocity); |
158 | > | else |
159 | > | cgColData.setStorageLayout(DataStorage::dslPosition); |
160 | > | |
161 | identsRow.resize(nAtomsInRow_); | |
162 | identsCol.resize(nAtomsInCol_); | |
163 | ||
# | Line 164 | Line 175 | namespace OpenMD { | |
175 | ||
176 | pot_row.resize(nAtomsInRow_); | |
177 | pot_col.resize(nAtomsInCol_); | |
178 | + | |
179 | + | expot_row.resize(nAtomsInRow_); |
180 | + | expot_col.resize(nAtomsInCol_); |
181 | ||
182 | AtomRowToGlobal.resize(nAtomsInRow_); | |
183 | AtomColToGlobal.resize(nAtomsInCol_); | |
# | Line 296 | Line 310 | namespace OpenMD { | |
310 | ||
311 | RealType tol = 1e-6; | |
312 | largestRcut_ = 0.0; | |
299 | – | RealType rc; |
313 | int atid; | |
314 | set<AtomType*> atypes = info_->getSimulatedAtomTypes(); | |
315 | ||
# | Line 381 | Line 394 | namespace OpenMD { | |
394 | } | |
395 | ||
396 | bool gTypeFound = false; | |
397 | < | for (int gt = 0; gt < gTypeCutoffs.size(); gt++) { |
397 | > | for (unsigned int gt = 0; gt < gTypeCutoffs.size(); gt++) { |
398 | if (abs(groupCutoff[cg1] - gTypeCutoffs[gt]) < tol) { | |
399 | groupToGtype[cg1] = gt; | |
400 | gTypeFound = true; | |
# | Line 406 | Line 419 | namespace OpenMD { | |
419 | ||
420 | RealType tradRcut = groupMax; | |
421 | ||
422 | < | for (int i = 0; i < gTypeCutoffs.size(); i++) { |
423 | < | for (int j = 0; j < gTypeCutoffs.size(); j++) { |
422 | > | for (unsigned int i = 0; i < gTypeCutoffs.size(); i++) { |
423 | > | for (unsigned int j = 0; j < gTypeCutoffs.size(); j++) { |
424 | RealType thisRcut; | |
425 | switch(cutoffPolicy_) { | |
426 | case TRADITIONAL: | |
# | Line 450 | Line 463 | namespace OpenMD { | |
463 | } | |
464 | } | |
465 | ||
453 | – | |
466 | groupCutoffs ForceMatrixDecomposition::getGroupCutoffs(int cg1, int cg2) { | |
467 | int i, j; | |
468 | #ifdef IS_MPI | |
# | Line 464 | Line 476 | namespace OpenMD { | |
476 | } | |
477 | ||
478 | int ForceMatrixDecomposition::getTopologicalDistance(int atom1, int atom2) { | |
479 | < | for (int j = 0; j < toposForAtom[atom1].size(); j++) { |
479 | > | for (unsigned int j = 0; j < toposForAtom[atom1].size(); j++) { |
480 | if (toposForAtom[atom1][j] == atom2) | |
481 | return topoDist[atom1][j]; | |
482 | } | |
# | Line 474 | Line 486 | namespace OpenMD { | |
486 | void ForceMatrixDecomposition::zeroWorkArrays() { | |
487 | pairwisePot = 0.0; | |
488 | embeddingPot = 0.0; | |
489 | + | excludedPot = 0.0; |
490 | + | excludedSelfPot = 0.0; |
491 | ||
492 | #ifdef IS_MPI | |
493 | if (storageLayout_ & DataStorage::dslForce) { | |
# | Line 492 | Line 506 | namespace OpenMD { | |
506 | fill(pot_col.begin(), pot_col.end(), | |
507 | Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); | |
508 | ||
509 | + | fill(expot_row.begin(), expot_row.end(), |
510 | + | Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); |
511 | + | |
512 | + | fill(expot_col.begin(), expot_col.end(), |
513 | + | Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); |
514 | + | |
515 | if (storageLayout_ & DataStorage::dslParticlePot) { | |
516 | fill(atomRowData.particlePot.begin(), atomRowData.particlePot.end(), | |
517 | 0.0); | |
# | Line 525 | Line 545 | namespace OpenMD { | |
545 | atomColData.skippedCharge.end(), 0.0); | |
546 | } | |
547 | ||
548 | + | if (storageLayout_ & DataStorage::dslFlucQForce) { |
549 | + | fill(atomRowData.flucQFrc.begin(), |
550 | + | atomRowData.flucQFrc.end(), 0.0); |
551 | + | fill(atomColData.flucQFrc.begin(), |
552 | + | atomColData.flucQFrc.end(), 0.0); |
553 | + | } |
554 | + | |
555 | + | if (storageLayout_ & DataStorage::dslElectricField) { |
556 | + | fill(atomRowData.electricField.begin(), |
557 | + | atomRowData.electricField.end(), V3Zero); |
558 | + | fill(atomColData.electricField.begin(), |
559 | + | atomColData.electricField.end(), V3Zero); |
560 | + | } |
561 | + | |
562 | + | if (storageLayout_ & DataStorage::dslFlucQForce) { |
563 | + | fill(atomRowData.flucQFrc.begin(), atomRowData.flucQFrc.end(), |
564 | + | 0.0); |
565 | + | fill(atomColData.flucQFrc.begin(), atomColData.flucQFrc.end(), |
566 | + | 0.0); |
567 | + | } |
568 | + | |
569 | #endif | |
570 | // even in parallel, we need to zero out the local arrays: | |
571 | ||
# | Line 537 | Line 578 | namespace OpenMD { | |
578 | fill(snap_->atomData.density.begin(), | |
579 | snap_->atomData.density.end(), 0.0); | |
580 | } | |
581 | + | |
582 | if (storageLayout_ & DataStorage::dslFunctional) { | |
583 | fill(snap_->atomData.functional.begin(), | |
584 | snap_->atomData.functional.end(), 0.0); | |
585 | } | |
586 | + | |
587 | if (storageLayout_ & DataStorage::dslFunctionalDerivative) { | |
588 | fill(snap_->atomData.functionalDerivative.begin(), | |
589 | snap_->atomData.functionalDerivative.end(), 0.0); | |
590 | } | |
591 | + | |
592 | if (storageLayout_ & DataStorage::dslSkippedCharge) { | |
593 | fill(snap_->atomData.skippedCharge.begin(), | |
594 | snap_->atomData.skippedCharge.end(), 0.0); | |
595 | } | |
596 | < | |
596 | > | |
597 | > | if (storageLayout_ & DataStorage::dslElectricField) { |
598 | > | fill(snap_->atomData.electricField.begin(), |
599 | > | snap_->atomData.electricField.end(), V3Zero); |
600 | > | } |
601 | } | |
602 | ||
603 | ||
# | Line 572 | Line 620 | namespace OpenMD { | |
620 | cgPlanVectorColumn->gather(snap_->cgData.position, | |
621 | cgColData.position); | |
622 | ||
623 | + | |
624 | + | |
625 | + | if (needVelocities_) { |
626 | + | // gather up the atomic velocities |
627 | + | AtomPlanVectorColumn->gather(snap_->atomData.velocity, |
628 | + | atomColData.velocity); |
629 | + | |
630 | + | cgPlanVectorColumn->gather(snap_->cgData.velocity, |
631 | + | cgColData.velocity); |
632 | + | } |
633 | + | |
634 | ||
635 | // if needed, gather the atomic rotation matrices | |
636 | if (storageLayout_ & DataStorage::dslAmat) { | |
# | Line 589 | Line 648 | namespace OpenMD { | |
648 | atomColData.electroFrame); | |
649 | } | |
650 | ||
651 | + | // if needed, gather the atomic fluctuating charge values |
652 | + | if (storageLayout_ & DataStorage::dslFlucQPosition) { |
653 | + | AtomPlanRealRow->gather(snap_->atomData.flucQPos, |
654 | + | atomRowData.flucQPos); |
655 | + | AtomPlanRealColumn->gather(snap_->atomData.flucQPos, |
656 | + | atomColData.flucQPos); |
657 | + | } |
658 | + | |
659 | #endif | |
660 | } | |
661 | ||
# | Line 611 | Line 678 | namespace OpenMD { | |
678 | for (int i = 0; i < n; i++) | |
679 | snap_->atomData.density[i] += rho_tmp[i]; | |
680 | } | |
681 | + | |
682 | + | if (storageLayout_ & DataStorage::dslElectricField) { |
683 | + | |
684 | + | AtomPlanVectorRow->scatter(atomRowData.electricField, |
685 | + | snap_->atomData.electricField); |
686 | + | |
687 | + | int n = snap_->atomData.electricField.size(); |
688 | + | vector<Vector3d> field_tmp(n, V3Zero); |
689 | + | AtomPlanVectorColumn->scatter(atomColData.electricField, field_tmp); |
690 | + | for (int i = 0; i < n; i++) |
691 | + | snap_->atomData.electricField[i] += field_tmp[i]; |
692 | + | } |
693 | #endif | |
694 | } | |
695 | ||
# | Line 690 | Line 769 | namespace OpenMD { | |
769 | ||
770 | } | |
771 | ||
772 | + | if (storageLayout_ & DataStorage::dslFlucQForce) { |
773 | + | |
774 | + | int nq = snap_->atomData.flucQFrc.size(); |
775 | + | vector<RealType> fqfrc_tmp(nq, 0.0); |
776 | + | |
777 | + | AtomPlanRealRow->scatter(atomRowData.flucQFrc, fqfrc_tmp); |
778 | + | for (int i = 0; i < nq; i++) { |
779 | + | snap_->atomData.flucQFrc[i] += fqfrc_tmp[i]; |
780 | + | fqfrc_tmp[i] = 0.0; |
781 | + | } |
782 | + | |
783 | + | AtomPlanRealColumn->scatter(atomColData.flucQFrc, fqfrc_tmp); |
784 | + | for (int i = 0; i < nq; i++) |
785 | + | snap_->atomData.flucQFrc[i] += fqfrc_tmp[i]; |
786 | + | |
787 | + | } |
788 | + | |
789 | nLocal_ = snap_->getNumberOfAtoms(); | |
790 | ||
791 | vector<potVec> pot_temp(nLocal_, | |
792 | Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); | |
793 | + | vector<potVec> expot_temp(nLocal_, |
794 | + | Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); |
795 | ||
796 | // scatter/gather pot_row into the members of my column | |
797 | ||
798 | AtomPlanPotRow->scatter(pot_row, pot_temp); | |
799 | + | AtomPlanPotRow->scatter(expot_row, expot_temp); |
800 | ||
801 | < | for (int ii = 0; ii < pot_temp.size(); ii++ ) |
801 | > | for (int ii = 0; ii < pot_temp.size(); ii++ ) |
802 | pairwisePot += pot_temp[ii]; | |
803 | < | |
803 | > | |
804 | > | for (int ii = 0; ii < expot_temp.size(); ii++ ) |
805 | > | excludedPot += expot_temp[ii]; |
806 | > | |
807 | > | if (storageLayout_ & DataStorage::dslParticlePot) { |
808 | > | // This is the pairwise contribution to the particle pot. The |
809 | > | // embedding contribution is added in each of the low level |
810 | > | // non-bonded routines. In single processor, this is done in |
811 | > | // unpackInteractionData, not in collectData. |
812 | > | for (int ii = 0; ii < N_INTERACTION_FAMILIES; ii++) { |
813 | > | for (int i = 0; i < nLocal_; i++) { |
814 | > | // factor of two is because the total potential terms are divided |
815 | > | // by 2 in parallel due to row/ column scatter |
816 | > | snap_->atomData.particlePot[i] += 2.0 * pot_temp[i](ii); |
817 | > | } |
818 | > | } |
819 | > | } |
820 | > | |
821 | fill(pot_temp.begin(), pot_temp.end(), | |
822 | Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); | |
823 | + | fill(expot_temp.begin(), expot_temp.end(), |
824 | + | Vector<RealType, N_INTERACTION_FAMILIES> (0.0)); |
825 | ||
826 | AtomPlanPotColumn->scatter(pot_col, pot_temp); | |
827 | + | AtomPlanPotColumn->scatter(expot_col, expot_temp); |
828 | ||
829 | for (int ii = 0; ii < pot_temp.size(); ii++ ) | |
830 | pairwisePot += pot_temp[ii]; | |
831 | + | |
832 | + | for (int ii = 0; ii < expot_temp.size(); ii++ ) |
833 | + | excludedPot += expot_temp[ii]; |
834 | + | |
835 | + | if (storageLayout_ & DataStorage::dslParticlePot) { |
836 | + | // This is the pairwise contribution to the particle pot. The |
837 | + | // embedding contribution is added in each of the low level |
838 | + | // non-bonded routines. In single processor, this is done in |
839 | + | // unpackInteractionData, not in collectData. |
840 | + | for (int ii = 0; ii < N_INTERACTION_FAMILIES; ii++) { |
841 | + | for (int i = 0; i < nLocal_; i++) { |
842 | + | // factor of two is because the total potential terms are divided |
843 | + | // by 2 in parallel due to row/ column scatter |
844 | + | snap_->atomData.particlePot[i] += 2.0 * pot_temp[i](ii); |
845 | + | } |
846 | + | } |
847 | + | } |
848 | ||
849 | + | if (storageLayout_ & DataStorage::dslParticlePot) { |
850 | + | int npp = snap_->atomData.particlePot.size(); |
851 | + | vector<RealType> ppot_temp(npp, 0.0); |
852 | + | |
853 | + | // This is the direct or embedding contribution to the particle |
854 | + | // pot. |
855 | + | |
856 | + | AtomPlanRealRow->scatter(atomRowData.particlePot, ppot_temp); |
857 | + | for (int i = 0; i < npp; i++) { |
858 | + | snap_->atomData.particlePot[i] += ppot_temp[i]; |
859 | + | } |
860 | + | |
861 | + | fill(ppot_temp.begin(), ppot_temp.end(), 0.0); |
862 | + | |
863 | + | AtomPlanRealColumn->scatter(atomColData.particlePot, ppot_temp); |
864 | + | for (int i = 0; i < npp; i++) { |
865 | + | snap_->atomData.particlePot[i] += ppot_temp[i]; |
866 | + | } |
867 | + | } |
868 | + | |
869 | for (int ii = 0; ii < N_INTERACTION_FAMILIES; ii++) { | |
870 | RealType ploc1 = pairwisePot[ii]; | |
871 | RealType ploc2 = 0.0; | |
# | Line 718 | Line 874 | namespace OpenMD { | |
874 | } | |
875 | ||
876 | for (int ii = 0; ii < N_INTERACTION_FAMILIES; ii++) { | |
877 | < | RealType ploc1 = embeddingPot[ii]; |
877 | > | RealType ploc1 = excludedPot[ii]; |
878 | RealType ploc2 = 0.0; | |
879 | MPI::COMM_WORLD.Allreduce(&ploc1, &ploc2, 1, MPI::REALTYPE, MPI::SUM); | |
880 | < | embeddingPot[ii] = ploc2; |
880 | > | excludedPot[ii] = ploc2; |
881 | } | |
882 | ||
883 | + | // Here be dragons. |
884 | + | MPI::Intracomm col = colComm.getComm(); |
885 | + | |
886 | + | col.Allreduce(MPI::IN_PLACE, |
887 | + | &snap_->frameData.conductiveHeatFlux[0], 3, |
888 | + | MPI::REALTYPE, MPI::SUM); |
889 | + | |
890 | + | |
891 | #endif | |
892 | ||
893 | } | |
894 | ||
895 | + | /** |
896 | + | * Collects information obtained during the post-pair (and embedding |
897 | + | * functional) loops onto local data structures. |
898 | + | */ |
899 | + | void ForceMatrixDecomposition::collectSelfData() { |
900 | + | snap_ = sman_->getCurrentSnapshot(); |
901 | + | storageLayout_ = sman_->getStorageLayout(); |
902 | + | |
903 | + | #ifdef IS_MPI |
904 | + | for (int ii = 0; ii < N_INTERACTION_FAMILIES; ii++) { |
905 | + | RealType ploc1 = embeddingPot[ii]; |
906 | + | RealType ploc2 = 0.0; |
907 | + | MPI::COMM_WORLD.Allreduce(&ploc1, &ploc2, 1, MPI::REALTYPE, MPI::SUM); |
908 | + | embeddingPot[ii] = ploc2; |
909 | + | } |
910 | + | for (int ii = 0; ii < N_INTERACTION_FAMILIES; ii++) { |
911 | + | RealType ploc1 = excludedSelfPot[ii]; |
912 | + | RealType ploc2 = 0.0; |
913 | + | MPI::COMM_WORLD.Allreduce(&ploc1, &ploc2, 1, MPI::REALTYPE, MPI::SUM); |
914 | + | excludedSelfPot[ii] = ploc2; |
915 | + | } |
916 | + | #endif |
917 | + | |
918 | + | } |
919 | + | |
920 | + | |
921 | + | |
922 | int ForceMatrixDecomposition::getNAtomsInRow() { | |
923 | #ifdef IS_MPI | |
924 | return nAtomsInRow_; | |
# | Line 768 | Line 959 | namespace OpenMD { | |
959 | return d; | |
960 | } | |
961 | ||
962 | + | Vector3d ForceMatrixDecomposition::getGroupVelocityColumn(int cg2){ |
963 | + | #ifdef IS_MPI |
964 | + | return cgColData.velocity[cg2]; |
965 | + | #else |
966 | + | return snap_->cgData.velocity[cg2]; |
967 | + | #endif |
968 | + | } |
969 | ||
970 | + | Vector3d ForceMatrixDecomposition::getAtomVelocityColumn(int atom2){ |
971 | + | #ifdef IS_MPI |
972 | + | return atomColData.velocity[atom2]; |
973 | + | #else |
974 | + | return snap_->atomData.velocity[atom2]; |
975 | + | #endif |
976 | + | } |
977 | + | |
978 | + | |
979 | Vector3d ForceMatrixDecomposition::getAtomToGroupVectorRow(int atom1, int cg1){ | |
980 | ||
981 | Vector3d d; | |
# | Line 834 | Line 1041 | namespace OpenMD { | |
1041 | * We need to exclude some overcounted interactions that result from | |
1042 | * the parallel decomposition. | |
1043 | */ | |
1044 | < | bool ForceMatrixDecomposition::skipAtomPair(int atom1, int atom2) { |
1044 | > | bool ForceMatrixDecomposition::skipAtomPair(int atom1, int atom2, int cg1, int cg2) { |
1045 | int unique_id_1, unique_id_2; | |
1046 | ||
1047 | #ifdef IS_MPI | |
1048 | // in MPI, we have to look up the unique IDs for each atom | |
1049 | unique_id_1 = AtomRowToGlobal[atom1]; | |
1050 | unique_id_2 = AtomColToGlobal[atom2]; | |
1051 | + | // group1 = cgRowToGlobal[cg1]; |
1052 | + | // group2 = cgColToGlobal[cg2]; |
1053 | #else | |
1054 | unique_id_1 = AtomLocalToGlobal[atom1]; | |
1055 | unique_id_2 = AtomLocalToGlobal[atom2]; | |
1056 | + | int group1 = cgLocalToGlobal[cg1]; |
1057 | + | int group2 = cgLocalToGlobal[cg2]; |
1058 | #endif | |
1059 | ||
1060 | if (unique_id_1 == unique_id_2) return true; | |
# | Line 855 | Line 1066 | namespace OpenMD { | |
1066 | } else { | |
1067 | if ((unique_id_1 + unique_id_2) % 2 == 1) return true; | |
1068 | } | |
1069 | + | #endif |
1070 | + | |
1071 | + | #ifndef IS_MPI |
1072 | + | if (group1 == group2) { |
1073 | + | if (unique_id_1 < unique_id_2) return true; |
1074 | + | } |
1075 | #endif | |
1076 | ||
1077 | return false; | |
# | Line 950 | Line 1167 | namespace OpenMD { | |
1167 | idat.skippedCharge2 = &(atomColData.skippedCharge[atom2]); | |
1168 | } | |
1169 | ||
1170 | < | #else |
1170 | > | if (storageLayout_ & DataStorage::dslFlucQPosition) { |
1171 | > | idat.flucQ1 = &(atomRowData.flucQPos[atom1]); |
1172 | > | idat.flucQ2 = &(atomColData.flucQPos[atom2]); |
1173 | > | } |
1174 | ||
1175 | + | #else |
1176 | + | |
1177 | idat.atypes = make_pair( atypesLocal[atom1], atypesLocal[atom2]); | |
956 | – | //idat.atypes = make_pair( ff_->getAtomType(idents[atom1]), |
957 | – | // ff_->getAtomType(idents[atom2]) ); |
1178 | ||
1179 | if (storageLayout_ & DataStorage::dslAmat) { | |
1180 | idat.A1 = &(snap_->atomData.aMat[atom1]); | |
# | Line 995 | Line 1215 | namespace OpenMD { | |
1215 | idat.skippedCharge1 = &(snap_->atomData.skippedCharge[atom1]); | |
1216 | idat.skippedCharge2 = &(snap_->atomData.skippedCharge[atom2]); | |
1217 | } | |
1218 | + | |
1219 | + | if (storageLayout_ & DataStorage::dslFlucQPosition) { |
1220 | + | idat.flucQ1 = &(snap_->atomData.flucQPos[atom1]); |
1221 | + | idat.flucQ2 = &(snap_->atomData.flucQPos[atom2]); |
1222 | + | } |
1223 | + | |
1224 | #endif | |
1225 | } | |
1226 | ||
# | Line 1003 | Line 1229 | namespace OpenMD { | |
1229 | #ifdef IS_MPI | |
1230 | pot_row[atom1] += RealType(0.5) * *(idat.pot); | |
1231 | pot_col[atom2] += RealType(0.5) * *(idat.pot); | |
1232 | + | expot_row[atom1] += RealType(0.5) * *(idat.excludedPot); |
1233 | + | expot_col[atom2] += RealType(0.5) * *(idat.excludedPot); |
1234 | ||
1235 | atomRowData.force[atom1] += *(idat.f1); | |
1236 | atomColData.force[atom2] -= *(idat.f1); | |
1237 | + | |
1238 | + | if (storageLayout_ & DataStorage::dslFlucQForce) { |
1239 | + | atomRowData.flucQFrc[atom1] -= *(idat.dVdFQ1); |
1240 | + | atomColData.flucQFrc[atom2] -= *(idat.dVdFQ2); |
1241 | + | } |
1242 | + | |
1243 | + | if (storageLayout_ & DataStorage::dslElectricField) { |
1244 | + | atomRowData.electricField[atom1] += *(idat.eField1); |
1245 | + | atomColData.electricField[atom2] += *(idat.eField2); |
1246 | + | } |
1247 | + | |
1248 | #else | |
1249 | pairwisePot += *(idat.pot); | |
1250 | + | excludedPot += *(idat.excludedPot); |
1251 | ||
1252 | snap_->atomData.force[atom1] += *(idat.f1); | |
1253 | snap_->atomData.force[atom2] -= *(idat.f1); | |
1254 | + | |
1255 | + | if (idat.doParticlePot) { |
1256 | + | // This is the pairwise contribution to the particle pot. The |
1257 | + | // embedding contribution is added in each of the low level |
1258 | + | // non-bonded routines. In parallel, this calculation is done |
1259 | + | // in collectData, not in unpackInteractionData. |
1260 | + | snap_->atomData.particlePot[atom1] += *(idat.vpair) * *(idat.sw); |
1261 | + | snap_->atomData.particlePot[atom2] += *(idat.vpair) * *(idat.sw); |
1262 | + | } |
1263 | + | |
1264 | + | if (storageLayout_ & DataStorage::dslFlucQForce) { |
1265 | + | snap_->atomData.flucQFrc[atom1] -= *(idat.dVdFQ1); |
1266 | + | snap_->atomData.flucQFrc[atom2] -= *(idat.dVdFQ2); |
1267 | + | } |
1268 | + | |
1269 | + | if (storageLayout_ & DataStorage::dslElectricField) { |
1270 | + | snap_->atomData.electricField[atom1] += *(idat.eField1); |
1271 | + | snap_->atomData.electricField[atom2] += *(idat.eField2); |
1272 | + | } |
1273 | + | |
1274 | #endif | |
1275 | ||
1276 | } | |
# | Line 1035 | Line 1295 | namespace OpenMD { | |
1295 | #endif | |
1296 | ||
1297 | RealType rList_ = (largestRcut_ + skinThickness_); | |
1038 | – | RealType rl2 = rList_ * rList_; |
1298 | Snapshot* snap_ = sman_->getCurrentSnapshot(); | |
1299 | Mat3x3d Hmat = snap_->getHmat(); | |
1300 | Vector3d Hx = Hmat.getColumn(0); | |
# | Line 1079 | Line 1338 | namespace OpenMD { | |
1338 | for (int j = 0; j < 3; j++) { | |
1339 | scaled[j] -= roundMe(scaled[j]); | |
1340 | scaled[j] += 0.5; | |
1341 | + | // Handle the special case when an object is exactly on the |
1342 | + | // boundary (a scaled coordinate of 1.0 is the same as |
1343 | + | // scaled coordinate of 0.0) |
1344 | + | if (scaled[j] >= 1.0) scaled[j] -= 1.0; |
1345 | } | |
1346 | ||
1347 | // find xyz-indices of cell that cutoffGroup is in. | |
# | Line 1103 | Line 1366 | namespace OpenMD { | |
1366 | for (int j = 0; j < 3; j++) { | |
1367 | scaled[j] -= roundMe(scaled[j]); | |
1368 | scaled[j] += 0.5; | |
1369 | + | // Handle the special case when an object is exactly on the |
1370 | + | // boundary (a scaled coordinate of 1.0 is the same as |
1371 | + | // scaled coordinate of 0.0) |
1372 | + | if (scaled[j] >= 1.0) scaled[j] -= 1.0; |
1373 | } | |
1374 | ||
1375 | // find xyz-indices of cell that cutoffGroup is in. | |
# | Line 1129 | Line 1396 | namespace OpenMD { | |
1396 | for (int j = 0; j < 3; j++) { | |
1397 | scaled[j] -= roundMe(scaled[j]); | |
1398 | scaled[j] += 0.5; | |
1399 | + | // Handle the special case when an object is exactly on the |
1400 | + | // boundary (a scaled coordinate of 1.0 is the same as |
1401 | + | // scaled coordinate of 0.0) |
1402 | + | if (scaled[j] >= 1.0) scaled[j] -= 1.0; |
1403 | } | |
1404 | ||
1405 | // find xyz-indices of cell that cutoffGroup is in. |
– | Removed lines |
+ | Added lines |
< | Changed lines |
> | Changed lines |