Merge branch 'multigrid_updateall' into 'master'

fixes for multigrid and changing to smart pointers See merge request jschoeberl/ngsolve!265
NGSolve · Nov 24, 2017 · 8bcaf8793b8f28d660d6b918639e1e504fc42145 · 8bcaf87
2 parents a5835d8 + bbdd678
commit 8bcaf8793b8f28d660d6b918639e1e504fc42145
diff --git a/comp/bddc.cpp b/comp/bddc.cpp
@@ -349,27 +349,27 @@ namespace ngcomp
 	  Flags flags;
 	  flags.SetFlag("eliminate_internal");
 	  flags.SetFlag("subassembled");
-	  cout << IM(1) << "call Create Smoothing Blocks of " << bfa->GetFESpace()->GetName() << endl;
+	  cout << IM(3) << "call Create Smoothing Blocks of " << bfa->GetFESpace()->GetName() << endl;
 	  shared_ptr<Table<int>> blocks = (bfa->GetFESpace()->CreateSmoothingBlocks(flags));
-	  cout << IM(1) << "has blocks" << endl << endl;
+	  cout << IM(3) << "has blocks" << endl << endl;
 	  // *testout << "blocks = " << endl << blocks << endl;
 	  // *testout << "pwbmat = " << endl << *pwbmat << endl;
-	  cout << IM(1) << "call block-jacobi inverse" << endl;
+	  cout << IM(3) << "call block-jacobi inverse" << endl;

 	  inv = dynamic_pointer_cast<BaseSparseMatrix> (pwbmat)->CreateBlockJacobiPrecond(blocks, 0, 0, 0);
 	  // inv = dynamic_cast<BaseSparseMatrix*> (pwbmat)->CreateJacobiPrecond(wb_free_dofs);

-	  cout << IM(1) << "has inverse" << endl << endl;
+	  cout << IM(3) << "has inverse" << endl << endl;
 	  // *testout << "blockjacobi = " << endl << *inv << endl;
 	  
 	  //Coarse Grid of Wirebasket
-	  cout << IM(1) << "call directsolverclusters inverse" << endl;
-	  Array<int> & clusters = *(bfa->GetFESpace()->CreateDirectSolverClusters(flags));
-	  cout << IM(1) << "has clusters" << endl << endl;
+	  cout << IM(3) << "call directsolverclusters inverse" << endl;
+	  auto clusters = bfa->GetFESpace()->CreateDirectSolverClusters(flags);
+	  cout << IM(3) << "has clusters" << endl << endl;
 	  
-	  cout << IM(1) << "call coarse wirebasket grid inverse" << endl;
-	  inv_coarse = pwbmat->InverseMatrix(&clusters);
-	  cout << IM(1) << "has inverse" << endl << endl;
+	  cout << IM(3) << "call coarse wirebasket grid inverse" << endl;
+	  inv_coarse = pwbmat->InverseMatrix(clusters);
+	  cout << IM(3) << "has inverse" << endl << endl;
 	  
 	  tmp = new VVector<>(ndof);
 	  tmp2 = new VVector<>(ndof);
@@ -412,18 +412,18 @@ namespace ngcomp

              if (coarse)
              {
-                cout << IM(1) << "call wirebasket preconditioner finalize ( with " << cntfreedofs
+                cout << IM(3) << "call wirebasket preconditioner finalize ( with " << cntfreedofs
                     << " free dofs out of " << pwbmat->Height() << " )" << endl;
                // throw Exception("combination of coarse and block not implemented! ");
                dynamic_pointer_cast<Preconditioner>(inv) -> FinalizeLevel(pwbmat.get());
              }
              else
              {
-                cout << IM(1) << "call wirebasket inverse ( with " << cntfreedofs
+                cout << IM(3) << "call wirebasket inverse ( with " << cntfreedofs
                     << " free dofs out of " << pwbmat->Height() << " )" << endl;
                inv = pwbmat->InverseMatrix(wb_free_dofs);
              }
-	      cout << IM(1) << "has inverse" << endl;
+	      cout << IM(3) << "has inverse" << endl;
 	      tmp = new VVector<TV>(ndof);
 	    }
 	}
diff --git a/comp/bilinearform.cpp b/comp/bilinearform.cpp
@@ -5127,11 +5127,8 @@ namespace ngcomp
            Restrict( *prolMat, &( dynamic_cast< BaseSparseMatrix& >
            ( GetMatrix( i-1 ) ) ) ) );
          */
-          mats[i-1] =
-            shared_ptr<BaseMatrix>
-            ( dynamic_cast< const BaseSparseMatrix& >(GetMatrix(i) ).
-              Restrict( *prolMat, &( dynamic_cast< BaseSparseMatrix& >
-                                     (GetMatrix( i-1 ) ) ) ) );
+          mats[i-1] = dynamic_cast< const BaseSparseMatrix& >(GetMatrix(i)).
+            Restrict(*prolMat,dynamic_pointer_cast<BaseSparseMatrix>(GetMatrixPtr(i-1)));
          
          delete prolMat;
        }
diff --git a/comp/facetfespace.cpp b/comp/facetfespace.cpp
@@ -855,28 +855,29 @@ namespace ngcomp


  
-  Array<int> * FacetFESpace :: CreateDirectSolverClusters (const Flags & precflags) const
+  shared_ptr<Array<int>> FacetFESpace :: CreateDirectSolverClusters (const Flags & precflags) const
  {
-		if (all_dofs_together)
-			throw Exception("FacetFESpace ::CreateDirectSolverClusters not implemented for case all_dofs_together!");
+    if (all_dofs_together)
+      throw Exception("FacetFESpace ::CreateDirectSolverClusters not implemented for case all_dofs_together!");

-		Array<int> & clusters = *new Array<int> (GetNDof());
+    auto spclusters = make_shared<Array<int>>(GetNDof());
+    Array<int> & clusters = *spclusters;

    clusters.SetSize(ndof);
    clusters = 0;
    
    for (int i = 0; i < nfa; i++)
      clusters[i] = 1;
    
-    return &clusters;
+    return spclusters;
    
    //
  
    for (int i=0; i<nfa-ncfa; i++)
      clusters[i]  = 1;
  
    // cout << "direct solver cluster = " << clusters << endl;
-    return & clusters;
+    return spclusters;
  }


@@ -982,12 +983,13 @@ namespace ngcomp
  HybridDGFESpace :: ~HybridDGFESpace () { ; }


-  Array<int> * HybridDGFESpace :: CreateDirectSolverClusters (const Flags & flags) const
+  shared_ptr<Array<int>> HybridDGFESpace :: CreateDirectSolverClusters (const Flags & flags) const
  {
    if (flags.GetDefineFlag("subassembled"))
      {
-	cout << "creating bddc-coarse grid(vertices)" << endl;
-	Array<int> & clusters = *new Array<int> (GetNDof());
+	cout << IM(3) << "creating bddc-coarse grid(vertices)" << endl;
+        auto spclusters = make_shared<Array<int>> (GetNDof());
+	Array<int> & clusters = *spclusters;
 	clusters = 0;
 	/*
          int nv = ma->GetNV();
@@ -1037,11 +1039,12 @@ namespace ngcomp
        }
 	*/

-	return &clusters;	
+	return spclusters;
      }
    else
      {
-	Array<int> & clusters = *new Array<int> (GetNDof());
+        auto spclusters = make_shared<Array<int>> (GetNDof());
+	Array<int> & clusters = *spclusters;
 	clusters = 0;

 	Array<DofId> dnums;
@@ -1061,7 +1064,7 @@ namespace ngcomp
 	for (int i = 0; i < freedofs.Size(); i++)
 	  if (!freedofs.Test(i)) clusters[i] = 0;
 	*testout << "Hybrid-FESpace, dsc = " << endl << clusters << endl;
-	return &clusters;
+	return spclusters;
      }
  }

diff --git a/comp/facetfespace.hpp b/comp/facetfespace.hpp
@@ -95,7 +95,7 @@ namespace ngcomp
    ///
    virtual shared_ptr<Table<int>> CreateSmoothingBlocks (const Flags & precflags) const override;
    ///
-    virtual Array<int> * CreateDirectSolverClusters (const Flags & precflags) const override;
+    virtual shared_ptr<Array<int>> CreateDirectSolverClusters (const Flags & precflags) const override;


    virtual INT<2> GetFacetOrder(int fnr) 
@@ -155,7 +155,7 @@ namespace ngcomp
    virtual ~HybridDGFESpace ();
    virtual string GetClassName () const { return "HybridDGFESpace"; }

-    virtual Array<int> * CreateDirectSolverClusters (const Flags & flags) const;
+    virtual shared_ptr<Array<int>> CreateDirectSolverClusters (const Flags & flags) const override;
    virtual shared_ptr<Table<int>> CreateSmoothingBlocks (const Flags & precflags) const;
  };

diff --git a/comp/fespace.cpp b/comp/fespace.cpp
@@ -1802,10 +1802,11 @@ lot of new non-zero entries in the matrix!\n" << endl;
  }

  
-  Array<int> * 
+  shared_ptr<Array<int>>
  NodalFESpace :: CreateDirectSolverClusters (const Flags & flags) const
  {
-    Array<int> & clusters = *new Array<int> (GetNDof());
+    auto spclusters = make_shared<Array<int>> (GetNDof());
+    Array<int> & clusters = *spclusters;
    clusters = 0;

    const int stdoffset = 1;
@@ -1829,11 +1830,10 @@ lot of new non-zero entries in the matrix!\n" << endl;
      if (clusters[i]) nonzero = true;
    if (!nonzero)
      {
-	delete &clusters;
-	return 0;
+	return nullptr;
      }

-    return &clusters;
+    return spclusters;
  }


diff --git a/comp/fespace.hpp b/comp/fespace.hpp
@@ -592,8 +592,8 @@ namespace ngcomp
    /// 
    virtual shared_ptr<Table<int>> CreateSmoothingBlocks (const Flags & flags) const;
    /// for anisotropic plane smoothing:
-    virtual Array<int> * CreateDirectSolverClusters (const Flags & flags) const
-    { return 0; }
+    virtual shared_ptr<Array<int>> CreateDirectSolverClusters (const Flags & flags) const
+    { return nullptr; }

    virtual void AddDirectSolverClusterDof(int dn) const
    { adddirectsolverdofs.Append(dn); }
@@ -967,7 +967,7 @@ namespace ngcomp
    virtual void GetFaceDofNrs (int fanr, Array<DofId> & dnums) const override;
    virtual void GetInnerDofNrs (int elnr, Array<DofId> & dnums) const override;

-    virtual Array<int> * CreateDirectSolverClusters (const Flags & flags) const override;
+    virtual shared_ptr<Array<int>> CreateDirectSolverClusters (const Flags & flags) const override;
  };


diff --git a/comp/h1hofespace.cpp b/comp/h1hofespace.cpp
@@ -986,7 +986,7 @@ namespace ngcomp
    size_t nfa = (ma->GetDimension() == 2) ? 0 : ma->GetNFaces();
    size_t ni = (eliminate_internal) ? 0 : ma->GetNE(); 
   
-    cout << " blocktype " << smoothing_type << endl; 
+    cout << IM(4) << " blocktype " << smoothing_type << endl;
    // cout << " Use H1-Block Smoother:  "; 

    FilteredTableCreator creator(GetFreeDofs().get());
@@ -1065,7 +1065,7 @@ namespace ngcomp
 		
            if (creator.GetMode() == 1)
              {
-                cout << " VE + F + I " << endl;
+                cout << IM(4) << " VE + F + I " << endl;
                creator.SetSize(nv+nfa+ni);
                break;
              }
@@ -1410,26 +1410,28 @@ namespace ngcomp
    


-  Array<int> * 
+  shared_ptr<Array<int>>
  H1HighOrderFESpace :: CreateDirectSolverClusters (const Flags & flags) const
  {
    if (flags.GetDefineFlag("subassembled"))
    {
 	cout << "creating bddc-coarse grid(vertices)" << endl;
-	Array<int> & clusters = *new Array<int> (GetNDof());
+        auto spclusters = make_shared<Array<int>> (GetNDof());
+	Array<int> & clusters = *spclusters;
 	clusters = 0;
 	int nv = ma->GetNV();
 	for (int i = 0; i < nv; i++)
 	  if (!IsDirichletVertex(i))
 	    clusters[i] = 1;		
-	return &clusters;	
+	return spclusters;
    }
    
    if (flags.NumFlagDefined ("ds_order"))
      {
 	int ds_order = int (flags.GetNumFlag ("ds_order", 1));

-	Array<int> & clusters = *new Array<int> (GetNDof());
+        auto spclusters = make_shared<Array<int>> (GetNDof());
+	Array<int> & clusters = *spclusters;
 	clusters = 0;
 	
 	int ned = ma->GetNEdges();
@@ -1475,7 +1477,7 @@ namespace ngcomp
 	  }
 	*/

-	return &clusters;
+	return spclusters;
      }


@@ -1486,7 +1488,8 @@ namespace ngcomp
    // int nv = ma->GetNV();
    // int nd = GetNDof();
    int ne = ma->GetNE();
-    Array<int> & clusters = *new Array<int> (GetNDof());
+    auto spclusters = make_shared<Array<int>> (GetNDof());
+    Array<int> & clusters = *spclusters;
    clusters = 0;

    // all vertices in global space
@@ -1606,11 +1609,10 @@ namespace ngcomp
      if (clusters[i]) nonzero = true;
    if (!nonzero)
      {
-	delete &clusters;
-	return 0;
+	return nullptr;
      }

-    return &clusters;
+    return spclusters;
  }

  template<>
diff --git a/comp/h1hofespace.hpp b/comp/h1hofespace.hpp
@@ -102,7 +102,7 @@ namespace ngcomp
    virtual shared_ptr<Table<int>> CreateSmoothingBlocks (const Flags & precflags) const override; 
    // virtual void CreateSmoothingBlocks2 (SmoothingBlocksCreator & sbc, const Flags & precflags) const; 
    ///
-    virtual Array<int> * CreateDirectSolverClusters (const Flags & flags) const override;
+    virtual shared_ptr<Array<int>> CreateDirectSolverClusters (const Flags & flags) const override;

    void UpdateDofTables ();
    ///
diff --git a/comp/hcurlhdivfes.cpp b/comp/hcurlhdivfes.cpp
@@ -2488,7 +2488,7 @@ namespace ngcomp
  }


-  Array<int> * 
+  shared_ptr<Array<int>>
  NedelecFESpace2 :: CreateDirectSolverClusters (const Flags & flags) const
  {
    (*testout) << "CreateDirectSolverClusters" << endl;
@@ -2498,7 +2498,8 @@ namespace ngcomp
    int ne = ma->GetNE();


-    Array<int> & clusters = *new Array<int> (nd);
+    auto spclusters = make_shared<Array<int>> (nd);
+    Array<int> & clusters = *spclusters;
    clusters = 0;


@@ -2624,7 +2625,7 @@ namespace ngcomp
    //     }
 	  
 	 
-    return &clusters;
+    return spclusters;
    

  }
diff --git a/comp/hcurlhdivfes.hpp b/comp/hcurlhdivfes.hpp
@@ -267,7 +267,7 @@ class NGS_DLL_HEADER NedelecFESpace2 : public FESpace
  SparseMatrix<double> * CreateGradient() const;

  
-  virtual Array<int> * CreateDirectSolverClusters (const Flags & flags) const override;
+  virtual shared_ptr<Array<int>> CreateDirectSolverClusters (const Flags & flags) const override;


  virtual void GetVertexDofNrs (int vnr, Array<DofId> & dnums) const override;