[ViewVC] Diff of: OpenMD/branches/devel_omp/src/parallel/ForceMatrixDecomposition.cpp

Comparing branches/devel_omp/src/parallel/ForceMatrixDecomposition.cpp (file contents):
Revision 1594 by chuckv, Tue Jul 19 16:45:30 2011 UTC vs.
Revision 1595 by chuckv, Tue Jul 19 18:50:04 2011 UTC

+#endif
-+
-+
+vector<vector<int> > ForceMatrixDecomposition::buildLayerBasedNeighborList() {
-+
+        printf("buildLayerBasedNeighborList; nGroups:%d\n", nGroups_);
-+
+        // Na = nGroups_
-+
+        /* cell occupancy counter */
-+
+        vector<int> k_c;
-+
+        /* c_i - has cell containing atom i (size Na) */
-+
+        vector<int> c;
-+
+        /* l_i - layer containing atom i (size Na) */
-+
+        vector<int> l;
-+
-+
+//      cellList_.clear();
-+
-+
+        RealType rList_ = (largestRcut_ + skinThickness_);
-+
+        Snapshot* snap_ = sman_->getCurrentSnapshot();
-+
+        Mat3x3d Hmat = snap_->getHmat();
-+
+        Vector3d Hx = Hmat.getColumn(0);
-+
+        Vector3d Hy = Hmat.getColumn(1);
-+
+        Vector3d Hz = Hmat.getColumn(2);
-+
-+
+        nCells_.x() = (int) (Hx.length()) / rList_;
-+
+        nCells_.y() = (int) (Hy.length()) / rList_;
-+
+        nCells_.z() = (int) (Hz.length()) / rList_;
-+
-+
+        Mat3x3d invHmat = snap_->getInvHmat();
-+
+        Vector3d rs, scaled, dr;
-+
+        Vector3i whichCell;
-+
+        int cellIndex;
-+
+        int nCtot = nCells_.x() * nCells_.y() * nCells_.z();
-+
-+
+//      cellList_.resize(nCtot);
-+
+        k_c = vector<int>(nCtot, 0);
-+
-+
+        for (int i = 0; i < nGroups_; i++)
-+
+        {
-+
+                rs = snap_->cgData.position[i];
-+
-+
+                // scaled positions relative to the box vectors
-+
+                scaled = invHmat * rs;
-+
-+
+                // wrap the vector back into the unit box by subtracting integer box
-+
+                // numbers
-+
+                for (int j = 0; j < 3; j++)
-+
+                {
-+
+                        scaled[j] -= roundMe(scaled[j]);
-+
+                        scaled[j] += 0.5;
-+
+                }
-+
-+
+                // find xyz-indices of cell that cutoffGroup is in.
-+
+                whichCell.x() = nCells_.x() * scaled.x();
-+
+                whichCell.y() = nCells_.y() * scaled.y();
-+
+                whichCell.z() = nCells_.z() * scaled.z();
-+
-+
+                // find single index of this cell:
-+
+                cellIndex = Vlinear(whichCell, nCells_);
-+
-+
+                c.push_back(cellIndex);
-+
-+
+//              // add this cutoff group to the list of groups in this cell;
-+
+//              cellList_[cellIndex].push_back(i);
-+
+        }
-+
-+
+        int k_c_curr;
-+
+        int k_c_max = 0;
-+
+        /* the cell-layer occupancy matrix */
-+
+        vector<vector<int> > H_c_l = vector<vector<int> >(nCtot);
-+
-+
+        for(int i = 0; i < nGroups_; ++i)
-+
+        {
-+
+                k_c_curr = ++k_c[c[i]];
-+
+                l.push_back(k_c_curr);
-+
-+
+                /* determines the number of layers in use */
-+
+                if(k_c_max < k_c_curr)
-+
+                {
-+
+                        k_c_max = k_c_curr;
-+
+                }
-+
-+
+                H_c_l[c[i]].push_back(/*l[*/i/*]*/);
-+
+        }
-+
+        int m;
-+
+        /* the neighbor matrix */
-+
+        vector<vector<int> >neighborMatW = vector<vector<int> >(nGroups_);
-+
-+
+//      vector<pair<int, int> > neighborList;
-+
+        groupCutoffs cuts;
-+
-+
+        /* loops over objects(atoms, rigidBodies, cutoffGroups, etc.) */
-+
+        for(int i = 0; i < nGroups_; ++i)
-+
+        {
-+
+                m = 0;
-+
+                /* c' */
-+
+                int c1 = c[i];
-+
+                Vector3i c1v = idxToV(c1, nCells_);
-+
-+
+                /* loops over the neighboring cells c'' */
-+
+                for (vector<Vector3i>::iterator os = cellOffsets_.begin(); os != cellOffsets_.end(); ++os)
-+
+                {
-+
+                        Vector3i c2v = c1v + (*os);
-+
-+
+                        if (c2v.x() >= nCells_.x())
-+
+                        {
-+
+                                c2v.x() = 0;
-+
+                        } else if (c2v.x() < 0)
-+
+                        {
-+
+                                c2v.x() = nCells_.x() - 1;
-+
+                        }
-+
-+
+                        if (c2v.y() >= nCells_.y())
-+
+                        {
-+
+                                c2v.y() = 0;
-+
+                        } else if (c2v.y() < 0)
-+
+                        {
-+
+                                c2v.y() = nCells_.y() - 1;
-+
+                        }
-+
-+
+                        if (c2v.z() >= nCells_.z())
-+
+                        {
-+
+                                c2v.z() = 0;
-+
+                        } else if (c2v.z() < 0)
-+
+                        {
-+
+                                c2v.z() = nCells_.z() - 1;
-+
+                        }
-+
-+
+                        int c2 = Vlinear(c2v, nCells_);
-+
+                        /* loops over layers l to access the neighbor atoms */
-+
+                        for (vector<int>::iterator j = H_c_l[c2].begin(); j != H_c_l[c2].end(); ++j)
-+
+                        {
-+
+//                              if i'' = 0 then break // doesn't apply to vector implementation of matrix
-+
+//                              if(i != *j)
-+
+                                if (c2 != c1 || (*j) < (i))
-+
+                                {
-+
+                                        dr = snap_->cgData.position[(*j)] - snap_->cgData.position[(i)];
-+
+                                        snap_->wrapVector(dr);
-+
+                                        cuts = getGroupCutoffs((i), (*j));
-+
+                                        if (dr.lengthSquare() < cuts.third)
-+
+                                        {
-+
+                                                ++m;
-+
+                                                /* transposed version of Rapaport W mat, to occupy successive memory locations on CPU */
-+
+                                                neighborMatW[i].push_back(*j);
-+
+//                                              neighborList.push_back(make_pair((i), (*j)));
-+
+                                        }
-+
+                                }
-+
+                        }
-+
+                }
-+
+        }
-+
-+
+        // save the local cutoff group positions for the check that is
-+
+        // done on each loop:
-+
+        saved_CG_positions_.clear();
-+
+        for (int i = 0; i < nGroups_; i++)
-+
+                saved_CG_positions_.push_back(snap_->cgData.position[i]);
-+
-+
+        return neighborMatW;
-+
+}
-+
+  /*
+   * buildNeighborList

Diff Legend

-–
+Removed lines
-+
+Added lines
-<
+Changed lines
->
+Changed lines

Comparing branches/devel_omp/src/parallel/ForceMatrixDecomposition.cpp (file contents): Revision 1594 by chuckv, Tue Jul 19 16:45:30 2011 UTC vs. Revision 1595 by chuckv, Tue Jul 19 18:50:04 2011 UTC

Diff Legend

Comparing branches/devel_omp/src/parallel/ForceMatrixDecomposition.cpp (file contents):
Revision 1594 by chuckv, Tue Jul 19 16:45:30 2011 UTC vs.
Revision 1595 by chuckv, Tue Jul 19 18:50:04 2011 UTC