--- branches/development/src/parallel/ForceDecomposition.hpp	2011/01/17 21:34:36	1540
+++ branches/development/src/parallel/ForceMatrixDecomposition.hpp	2011/05/24 21:24:45	1567
@@ -39,47 +39,114 @@
  * [4]  Vardeman & Gezelter, in progress (2009).                        
  */
  
-#ifndef PARALLEL_FORCEDECOMPOSITION_HPP
-#define PARALLEL_FORCEDECOMPOSITION_HPP
+#ifndef PARALLEL_FORCEMATRIXDECOMPOSITION_HPP
+#define PARALLEL_FORCEMATRIXDECOMPOSITION_HPP
 
-#include "Parallel/Decomposition.hpp"
-#include "Parallel/Communicator.hpp"
+#include "parallel/ForceDecomposition.hpp"
 #include "math/SquareMatrix3.hpp"
+#include "brains/Snapshot.hpp"
 
+#ifdef IS_MPI
+#include "parallel/Communicator.hpp"
+#endif
+
+using namespace std;
 namespace OpenMD {
   
-  class ForceDecomposition : public Decomposition {
+  class ForceMatrixDecomposition : public ForceDecomposition {
   public:
-    ForceDecomposition(Snapshot* sman);
+    ForceMatrixDecomposition(SimInfo* info) : ForceDecomposition(info) {sman_ = info_->getSnapshotManager();}
     void distributeInitialData();
     void distributeData();
     void collectIntermediateData();
     void distributeIntermediateData();
     void collectData();
 
-    unsigned int getNcutoffGroupsI();
-    unsigned int getNcutoffGroupsJ();
+    // neighbor list routines
+    bool checkNeighborList();
+    vector<pair<int, int> >  buildNeighborList();
 
-    vector<int> getAtomsInGroupI(int whichCGI);
-    vector<int> getAtomsInGroupJ(int whichCGJ);
+    // group bookkeeping
+    pair<int, int> getGroupTypes(int cg1, int cg2);
 
-    AtomType* getAtomTypeI(int whichAtomI);
-    AtomType* getAtomTypeJ(int whichAtomJ);   
+    // Group->atom bookkeeping
+    vector<int> getAtomsInGroupRow(int cg1);
+    vector<int> getAtomsInGroupColumn(int cg2);
+    Vector3d getAtomToGroupVectorRow(int atom1, int cg1);
+    Vector3d getAtomToGroupVectorColumn(int atom2, int cg2);
+    RealType getMfactRow(int atom1);
+    RealType getMfactColumn(int atom2);
 
-#ifdef IS_MPI
-    
-    Communicator<Row, RealType>* AtomCommRealI; 
-    Communicator<Row, Vector3d>* AtomCommVectorI; 
-    Communicator<Row, Mat3x3d>*  AtomCommMatrixI; 
+    // spatial data
+    Vector3d getIntergroupVector(int cg1, int cg2);
+    Vector3d getInteratomicVector(int atom1, int atom2);
+       
+    // atom bookkeeping
+    vector<int> getAtomList();
+    vector<int> getSkipsForAtom(int atom1);
+    bool skipAtomPair(int atom1, int atom2);
+    void addForceToAtomRow(int atom1, Vector3d fg);
+    void addForceToAtomColumn(int atom2, Vector3d fg);
 
-    Communicator<Column, RealType>* AtomCommRealJ; 
-    Communicator<Column, Vector3d>* AtomCommVectorJ; 
-    Communicator<Column, Mat3x3d>*  AtomCommMatrixJ; 
+    // filling interaction blocks with pointers
+    InteractionData fillInteractionData(int atom1, int atom2);
+    InteractionData fillSkipData(int atom1, int atom2);
+    SelfData fillSelfData(int atom1);
 
-    Communicator<Row, Vector3d>* cgCommVectorI; 
-    Communicator<Column, Vector3d>* cgCommVectorJ; 
-    
+  private: 
+    SnapshotManager* sman_;    
+    Snapshot* snap_;
+    int storageLayout_;
+    vector<Vector3i> Cells;
+
+    int nLocal_;
+    int nGroups_;
+
+#ifdef IS_MPI    
+    DataStorage atomRowData;
+    DataStorage atomColData;
+    DataStorage cgRowData;
+    DataStorage cgColData;
+
+    int nAtomsInRow_; 
+    int nAtomsInCol_;
+    int nGroupsInRow_;
+    int nGroupsInCol_;
+
+    Communicator<Row, int>* AtomCommIntRow;
+    Communicator<Row, RealType>* AtomCommRealRow; 
+    Communicator<Row, Vector3d>* AtomCommVectorRow; 
+    Communicator<Row, Mat3x3d>*  AtomCommMatrixRow; 
+
+    Communicator<Column, int>* AtomCommIntColumn;
+    Communicator<Column, RealType>* AtomCommRealColumn; 
+    Communicator<Column, Vector3d>* AtomCommVectorColumn; 
+    Communicator<Column, Mat3x3d>*  AtomCommMatrixColumn; 
+
+    Communicator<Row, int>* cgCommIntRow;
+    Communicator<Row, Vector3d>* cgCommVectorRow; 
+    Communicator<Column, int>* cgCommIntColumn;
+    Communicator<Column, Vector3d>* cgCommVectorColumn; 
+
+    vector<vector<RealType> > pot_row;
+    vector<vector<RealType> > pot_col;
+    vector<int> identsRow;
+    vector<int> identsCol;
+
+    vector<int> AtomLocalToGlobal;
+    vector<int> AtomRowToGlobal;
+    vector<int> AtomColToGlobal;
+    vector<int> cgLocalToGlobal;
+    vector<int> cgRowToGlobal;
+    vector<int> cgColToGlobal;
+
+    vector<vector<int> > CellListRow;
+    vector<vector<int> > CellListCol;
+#else 
+    vector<vector<int> > CellList;
 #endif
+    vector<RealType> pot_local;
+
   };
 
 }