Bug fix in FMM_Pts, FMM_Tree

Vector for [up/dn]wd_[equiv/check]_surface coordinates was a member of
FMM_Pts, however the vector points to memory in setup_data which belongs
to FMM_Tree.  This caused access to freed memory if the tree was
destroyed.  The fix is to move this coordinate data to the tree.

include/fmm_cheb.hpp

@@ -28,6 +28,7 @@ class FMM_Cheb: public FMM_Pts<FMMNode>{
 
  public:
 
+  typedef typename FMM_Pts<FMMNode>::FMMTree_t FMMTree_t;
   typedef typename FMMNode::Real_t Real_t;
   typedef FMMNode FMMNode_t;
 
@@ -69,37 +70,37 @@ class FMM_Cheb: public FMM_Pts<FMMNode>{
    */
   int& ChebDeg(){return cheb_deg;}
 
-  virtual void CollectNodeData(std::vector<FMMNode*>& nodes, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, std::vector<std::vector<Vector<Real_t>* > > vec_list = std::vector<std::vector<Vector<Real_t>* > >(0));
+  virtual void CollectNodeData(FMMTree_t* tree, std::vector<FMMNode*>& nodes, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, std::vector<std::vector<Vector<Real_t>* > > vec_list = std::vector<std::vector<Vector<Real_t>* > >(0));
 
   /**
    * \brief Initialize multipole expansions for the given array of leaf nodes
    * at a given level.
    */
-  virtual void Source2UpSetup(SetupData<Real_t>& setup_data, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
+  virtual void Source2UpSetup(SetupData<Real_t>& setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
   virtual void Source2Up     (SetupData<Real_t>& setup_data, bool device=false);
 
   /**
    * \brief Compute X-List intractions.
    */
-  virtual void X_ListSetup(SetupData<Real_t>& setup_data, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
+  virtual void X_ListSetup(SetupData<Real_t>& setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
   virtual void X_List     (SetupData<Real_t>& setup_data, bool device=false);
 
   /**
    * \brief Compute target potential from the local expansion.
    */
-  virtual void Down2TargetSetup(SetupData<Real_t>& setup_data, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
+  virtual void Down2TargetSetup(SetupData<Real_t>& setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
   virtual void Down2Target     (SetupData<Real_t>& setup_data, bool device=false);
 
   /**
    * \brief Compute W-List intractions.
    */
-  virtual void W_ListSetup(SetupData<Real_t>& setup_data, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
+  virtual void W_ListSetup(SetupData<Real_t>& setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
   virtual void W_List     (SetupData<Real_t>& setup_data, bool device=false);
 
   /**
    * \brief Compute U-List intractions.
    */
-  virtual void U_ListSetup(SetupData<Real_t>& setup_data, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
+  virtual void U_ListSetup(SetupData<Real_t>& setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
   virtual void U_List     (SetupData<Real_t>& setup_data, bool device=false);
 
   virtual void PostProcessing(std::vector<FMMNode_t*>& nodes);

include/fmm_cheb.txx

@@ -817,7 +817,7 @@ Matrix<typename FMMNode::Real_t>& FMM_Cheb<FMMNode>::Precomp(int level, Mat_Type
 
 
 template <class FMMNode>
-void FMM_Cheb<FMMNode>::CollectNodeData(std::vector<FMMNode*>& node, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, std::vector<std::vector<Vector<Real_t>* > > vec_list){
+void FMM_Cheb<FMMNode>::CollectNodeData(FMMTree_t* tree, std::vector<FMMNode*>& node, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, std::vector<std::vector<Vector<Real_t>* > > vec_list){
   if(vec_list.size()<6) vec_list.resize(6);
   size_t n_coeff=(cheb_deg+1)*(cheb_deg+2)*(cheb_deg+3)/6;
   if(node.size()==0) return;
@@ -843,14 +843,14 @@ void FMM_Cheb<FMMNode>::CollectNodeData(std::vector<FMMNode*>& node, std::vector
       }
     }
   }
-  FMM_Pts<FMMNode_t>::CollectNodeData(node, buff, n_list, vec_list);
+  FMM_Pts<FMMNode_t>::CollectNodeData(tree, node, buff, n_list, vec_list);
 }
 
 
 template <class FMMNode>
-void FMM_Cheb<FMMNode>::Source2UpSetup(SetupData<Real_t>& setup_data, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
+void FMM_Cheb<FMMNode>::Source2UpSetup(SetupData<Real_t>& setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
   if(this->MultipoleOrder()==0) return;
-  FMM_Pts<FMMNode>::Source2UpSetup(setup_data, buff, n_list, level, device);
+  FMM_Pts<FMMNode>::Source2UpSetup(setup_data, tree, buff, n_list, level, device);
 
   { // Set setup_data
     setup_data.level=level;
@@ -887,9 +887,9 @@ void FMM_Cheb<FMMNode>::Source2Up     (SetupData<Real_t>& setup_data, bool devic
 
 
 template <class FMMNode>
-void FMM_Cheb<FMMNode>::X_ListSetup(SetupData<Real_t>& setup_data, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
+void FMM_Cheb<FMMNode>::X_ListSetup(SetupData<Real_t>& setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
   if(this->MultipoleOrder()==0) return;
-  FMM_Pts<FMMNode>::X_ListSetup(setup_data, buff, n_list, level, device);
+  FMM_Pts<FMMNode>::X_ListSetup(setup_data, tree, buff, n_list, level, device);
 
   { // Set setup_data
     setup_data.level=level;
@@ -932,7 +932,7 @@ void FMM_Cheb<FMMNode>::X_List     (SetupData<Real_t>& setup_data, bool device){
 
 
 template <class FMMNode>
-void FMM_Cheb<FMMNode>::W_ListSetup(SetupData<Real_t>& setup_data, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
+void FMM_Cheb<FMMNode>::W_ListSetup(SetupData<Real_t>& setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
   if(this->MultipoleOrder()==0) return;
   { // Set setup_data
     setup_data.level=level;
@@ -974,8 +974,8 @@ void FMM_Cheb<FMMNode>::W_List     (SetupData<Real_t>& setup_data, bool device){
 
 
 template <class FMMNode>
-void FMM_Cheb<FMMNode>::U_ListSetup(SetupData<Real_t>& setup_data, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
-  FMM_Pts<FMMNode>::U_ListSetup(setup_data, buff, n_list, level, device);
+void FMM_Cheb<FMMNode>::U_ListSetup(SetupData<Real_t>& setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
+  FMM_Pts<FMMNode>::U_ListSetup(setup_data, tree, buff, n_list, level, device);
 
   { // Set setup_data
     setup_data.level=level;
@@ -1020,7 +1020,7 @@ void FMM_Cheb<FMMNode>::U_List     (SetupData<Real_t>& setup_data, bool device){
 
 
 template <class FMMNode>
-void FMM_Cheb<FMMNode>::Down2TargetSetup(SetupData<Real_t>& setup_data, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
+void FMM_Cheb<FMMNode>::Down2TargetSetup(SetupData<Real_t>& setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
   if(this->MultipoleOrder()==0) return;
   { // Set setup_data
     setup_data.level=level;

include/fmm_pts.hpp

@@ -86,12 +86,15 @@ struct SetupData{
   Matrix<Real_t>* output_data;
 };
 
+template <class FMM_Mat_t>
+class FMM_Tree;
 
 template <class FMMNode>
 class FMM_Pts{
 
  public:
 
+  typedef FMM_Tree<FMM_Pts<FMMNode> > FMMTree_t;
   typedef typename FMMNode::Real_t Real_t;
   typedef FMMNode FMMNode_t;
 
@@ -133,7 +136,7 @@ class FMM_Pts{
    */
   bool Homogen(){return kernel->homogen;}
 
-  virtual void CollectNodeData(std::vector<FMMNode*>& nodes, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, std::vector<std::vector<Vector<Real_t>* > > vec_list = std::vector<std::vector<Vector<Real_t>* > >(0));
+  virtual void CollectNodeData(FMMTree_t* tree, std::vector<FMMNode*>& nodes, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, std::vector<std::vector<Vector<Real_t>* > > vec_list = std::vector<std::vector<Vector<Real_t>* > >(0));
 
   void SetupPrecomp(SetupData<Real_t>& setup_data, bool device=false);
   void SetupInterac(SetupData<Real_t>& setup_data, bool device=false);
@@ -148,14 +151,14 @@ class FMM_Pts{
    * \brief Initialize multipole expansions for the given array of leaf nodes
    * at a given level.
    */
-  virtual void Source2UpSetup(SetupData<Real_t>&  setup_data, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
+  virtual void Source2UpSetup(SetupData<Real_t>&  setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
   virtual void Source2Up     (SetupData<Real_t>&  setup_data, bool device=false);
 
   /**
    * \brief Initialize multipole expansions for the given array of non-leaf
    * nodes from that of its children.
    */
-  virtual void Up2UpSetup(SetupData<Real_t>&  setup_data, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
+  virtual void Up2UpSetup(SetupData<Real_t>&  setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
   virtual void Up2Up     (SetupData<Real_t>&  setup_data, bool device=false);
 
   virtual void PeriodicBC(FMMNode* node);
@@ -163,37 +166,37 @@ class FMM_Pts{
   /**
    * \brief Compute V-List intractions.
    */
-  virtual void V_ListSetup(SetupData<Real_t>&  setup_data, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
+  virtual void V_ListSetup(SetupData<Real_t>&  setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
   virtual void V_List     (SetupData<Real_t>&  setup_data, bool device=false);
 
   /**
    * \brief Compute X-List intractions.
    */
-  virtual void X_ListSetup(SetupData<Real_t>&  setup_data, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
+  virtual void X_ListSetup(SetupData<Real_t>&  setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
   virtual void X_List     (SetupData<Real_t>&  setup_data, bool device=false);
 
   /**
    * \brief Compute contribution of local expansion from the parent.
    */
-  virtual void Down2DownSetup(SetupData<Real_t>&  setup_data, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
+  virtual void Down2DownSetup(SetupData<Real_t>&  setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
   virtual void Down2Down     (SetupData<Real_t>&  setup_data, bool device=false);
 
   /**
    * \brief Compute target potential from the local expansion.
    */
-  virtual void Down2TargetSetup(SetupData<Real_t>&  setup_data, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
+  virtual void Down2TargetSetup(SetupData<Real_t>&  setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
   virtual void Down2Target     (SetupData<Real_t>&  setup_data, bool device=false);
 
   /**
    * \brief Compute W-List intractions.
    */
-  virtual void W_ListSetup(SetupData<Real_t>&  setup_data, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
+  virtual void W_ListSetup(SetupData<Real_t>&  setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
   virtual void W_List     (SetupData<Real_t>&  setup_data, bool device=false);
 
   /**
    * \brief Compute U-List intractions.
    */
-  virtual void U_ListSetup(SetupData<Real_t>&  setup_data, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
+  virtual void U_ListSetup(SetupData<Real_t>&  setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& node_data, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device);
   virtual void U_List     (SetupData<Real_t>&  setup_data, bool device=false);
 
   virtual void PostProcessing(std::vector<FMMNode_t*>& nodes);
@@ -223,11 +226,6 @@ class FMM_Pts{
       Vector<Real_t>& input_data, Vector<Real_t>& output_data, Vector<Real_t>& buffer_, Matrix<Real_t>& M);
   typename FFTW_t<Real_t>::plan vlist_ifftplan; bool vlist_ifft_flag;
 
-  std::vector<Vector<Real_t> > upwd_check_surf;
-  std::vector<Vector<Real_t> > upwd_equiv_surf;
-  std::vector<Vector<Real_t> > dnwd_check_surf;
-  std::vector<Vector<Real_t> > dnwd_equiv_surf;
-
   mem::MemoryManager* mem_mgr;
   InteracList<FMMNode> interac_list;
   const Kernel<Real_t>* kernel;    //The kernel function.

include/fmm_pts.txx

@@ -1025,7 +1025,7 @@ void FMM_Pts<FMMNode>::PrecompAll(Mat_Type type, int level){
 }
 
 template <class FMMNode>
-void FMM_Pts<FMMNode>::CollectNodeData(std::vector<FMMNode*>& node, std::vector<Matrix<Real_t> >& buff_list, std::vector<Vector<FMMNode_t*> >& n_list, std::vector<std::vector<Vector<Real_t>* > > vec_list){
+void FMM_Pts<FMMNode>::CollectNodeData(FMMTree_t* tree, std::vector<FMMNode*>& node, std::vector<Matrix<Real_t> >& buff_list, std::vector<Vector<FMMNode_t*> >& n_list, std::vector<std::vector<Vector<Real_t>* > > vec_list){
   if(buff_list.size()<7) buff_list.resize(7);
   if(   n_list.size()<7)    n_list.resize(7);
   if( vec_list.size()<7)  vec_list.resize(7);
@@ -1221,29 +1221,29 @@ void FMM_Pts<FMMNode>::CollectNodeData(std::vector<FMMNode*>& node, std::vector<
       }
     }
     { // check and equiv surfaces.
-      if(upwd_check_surf.size()==0){
+      if(tree->upwd_check_surf.size()==0){
         size_t m=MultipoleOrder();
-        upwd_check_surf.resize(MAX_DEPTH);
-        upwd_equiv_surf.resize(MAX_DEPTH);
-        dnwd_check_surf.resize(MAX_DEPTH);
-        dnwd_equiv_surf.resize(MAX_DEPTH);
+        tree->upwd_check_surf.resize(MAX_DEPTH);
+        tree->upwd_equiv_surf.resize(MAX_DEPTH);
+        tree->dnwd_check_surf.resize(MAX_DEPTH);
+        tree->dnwd_equiv_surf.resize(MAX_DEPTH);
         for(size_t depth=0;depth<MAX_DEPTH;depth++){
           Real_t c[3]={0.0,0.0,0.0};
-          upwd_check_surf[depth].ReInit((6*(m-1)*(m-1)+2)*COORD_DIM);
-          upwd_equiv_surf[depth].ReInit((6*(m-1)*(m-1)+2)*COORD_DIM);
-          dnwd_check_surf[depth].ReInit((6*(m-1)*(m-1)+2)*COORD_DIM);
-          dnwd_equiv_surf[depth].ReInit((6*(m-1)*(m-1)+2)*COORD_DIM);
-          upwd_check_surf[depth]=u_check_surf(m,c,depth);
-          upwd_equiv_surf[depth]=u_equiv_surf(m,c,depth);
-          dnwd_check_surf[depth]=d_check_surf(m,c,depth);
-          dnwd_equiv_surf[depth]=d_equiv_surf(m,c,depth);
+          tree->upwd_check_surf[depth].ReInit((6*(m-1)*(m-1)+2)*COORD_DIM);
+          tree->upwd_equiv_surf[depth].ReInit((6*(m-1)*(m-1)+2)*COORD_DIM);
+          tree->dnwd_check_surf[depth].ReInit((6*(m-1)*(m-1)+2)*COORD_DIM);
+          tree->dnwd_equiv_surf[depth].ReInit((6*(m-1)*(m-1)+2)*COORD_DIM);
+          tree->upwd_check_surf[depth]=u_check_surf(m,c,depth);
+          tree->upwd_equiv_surf[depth]=u_equiv_surf(m,c,depth);
+          tree->dnwd_check_surf[depth]=d_check_surf(m,c,depth);
+          tree->dnwd_equiv_surf[depth]=d_equiv_surf(m,c,depth);
         }
       }
       for(size_t depth=0;depth<MAX_DEPTH;depth++){
-        vec_lst.push_back(&upwd_check_surf[depth]);
-        vec_lst.push_back(&upwd_equiv_surf[depth]);
-        vec_lst.push_back(&dnwd_check_surf[depth]);
-        vec_lst.push_back(&dnwd_equiv_surf[depth]);
+        vec_lst.push_back(&tree->upwd_check_surf[depth]);
+        vec_lst.push_back(&tree->upwd_equiv_surf[depth]);
+        vec_lst.push_back(&tree->dnwd_check_surf[depth]);
+        vec_lst.push_back(&tree->dnwd_equiv_surf[depth]);
       }
     }
 
@@ -1931,7 +1931,7 @@ void FMM_Pts<FMMNode>::EvalList(SetupData<Real_t>& setup_data, bool device){
 
 
 template <class FMMNode>
-void FMM_Pts<FMMNode>::Source2UpSetup(SetupData<Real_t>&  setup_data, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
+void FMM_Pts<FMMNode>::Source2UpSetup(SetupData<Real_t>&  setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
   if(this->MultipoleOrder()==0) return;
   { // Set setup_data
     setup_data.level=level;
@@ -1962,7 +1962,7 @@ void FMM_Pts<FMMNode>::Source2UpSetup(SetupData<Real_t>&  setup_data, std::vecto
     input_vector .push_back(&((FMMNode*)nodes_in [i])->surf_value);
   }
   for(size_t i=0;i<nodes_out.size();i++){
-    output_vector.push_back(&upwd_check_surf[((FMMNode*)nodes_out[i])->Depth()]);
+    output_vector.push_back(&tree->upwd_check_surf[((FMMNode*)nodes_out[i])->Depth()]);
     output_vector.push_back(&((FMMData*)((FMMNode*)nodes_out[i])->FMMData())->upward_equiv);
   }
 
@@ -1983,7 +1983,7 @@ void FMM_Pts<FMMNode>::Source2Up(SetupData<Real_t>&  setup_data, bool device){
 
 
 template <class FMMNode>
-void FMM_Pts<FMMNode>::Up2UpSetup(SetupData<Real_t>& setup_data, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
+void FMM_Pts<FMMNode>::Up2UpSetup(SetupData<Real_t>& setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
   if(this->MultipoleOrder()==0) return;
   { // Set setup_data
     setup_data.level=level;
@@ -2716,7 +2716,7 @@ void VListHadamard(size_t dof, size_t M_dim, size_t ker_dim0, size_t ker_dim1, V
 }
 
 template <class FMMNode>
-void FMM_Pts<FMMNode>::V_ListSetup(SetupData<Real_t>&  setup_data, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
+void FMM_Pts<FMMNode>::V_ListSetup(SetupData<Real_t>&  setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
   if(this->MultipoleOrder()==0) return;
   if(level==0) return;
   { // Set setup_data
@@ -3135,7 +3135,7 @@ void FMM_Pts<FMMNode>::V_List     (SetupData<Real_t>&  setup_data, bool device){
 
 
 template <class FMMNode>
-void FMM_Pts<FMMNode>::Down2DownSetup(SetupData<Real_t>& setup_data, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
+void FMM_Pts<FMMNode>::Down2DownSetup(SetupData<Real_t>& setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
   if(this->MultipoleOrder()==0) return;
   { // Set setup_data
     setup_data.level=level;
@@ -3626,7 +3626,7 @@ void FMM_Pts<FMMNode>::EvalListPts(SetupData<Real_t>& setup_data, bool device){
 
 
 template <class FMMNode>
-void FMM_Pts<FMMNode>::X_ListSetup(SetupData<Real_t>&  setup_data, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
+void FMM_Pts<FMMNode>::X_ListSetup(SetupData<Real_t>&  setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
   if(this->MultipoleOrder()==0) return;
   { // Set setup_data
     setup_data.level=level;
@@ -3657,7 +3657,7 @@ void FMM_Pts<FMMNode>::X_ListSetup(SetupData<Real_t>&  setup_data, std::vector<M
     input_vector .push_back(&((FMMNode*)nodes_in [i])->surf_value);
   }
   for(size_t i=0;i<nodes_out.size();i++){
-    output_vector.push_back(&dnwd_check_surf[((FMMNode*)nodes_out[i])->Depth()]);
+    output_vector.push_back(&tree->dnwd_check_surf[((FMMNode*)nodes_out[i])->Depth()]);
     output_vector.push_back(&((FMMData*)((FMMNode*)nodes_out[i])->FMMData())->dnward_equiv);
   }
 
@@ -3679,7 +3679,7 @@ void FMM_Pts<FMMNode>::X_List     (SetupData<Real_t>&  setup_data, bool device){
 
 
 template <class FMMNode>
-void FMM_Pts<FMMNode>::W_ListSetup(SetupData<Real_t>&  setup_data, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
+void FMM_Pts<FMMNode>::W_ListSetup(SetupData<Real_t>&  setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
   if(this->MultipoleOrder()==0) return;
   { // Set setup_data
     setup_data.level=level;
@@ -3704,7 +3704,7 @@ void FMM_Pts<FMMNode>::W_ListSetup(SetupData<Real_t>&  setup_data, std::vector<M
   std::vector<Vector<Real_t>*>&  input_vector=setup_data. input_vector;  input_vector.clear();
   std::vector<Vector<Real_t>*>& output_vector=setup_data.output_vector; output_vector.clear();
   for(size_t i=0;i<nodes_in .size();i++){
-    input_vector .push_back(&upwd_equiv_surf[((FMMNode*)nodes_in [i])->Depth()]);
+    input_vector .push_back(&tree->upwd_equiv_surf[((FMMNode*)nodes_in [i])->Depth()]);
     input_vector .push_back(&((FMMData*)((FMMNode*)nodes_in [i])->FMMData())->upward_equiv);
     input_vector .push_back(NULL);
     input_vector .push_back(NULL);
@@ -3730,7 +3730,7 @@ void FMM_Pts<FMMNode>::W_List     (SetupData<Real_t>&  setup_data, bool device){
 
 
 template <class FMMNode>
-void FMM_Pts<FMMNode>::U_ListSetup(SetupData<Real_t>& setup_data, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
+void FMM_Pts<FMMNode>::U_ListSetup(SetupData<Real_t>& setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
   { // Set setup_data
     setup_data.level=level;
     setup_data.kernel=kernel->k_s2t;
@@ -3782,7 +3782,7 @@ void FMM_Pts<FMMNode>::U_List     (SetupData<Real_t>&  setup_data, bool device){
 
 
 template <class FMMNode>
-void FMM_Pts<FMMNode>::Down2TargetSetup(SetupData<Real_t>&  setup_data, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
+void FMM_Pts<FMMNode>::Down2TargetSetup(SetupData<Real_t>&  setup_data, FMMTree_t* tree, std::vector<Matrix<Real_t> >& buff, std::vector<Vector<FMMNode_t*> >& n_list, int level, bool device){
   if(this->MultipoleOrder()==0) return;
   { // Set setup_data
     setup_data.level=level;
@@ -3807,7 +3807,7 @@ void FMM_Pts<FMMNode>::Down2TargetSetup(SetupData<Real_t>&  setup_data, std::vec
   std::vector<Vector<Real_t>*>&  input_vector=setup_data. input_vector;  input_vector.clear();
   std::vector<Vector<Real_t>*>& output_vector=setup_data.output_vector; output_vector.clear();
   for(size_t i=0;i<nodes_in .size();i++){
-    input_vector .push_back(&dnwd_equiv_surf[((FMMNode*)nodes_in [i])->Depth()]);
+    input_vector .push_back(&tree->dnwd_equiv_surf[((FMMNode*)nodes_in [i])->Depth()]);
     input_vector .push_back(&((FMMData*)((FMMNode*)nodes_in [i])->FMMData())->dnward_equiv);
     input_vector .push_back(NULL);
     input_vector .push_back(NULL);

include/fmm_tree.hpp

@@ -24,6 +24,7 @@ namespace pvfmm{
  */
 template <class FMM_Mat_t>
 class FMM_Tree: public MPI_Tree<typename FMM_Mat_t::FMMNode_t>{
+  friend FMM_Mat_t;
 
  public:
 
@@ -91,8 +92,7 @@ class FMM_Tree: public MPI_Tree<typename FMM_Mat_t::FMMNode_t>{
    */
   void Copy_FMMOutput();
 
- private:
-
+ protected:
 
   std::vector<Matrix<Real_t> > node_data_buff;
   InteracList<Node_t> interac_list;
@@ -101,6 +101,11 @@ class FMM_Tree: public MPI_Tree<typename FMM_Mat_t::FMMNode_t>{
 
   std::vector<Matrix<char> > precomp_lst; //Precomputed data for each interaction type.
   std::vector<SetupData<Real_t> > setup_data;
+
+  std::vector<Vector<Real_t> > upwd_check_surf;
+  std::vector<Vector<Real_t> > upwd_equiv_surf;
+  std::vector<Vector<Real_t> > dnwd_check_surf;
+  std::vector<Vector<Real_t> > dnwd_equiv_surf;
 };
 
 }//end namespace

include/fmm_tree.txx

@@ -72,6 +72,7 @@ void FMM_Tree<FMM_Mat_t>::InitFMM_Tree(bool refine, BoundaryType bndry_) {
 template <class FMM_Mat_t>
 void FMM_Tree<FMM_Mat_t>::SetupFMM(FMM_Mat_t* fmm_mat_) {
   Profile::Tic("SetupFMM",this->Comm(),true);{
+  typedef typename FMM_Mat_t::FMMTree_t MatTree_t;
   bool device=true;
 
   #ifdef __INTEL_OFFLOAD
@@ -111,7 +112,7 @@ void FMM_Tree<FMM_Mat_t>::SetupFMM(FMM_Mat_t* fmm_mat_) {
   }
   //Collect node data into continuous array.
   std::vector<Vector<Node_t*> > node_lists; // TODO: Remove this parameter, not really needed
-  fmm_mat->CollectNodeData(all_nodes, node_data_buff, node_lists);
+  fmm_mat->CollectNodeData((MatTree_t*)this,all_nodes, node_data_buff, node_lists);
   Profile::Toc();
 
   setup_data.resize(8*MAX_DEPTH);
@@ -120,50 +121,50 @@ void FMM_Tree<FMM_Mat_t>::SetupFMM(FMM_Mat_t* fmm_mat_) {
   Profile::Tic("UListSetup",this->Comm(),false,3);
   for(size_t i=0;i<MAX_DEPTH;i++){
     setup_data[i+MAX_DEPTH*0].precomp_data=&precomp_lst[0];
-    fmm_mat->U_ListSetup(setup_data[i+MAX_DEPTH*0],node_data_buff,node_lists,fmm_mat->Homogen()?(i==0?-1:MAX_DEPTH+1):i, device);
+    fmm_mat->U_ListSetup(setup_data[i+MAX_DEPTH*0],(MatTree_t*)this,node_data_buff,node_lists,fmm_mat->Homogen()?(i==0?-1:MAX_DEPTH+1):i, device);
   }
   Profile::Toc();
   Profile::Tic("WListSetup",this->Comm(),false,3);
   for(size_t i=0;i<MAX_DEPTH;i++){
     setup_data[i+MAX_DEPTH*1].precomp_data=&precomp_lst[1];
-    fmm_mat->W_ListSetup(setup_data[i+MAX_DEPTH*1],node_data_buff,node_lists,fmm_mat->Homogen()?(i==0?-1:MAX_DEPTH+1):i, device);
+    fmm_mat->W_ListSetup(setup_data[i+MAX_DEPTH*1],(MatTree_t*)this,node_data_buff,node_lists,fmm_mat->Homogen()?(i==0?-1:MAX_DEPTH+1):i, device);
   }
   Profile::Toc();
   Profile::Tic("XListSetup",this->Comm(),false,3);
   for(size_t i=0;i<MAX_DEPTH;i++){
     setup_data[i+MAX_DEPTH*2].precomp_data=&precomp_lst[2];
-    fmm_mat->X_ListSetup(setup_data[i+MAX_DEPTH*2],node_data_buff,node_lists,fmm_mat->Homogen()?(i==0?-1:MAX_DEPTH+1):i, device);
+    fmm_mat->X_ListSetup(setup_data[i+MAX_DEPTH*2],(MatTree_t*)this,node_data_buff,node_lists,fmm_mat->Homogen()?(i==0?-1:MAX_DEPTH+1):i, device);
   }
   Profile::Toc();
   Profile::Tic("VListSetup",this->Comm(),false,3);
   for(size_t i=0;i<MAX_DEPTH;i++){
     setup_data[i+MAX_DEPTH*3].precomp_data=&precomp_lst[3];
-    fmm_mat->V_ListSetup(setup_data[i+MAX_DEPTH*3],node_data_buff,node_lists,fmm_mat->Homogen()?(i==0?-1:MAX_DEPTH+1):i, /*device*/ false);
+    fmm_mat->V_ListSetup(setup_data[i+MAX_DEPTH*3],(MatTree_t*)this,node_data_buff,node_lists,fmm_mat->Homogen()?(i==0?-1:MAX_DEPTH+1):i, /*device*/ false);
   }
   Profile::Toc();
   Profile::Tic("D2DSetup",this->Comm(),false,3);
   for(size_t i=0;i<MAX_DEPTH;i++){
     setup_data[i+MAX_DEPTH*4].precomp_data=&precomp_lst[4];
-    fmm_mat->Down2DownSetup(setup_data[i+MAX_DEPTH*4],node_data_buff,node_lists,i, /*device*/ false);
+    fmm_mat->Down2DownSetup(setup_data[i+MAX_DEPTH*4],(MatTree_t*)this,node_data_buff,node_lists,i, /*device*/ false);
   }
   Profile::Toc();
   Profile::Tic("D2TSetup",this->Comm(),false,3);
   for(size_t i=0;i<MAX_DEPTH;i++){
     setup_data[i+MAX_DEPTH*5].precomp_data=&precomp_lst[5];
-    fmm_mat->Down2TargetSetup(setup_data[i+MAX_DEPTH*5],node_data_buff,node_lists,fmm_mat->Homogen()?(i==0?-1:MAX_DEPTH+1):i, /*device*/ false);
+    fmm_mat->Down2TargetSetup(setup_data[i+MAX_DEPTH*5],(MatTree_t*)this,node_data_buff,node_lists,fmm_mat->Homogen()?(i==0?-1:MAX_DEPTH+1):i, /*device*/ false);
   }
   Profile::Toc();
 
   Profile::Tic("S2USetup",this->Comm(),false,3);
   for(size_t i=0;i<MAX_DEPTH;i++){
     setup_data[i+MAX_DEPTH*6].precomp_data=&precomp_lst[6];
-    fmm_mat->Source2UpSetup(setup_data[i+MAX_DEPTH*6],node_data_buff,node_lists,fmm_mat->Homogen()?(i==0?-1:MAX_DEPTH+1):i, /*device*/ false);
+    fmm_mat->Source2UpSetup(setup_data[i+MAX_DEPTH*6],(MatTree_t*)this,node_data_buff,node_lists,fmm_mat->Homogen()?(i==0?-1:MAX_DEPTH+1):i, /*device*/ false);
   }
   Profile::Toc();
   Profile::Tic("U2USetup",this->Comm(),false,3);
   for(size_t i=0;i<MAX_DEPTH;i++){
     setup_data[i+MAX_DEPTH*7].precomp_data=&precomp_lst[7];
-    fmm_mat->Up2UpSetup(setup_data[i+MAX_DEPTH*7],node_data_buff,node_lists,i, /*device*/ false);
+    fmm_mat->Up2UpSetup(setup_data[i+MAX_DEPTH*7],(MatTree_t*)this,node_data_buff,node_lists,i, /*device*/ false);
   }
   Profile::Toc();