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

Comparing branches/devel_omp/src/brains/ForceManager.cpp (file contents):
Revision 1598 by mciznick, Wed Jul 27 14:26:53 2011 UTC vs.
Revision 1614 by mciznick, Tue Aug 23 20:55:51 2011 UTC

# Line 63 | Line 63
63   #include <iomanip>
64  
65   #include <omp.h>
66 + //#include <time.h>
67 + #include <sys/time.h>
68  
69   using namespace std;
70   namespace OpenMD {
71 +
72 + //long int elapsedTime = 0;
73  
74   ForceManager::ForceManager(SimInfo * info) :
75          info_(info) {
# Line 366 | Line 370 | void ForceManager::calcForces() {
370          preCalculation();
371          shortRangeInteractions();
372          //    longRangeInteractions();
373 <        longRangeInteractionsRapaport();
373 >        //      longRangeInteractionsRapaport();
374 >        longRangeInteractionsParallel();
375          postCalculation();
376   }
377  
# Line 529 | Line 534 | void ForceManager::shortRangeInteractions() {
534          curSnapshot->statData[Stats::DIHEDRAL_POTENTIAL] = torsionPotential;
535          curSnapshot->statData[Stats::INVERSION_POTENTIAL] = inversionPotential;
536   }
537 +
538 + void ForceManager::longRangeInteractionsParallel() {
539 +
540 +        Snapshot* curSnapshot = info_->getSnapshotManager()->getCurrentSnapshot();
541 +        DataStorage* config = &(curSnapshot->atomData);
542 +        DataStorage* cgConfig = &(curSnapshot->cgData);
543  
544 +        //calculate the center of mass of cutoff group
545 +
546 +        SimInfo::MoleculeIterator mi;
547 +        Molecule* mol;
548 +        Molecule::CutoffGroupIterator ci;
549 +        CutoffGroup* cg;
550 +
551 +        if (info_->getNCutoffGroups() > 0)
552 +        {
553 +                for (mol = info_->beginMolecule(mi); mol != NULL; mol = info_->nextMolecule(mi))
554 +                {
555 +                        for (cg = mol->beginCutoffGroup(ci); cg != NULL; cg = mol->nextCutoffGroup(ci))
556 +                        {
557 +                                //                              cerr << "branch1\n";
558 +                                //                              cerr << "globind = " << cg->getGlobalIndex() << ":" << __LINE__ << "\n";
559 +                                cg->updateCOM();
560 +
561 +                                //                              cerr << "gbI: " << cg->getGlobalIndex() << " locI: " << cg->getLocalIndex() << " x: "
562 +                                //                                              << cgConfig->position[cg->getLocalIndex()].x() << " y: " << cgConfig->position[cg->getLocalIndex()].y()
563 +                                //                                              << " z: " << cgConfig->position[cg->getLocalIndex()].z() << "\n";
564 +                        }
565 +                }
566 +        } else
567 +        {
568 +                // center of mass of the group is the same as position of the atom
569 +                // if cutoff group does not exist
570 +                //              cerr << ":" << __LINE__ << "branch2\n";
571 +                cgConfig->position = config->position;
572 +        }
573 +
574 +        fDecomp_->zeroWorkArrays();
575 +        fDecomp_->distributeData();
576 +
577 +        int atom1, atom2, topoDist;
578 +        Vector3d d_grp, dag, d;
579 +        RealType rgrpsq, rgrp, r2, r;
580 +        RealType electroMult, vdwMult;
581 +        RealType vij;
582 +        Vector3d fij, fg, f1;
583 +        tuple3<RealType, RealType, RealType> cuts;
584 +        RealType rCutSq;
585 +        bool in_switching_region;
586 +        RealType sw, dswdr, swderiv;
587 +        vector<int> atomListColumn, atomListRow, atomListLocal;
588 +
589 +        /* Defines local interaction data to each thread */
590 +        InteractionDataPrv idatPrv;
591 +
592 +        SelfData sdat;
593 +        RealType mf;
594 +        RealType lrPot;
595 +        RealType vpair;
596 +        potVec longRangePotential(0.0);
597 +        potVec workPot(0.0);
598 +
599 +        int loopStart, loopEnd;
600 +        sdat.pot = fDecomp_->getEmbeddingPotential();
601 +
602 +        vector<CutoffGroup *> cgs;
603 +
604 +        for (mol = info_->beginMolecule(mi); mol != NULL; mol = info_->nextMolecule(mi))
605 +        {
606 +                for (cg = mol->beginCutoffGroup(ci); cg != NULL; cg = mol->nextCutoffGroup(ci))
607 +                {
608 +                        cgs.push_back(cg);
609 +                }
610 +        }
611 +
612 +        loopEnd = PAIR_LOOP;
613 +        if (info_->requiresPrepair())
614 +        {
615 +                loopStart = PREPAIR_LOOP;
616 +        } else
617 +        {
618 +                loopStart = PAIR_LOOP;
619 +        }
620 +
621 +        for (int iLoop = loopStart; iLoop <= loopEnd; iLoop++)
622 +        {
623 +
624 +                if (iLoop == loopStart)
625 +                {
626 +                        bool update_nlist = fDecomp_->checkNeighborList();
627 +                        if (update_nlist)
628 +                                neighborMatW = fDecomp_->buildLayerBasedNeighborList();
629 +                }
630 +
631 +                /* Eager initialization */
632 +                /* Initializes InteractionManager before force calculations */
633 +                interactionMan_->initializeOMP();
634 +                /* Initializes forces used in simulation before force calculations */
635 +                interactionMan_->initNonbondedForces();
636 +
637 +                vector<CutoffGroup *>::iterator cg1;
638 +                vector<CutoffGroup *>::iterator cg2;
639 +
640 +                /* Defines local accumulators for each thread */
641 +                vector<Vector3d> forceLcl;
642 +                Vector3d fatom1Lcl;
643 +                Mat3x3d tauLcl;
644 +                potVec potLcl;
645 +
646 +                /* Defines the size of chunks into which the loop iterations will be split */
647 +                int chunkSize = cgs.size() / (omp_get_max_threads() * 5);
648 +
649 +                /*struct timeval tv, tv2;
650 +
651 +                gettimeofday(&tv, NULL);*/
652 +
653 +                /* Defines the parallel region and the list of variables that should be shared and private to each thread */
654 +                #pragma omp parallel default(none) shared(curSnapshot, iLoop, cgs, chunkSize, config) \
655 +                        private(cg1, cg2, cuts, d_grp, rgrpsq, rCutSq, idatPrv, vij, fij, in_switching_region, \
656 +                        dswdr, rgrp, atomListRow, atomListColumn, atom1, atom2, topoDist, d, r2, swderiv, fg, mf, \
657 +                        dag, tauLcl, forceLcl, fatom1Lcl, potLcl)
658 +                {
659 +                        idatPrv.shiftedPot = (cutoffMethod_ == SHIFTED_POTENTIAL) ? true : false;
660 +                        idatPrv.shiftedForce = (cutoffMethod_ == SHIFTED_FORCE) ? true : false;
661 +                        forceLcl = vector<Vector3d>(config->force.size());
662 +                        tauLcl *= 0;
663 +
664 +                        /* Spreads force calculations between threads. Each thread receives chunkSize portion of the CutoffGroups. */
665 +                        #pragma omp for schedule(dynamic, chunkSize)
666 +                        for (cg1 = cgs.begin(); cg1 < cgs.end(); ++cg1)
667 +                        {
668 +                                /* Iterates between neighbors of the CutoffGroup */
669 +                                for (cg2 = neighborMatW[(*cg1)->getGlobalIndex()].begin(); cg2 < neighborMatW[(*cg1)->getGlobalIndex()].end(); ++cg2)
670 +                                {
671 +
672 +                                        cuts = fDecomp_->getGroupCutoffs((*cg1)->getGlobalIndex(), (*cg2)->getGlobalIndex());
673 +
674 +                                        d_grp = fDecomp_->getIntergroupVector((*cg1), (*cg2));
675 +                                        curSnapshot->wrapVector(d_grp);
676 +                                        rgrpsq = d_grp.lengthSquare();
677 +
678 +                                        rCutSq = cuts.second;
679 +
680 + //                                      printf("Thread %d\tcg1:%d\tcg2:%d d_grp\tx:%f\ty:%f\tz:%f\trgrpsq:%f\n", omp_get_thread_num(), (*cg1)->getGlobalIndex(), (*cg2)->getGlobalIndex(), d_grp.x(), d_grp.y(), d_grp.z(), rgrpsq);
681 +
682 +                                        if (rgrpsq < rCutSq)
683 +                                        {
684 +                                                idatPrv.rcut = cuts.first;
685 +                                                if (iLoop == PAIR_LOOP)
686 +                                                {
687 +                                                        vij = 0.0;
688 +                                                        fij = V3Zero;
689 +                                                }
690 +
691 +                                                in_switching_region = switcher_->getSwitch(rgrpsq, idatPrv.sw, dswdr, rgrp);
692 +
693 + //                                              printf("in_switching_region:%d\trgrpsq:%f\t*idatPrv.sw:%f\tdswdr:%f\trgrp:%f\n", (in_switching_region == false ? 0 : 1), rgrpsq, idatPrv.sw, dswdr, rgrp);
694 +
695 +                                                atomListRow = fDecomp_->getAtomsInGroupRow((*cg1)->getGlobalIndex());
696 +                                                atomListColumn = fDecomp_->getAtomsInGroupColumn((*cg2)->getGlobalIndex());
697 +
698 +                                                for (vector<int>::iterator ia = atomListRow.begin(); ia != atomListRow.end(); ++ia)
699 +                                                {
700 +                                                        atom1 = (*ia);
701 +                                                        fatom1Lcl = V3Zero;
702 +
703 +                                                        for (vector<int>::iterator jb = atomListColumn.begin(); jb != atomListColumn.end(); ++jb)
704 +                                                        {
705 +                                                                atom2 = (*jb);
706 +
707 +                                                                if (!fDecomp_->skipAtomPair(atom1, atom2))
708 +                                                                {
709 +                                                                        idatPrv.vpair = 0.0;
710 +                                                                        idatPrv.pot = 0.0;
711 +                                                                        idatPrv.f1 = V3Zero;
712 +
713 +                                                                        fDecomp_->fillInteractionDataOMP(idatPrv, atom1, atom2);
714 +
715 +                                                                        topoDist = fDecomp_->getTopologicalDistance(atom1, atom2);
716 +                                                                        idatPrv.vdwMult = vdwScale_[topoDist];
717 +                                                                        idatPrv.electroMult = electrostaticScale_[topoDist];
718 +
719 + //                                                                      printf("topoDist:%d\tidatPrv.vdwMult:%f\tidatPrv.electroMult:%f\n", topoDist, idatPrv.vdwMult, idatPrv.electroMult);
720 +
721 +                                                                        if (atomListRow.size() == 1 && atomListColumn.size() == 1)
722 +                                                                        {
723 +                                                                                idatPrv.d = d_grp;
724 +                                                                                idatPrv.r2 = rgrpsq;
725 +                                                                                //cerr << "dgrp = " << d_grp << ":" << __LINE__ << "\n";
726 +                                                                        } else
727 +                                                                        {
728 +                                                                                d = fDecomp_->getInteratomicVector(atom1, atom2);
729 +                                                                                curSnapshot->wrapVector(d);
730 +                                                                                r2 = d.lengthSquare();
731 +                                                                                //cerr << "datm = " << d << ":" << __LINE__ << "\n";
732 +                                                                                idatPrv.d = d;
733 +                                                                                idatPrv.r2 = r2;
734 +                                                                        }
735 +
736 +                                                                        //printf("idatPrv.d x:%f\ty:%f\tz:%f\tidatPrv.r2:%f\n", (idatPrv.d).x(), (idatPrv.d).y(), (idatPrv.d).z(), idatPrv.r2);
737 +                                                                        //cerr << "idat.d = " << *(idat.d) << ":" << __LINE__ << "\n";
738 +                                                                        idatPrv.rij = sqrt((idatPrv.r2));
739 +                                                                        //cerr << "idat.rij = " << *(idat.rij) << "\n";
740 +
741 +                                                                        if (iLoop == PREPAIR_LOOP)
742 +                                                                        {
743 +                                                                                interactionMan_->doPrePairOMP(idatPrv);
744 +                                                                        } else
745 +                                                                        {
746 +                                                                                interactionMan_->doPairOMP(idatPrv);
747 +
748 +                                                                                /* Accumulates potential and forces in local arrays */
749 +                                                                                potLcl += idatPrv.pot;
750 +                                                                                fatom1Lcl += idatPrv.f1;
751 +                                                                                forceLcl[atom2] -= idatPrv.f1;
752 +                                                                                //cerr << "d = " << *(idat.d) << "\tv=" << vpair << "\tf=" << f1 << ":" << __LINE__ << "\n";
753 +                                                                                //printf("d x:%f y:%f z:%f vpair:%f f1 x:%f y:%f z:%f\n", idatPrv.d.x(), idatPrv.d.y(), idatPrv.d.z(), idatPrv.vpair, idatPrv.f1.x(), idatPrv.f1.y(), idatPrv.f1.z());
754 +                                                                                vij += idatPrv.vpair;
755 +                                                                                fij += idatPrv.f1;
756 +                                                                                tauLcl -= outProduct(idatPrv.d, idatPrv.f1);
757 +                                                                                //printf("vij:%f fij x:%f y:%f z:%f\n", vij, fij.x(), fij.y(), fij.z());
758 +                                                                        }
759 +                                                                }
760 +                                                        }
761 +                                                        /* We can accumulate force for current CutoffGroup in shared memory without worring about read after write bugs*/
762 +                                                        config->force[atom1] += fatom1Lcl;
763 +                                                }
764 +
765 +                                                if (iLoop == PAIR_LOOP)
766 +                                                {
767 +                                                        if (in_switching_region)
768 +                                                        {
769 +                                                                swderiv = vij * dswdr / rgrp;
770 +                                                                fg = swderiv * d_grp;
771 +                                                                fij += fg;
772 +
773 +                                                                if (atomListRow.size() == 1 && atomListColumn.size() == 1)
774 +                                                                {
775 +                                                                                tauLcl -= outProduct(idatPrv.d, fg);
776 +                                                                }
777 +
778 +                                                                for (vector<int>::iterator ia = atomListRow.begin(); ia != atomListRow.end(); ++ia)
779 +                                                                {
780 +                                                                        atom1 = (*ia);
781 +                                                                        mf = fDecomp_->getMassFactorRow(atom1);
782 +                                                                        // fg is the force on atom ia due to cutoff group's
783 +                                                                        // presence in switching region
784 +                                                                        fg = swderiv * d_grp * mf;
785 +                                                                        #pragma omp critical (forceLck)
786 +                                                                        {
787 +                                                                                fDecomp_->addForceToAtomRow(atom1, fg);
788 +                                                                        }
789 +
790 +                                                                        if (atomListRow.size() > 1)
791 +                                                                        {
792 +                                                                                if (info_->usesAtomicVirial())
793 +                                                                                {
794 +                                                                                        // find the distance between the atom
795 +                                                                                        // and the center of the cutoff group:
796 +                                                                                        dag = fDecomp_->getAtomToGroupVectorRow(atom1, (*cg1)->getGlobalIndex());
797 +                                                                                        tauLcl -= outProduct(dag, fg);
798 +                                                                                }
799 +                                                                        }
800 +                                                                }
801 +                                                                for (vector<int>::iterator jb = atomListColumn.begin(); jb != atomListColumn.end(); ++jb)
802 +                                                                {
803 +                                                                        atom2 = (*jb);
804 +                                                                        mf = fDecomp_->getMassFactorColumn(atom2);
805 +                                                                        // fg is the force on atom jb due to cutoff group's
806 +                                                                        // presence in switching region
807 +                                                                        fg = -swderiv * d_grp * mf;
808 +                                                                        #pragma omp critical (forceLck)
809 +                                                                        {
810 +                                                                                fDecomp_->addForceToAtomColumn(atom2, fg);
811 +                                                                        }
812 +
813 +                                                                        if (atomListColumn.size() > 1)
814 +                                                                        {
815 +                                                                                if (info_->usesAtomicVirial())
816 +                                                                                {
817 +                                                                                        // find the distance between the atom
818 +                                                                                        // and the center of the cutoff group:
819 +                                                                                        dag = fDecomp_->getAtomToGroupVectorColumn(atom2, (*cg2)->getGlobalIndex());
820 +                                                                                        tauLcl -= outProduct(dag, fg);
821 +                                                                                }
822 +                                                                        }
823 +                                                                }
824 +                                                        }
825 +                                                        //if (!SIM_uses_AtomicVirial) {
826 +                                                        //  tau -= outProduct(d_grp, fij);
827 +                                                        //}
828 +                                                }
829 +                                        }
830 +                                }
831 +                        }// END: omp for loop
832 +                        /* Critical region which accumulates forces from local to shared arrays  */
833 +                        #pragma omp critical (forceAdd)
834 +                        {
835 +                                for (int currAtom = 0; currAtom < config->force.size(); ++currAtom)
836 +                                {
837 +                                        config->force[currAtom] += forceLcl[currAtom];
838 +                                }
839 +
840 +                                tau -= tauLcl;
841 +                                *(fDecomp_->getPairwisePotential()) += potLcl;
842 +                        }
843 +                }// END: omp parallel region
844 +
845 +                /*gettimeofday(&tv2, NULL);
846 +
847 +                elapsedTime += 1000000 * (tv2.tv_sec - tv.tv_sec)
848 +                                                + (tv2.tv_usec - tv.tv_usec);
849 +
850 +                Globals *simParams_ = info_->getSimParams();
851 +
852 +                if(curSnapshot->getTime() >= simParams_->getRunTime() - 1)
853 +                {
854 +                        printf("Force calculation time: %ld [us]\n", elapsedTime);
855 +                }*/
856 +
857 +                if (iLoop == PREPAIR_LOOP)
858 +                {
859 +                        if (info_->requiresPrepair())
860 +                        {
861 +
862 +                                fDecomp_->collectIntermediateData();
863 +
864 +                                for (int atom1 = 0; atom1 < info_->getNAtoms(); atom1++)
865 +                                {
866 +                                        fDecomp_->fillSelfData(sdat, atom1);
867 +                                        interactionMan_->doPreForce(sdat);
868 +                                }
869 +
870 +                                fDecomp_->distributeIntermediateData();
871 +
872 +                        }
873 +                }
874 +        }
875 +
876 +        fDecomp_->collectData();
877 +
878 +        if (info_->requiresSelfCorrection())
879 +        {
880 +
881 +                for (int atom1 = 0; atom1 < info_->getNAtoms(); atom1++)
882 +                {
883 +                        fDecomp_->fillSelfData(sdat, atom1);
884 +                        interactionMan_->doSelfCorrection(sdat);
885 +                }
886 +
887 +        }
888 +
889 +        longRangePotential = *(fDecomp_->getEmbeddingPotential()) + *(fDecomp_->getPairwisePotential());
890 +
891 +        lrPot = longRangePotential.sum();
892 +
893 +        //store the tau and long range potential
894 +        curSnapshot->statData[Stats::LONG_RANGE_POTENTIAL] = lrPot;
895 +        curSnapshot->statData[Stats::VANDERWAALS_POTENTIAL] = longRangePotential[VANDERWAALS_FAMILY];
896 +        curSnapshot->statData[Stats::ELECTROSTATIC_POTENTIAL] = longRangePotential[ELECTROSTATIC_FAMILY];
897 + }
898 +
899   void ForceManager::longRangeInteractionsRapaport() {
900  
901          Snapshot* curSnapshot = info_->getSnapshotManager()->getCurrentSnapshot();
# Line 549 | Line 915 | void ForceManager::longRangeInteractionsRapaport() {
915                  {
916                          for (cg = mol->beginCutoffGroup(ci); cg != NULL; cg = mol->nextCutoffGroup(ci))
917                          {
918 < //                              cerr << "branch1\n";
919 < //                              cerr << "globind = " << cg->getGlobalIndex() << ":" << __LINE__ << "\n";
918 >                                //                              cerr << "branch1\n";
919 >                                //                              cerr << "globind = " << cg->getGlobalIndex() << ":" << __LINE__ << "\n";
920                                  cg->updateCOM();
921  
922 < //                              cerr << "gbI: " << cg->getGlobalIndex() << " locI: " << cg->getLocalIndex() << " x: "
923 < //                                              << cgConfig->position[cg->getLocalIndex()].x() << " y: " << cgConfig->position[cg->getLocalIndex()].y()
924 < //                                              << " z: " << cgConfig->position[cg->getLocalIndex()].z() << "\n";
922 >                                //                              cerr << "gbI: " << cg->getGlobalIndex() << " locI: " << cg->getLocalIndex() << " x: "
923 >                                //                                              << cgConfig->position[cg->getLocalIndex()].x() << " y: " << cgConfig->position[cg->getLocalIndex()].y()
924 >                                //                                              << " z: " << cgConfig->position[cg->getLocalIndex()].z() << "\n";
925                          }
926                  }
927          } else
928          {
929                  // center of mass of the group is the same as position of the atom
930                  // if cutoff group does not exist
931 < //              cerr << ":" << __LINE__ << "branch2\n";
931 >                //              cerr << ":" << __LINE__ << "branch2\n";
932                  cgConfig->position = config->position;
933          }
934  
# Line 620 | Line 986 | void ForceManager::longRangeInteractionsRapaport() {
986                                  neighborMatW = fDecomp_->buildLayerBasedNeighborList();
987                  }
988  
623                //              printf("before omp loop\n");
624                //#pragma omp parallel for num_threads(3) default(none) shared(curSnapshot, idat, iLoop, sw, cerr) \
625        private(i, j, cg1, cg2, cuts, d_grp, rgrpsq, rCutSq, vij, fij, in_switching_region, rgrp, dswdr, atomListRow, atomListColumn, atom1, atom2, vpair, workPot, f1, topoDist, vdwMult, electroMult, d, r2, r, swderiv, fg, mf, dag)
989                  for (mol = info_->beginMolecule(mi); mol != NULL; mol = info_->nextMolecule(mi))
990                  {
991                          for (cg1 = mol->beginCutoffGroup(ci); cg1 != NULL; cg1 = mol->nextCutoffGroup(ci))
992                          {
993                                  //                      printf("Thread %d executes loop iteration %d\n", omp_get_thread_num(), i);
994 <                                for (vector<CutoffGroup *>::iterator cg2 = neighborMatW[cg1->getGlobalIndex()].begin(); cg2 != neighborMatW[cg1->getGlobalIndex()].end(); ++cg2)
994 >                                for (vector<CutoffGroup *>::iterator cg2 = neighborMatW[cg1->getGlobalIndex()].begin(); cg2
995 >                                                != neighborMatW[cg1->getGlobalIndex()].end(); ++cg2)
996                                  {
997  
998                                          cuts = fDecomp_->getGroupCutoffs(cg1->getGlobalIndex(), (*cg2)->getGlobalIndex());
# Line 677 | Line 1041 | void ForceManager::longRangeInteractionsRapaport() {
1041                                                                          {
1042                                                                                  idat.d = &d_grp;
1043                                                                                  idat.r2 = &rgrpsq;
1044 < //                                                                              cerr << "dgrp = " << d_grp << ":" << __LINE__ << "\n";
1044 >                                                                                //                                                                              cerr << "dgrp = " << d_grp << ":" << __LINE__ << "\n";
1045                                                                          } else
1046                                                                          {
1047                                                                                  d = fDecomp_->getInteratomicVector(atom1, atom2);
1048                                                                                  curSnapshot->wrapVector(d);
1049                                                                                  r2 = d.lengthSquare();
1050 < //                                                                              cerr << "datm = " << d << ":" << __LINE__ << "\n";
1050 >                                                                                //                                                                              cerr << "datm = " << d << ":" << __LINE__ << "\n";
1051                                                                                  idat.d = &d;
1052                                                                                  idat.r2 = &r2;
1053                                                                          }
1054  
1055 < //                                                                      cerr << "idat.d = " << *(idat.d) << ":" << __LINE__ << "\n";
1055 >                                                                        //                                                                      cerr << "idat.d = " << *(idat.d) << ":" << __LINE__ << "\n";
1056                                                                          r = sqrt(*(idat.r2));
1057                                                                          idat.rij = &r;
1058 < //                                                                      cerr << "idat.rij = " << *(idat.rij) << "\n";
1058 >                                                                        //                                                                      cerr << "idat.rij = " << *(idat.rij) << "\n";
1059  
1060                                                                          if (iLoop == PREPAIR_LOOP)
1061                                                                          {
# Line 701 | Line 1065 | void ForceManager::longRangeInteractionsRapaport() {
1065                                                                                  interactionMan_->doPair(idat);
1066                                                                                  fDecomp_->unpackInteractionData(idat, atom1, atom2);
1067  
1068 < //                                                                              cerr << "d = " << *(idat.d) << "\tv=" << vpair << "\tf=" << f1 << ":" << __LINE__ << "\n";
1068 >                                                                                //                                                                              cerr << "d = " << *(idat.d) << "\tv=" << vpair << "\tf=" << f1 << ":" << __LINE__ << "\n";
1069                                                                                  vij += vpair;
1070                                                                                  fij += f1;
1071                                                                                  tau -= outProduct(*(idat.d), f1);
# Line 771 | Line 1135 | void ForceManager::longRangeInteractionsRapaport() {
1135                                          }
1136                                  }
1137                          }
1138 <                }// END: omp for loop
775 <                //              printf("after omp loop\n");
1138 >                }
1139  
1140                  if (iLoop == PREPAIR_LOOP)
1141                  {

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines