]> https://gitweb.dealii.org/ - dealii.git/commitdiff
new multilevel tests
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Tue, 8 Mar 2005 15:44:26 +0000 (15:44 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Tue, 8 Mar 2005 15:44:26 +0000 (15:44 +0000)
git-svn-id: https://svn.dealii.org/trunk@10044 0785d39b-7218-0410-832d-ea1e28bc413d

tests/multigrid/Makefile
tests/multigrid/cycles.cc [new file with mode: 0644]
tests/multigrid/mg1.cc [deleted file]
tests/multigrid/mg1.checked [deleted file]
tests/multigrid/transfer.cc

index fb01f0f972a639ae07b2d6be03c5f5804d77d66e..1c701056c03205dc65a3068c444f612b11e6436a 100644 (file)
@@ -1,7 +1,6 @@
-# Generated automatically from Makefile.in by configure.
 ############################################################
 # $Id$
-# Copyright (C) 2000, 2001, 2002, 2003 by the deal.II authors
+# Copyright (C) 2000, 2001, 2002, 2003, 2005 by the deal.II authors
 ############################################################
 
 ############################################################
@@ -23,9 +22,9 @@ default: run-tests
 ############################################################
 
 transfer.exe      : transfer.g.$(OBJEXT) $(lib-2d)           $(libraries)
-mg1.exe           : mg1.g.$(OBJEXT)      $(lib-2d)           $(libraries)
+cycles.exe        : cycles.g.$(OBJEXT)   $(lib-2d)           $(libraries)
 
-tests = transfer
+tests = cycles transfer
 
 ############################################################
 
diff --git a/tests/multigrid/cycles.cc b/tests/multigrid/cycles.cc
new file mode 100644 (file)
index 0000000..15b330a
--- /dev/null
@@ -0,0 +1,88 @@
+//----------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2000 - 2005 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------------------------------------------------------
+
+/*
+ * Test the general behavior of the multigrid cycles without any
+ * numerics. Therefore, all transfer operators are void and we use the
+ * same matrix on each level.
+ */
+
+#include <base/logstream.h>
+#include <lac/vector.h>
+#include <lac/full_matrix.h>
+#include <multigrid/mg_base.h>
+#include <multigrid/multigrid.h>
+#include <multigrid/mg_level_object.h>
+#include <multigrid/mg_matrix.h>
+
+#include <fstream>
+
+#define N 3
+typedef Vector<double> VECTOR;
+
+class MGAll
+  :
+  public MGSmootherBase<VECTOR>,
+  public MGTransferBase<VECTOR>,
+  public MGCoarseGridBase<VECTOR>
+{
+  public:
+    virtual ~MGAll()
+      {}
+    
+    virtual void smooth (const unsigned int,
+                        VECTOR &, const VECTOR &) const
+      {}
+
+    virtual void prolongate (const unsigned int,
+                            VECTOR &, const VECTOR &) const
+      {}
+    
+    virtual void restrict_and_add (const unsigned int,
+                                  VECTOR &, const VECTOR &) const
+      {}
+    
+    virtual void operator() (const unsigned int,
+                            VECTOR &, const VECTOR &) const
+      {}
+};
+
+void test_cycles(unsigned int minlevel, unsigned int maxlevel)
+{
+  MGAll all;
+  MGLevelObject<FullMatrix<double> > level_matrices(0, maxlevel);
+  for (unsigned int i=0;i<=maxlevel;++i)
+    level_matrices[i].reinit(N, N);
+  MGMatrix<FullMatrix<double>, VECTOR> mgmatrix(&level_matrices);
+  
+  Multigrid<VECTOR> mg1(minlevel, maxlevel, mgmatrix, all, all, all, all,
+                       Multigrid<VECTOR>::v_cycle);
+  mg1.set_debug(3);
+  for (unsigned int i=minlevel;i<=maxlevel;++i)
+    mg1.defect[i].reinit(N);
+  mg1.cycle();
+  deallog << std::endl;
+  
+  mg1.set_cycle(Multigrid<VECTOR>::w_cycle);
+  mg1.cycle();
+}
+
+int main()
+{
+  std::ofstream logfile("cycles.output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  
+  test_cycles (0,4);
+  test_cycles (2,5);
+}
diff --git a/tests/multigrid/mg1.cc b/tests/multigrid/mg1.cc
deleted file mode 100644 (file)
index 5dd8394..0000000
+++ /dev/null
@@ -1,205 +0,0 @@
-//----------------------------  mg.cc  ---------------------------
-//    $Id$
-//    Version: $Name$
-//
-//    Copyright (C) 1998, 1999, 2000, 2001 by the deal.II authors
-//
-//    This file is subject to QPL and may not be  distributed
-//    without copyright and license information. Please refer
-//    to the file deal.II/doc/license.html for the  text  and
-//    further information on this license.
-//
-//----------------------------  mg.cc  ---------------------------
-
-
-#include <base/function.h>
-#include <base/quadrature_lib.h>
-#include <base/logstream.h>
-#include <lac/vector.h>
-#include <lac/sparse_matrix.h>
-#include <lac/solver_cg.h>
-#include <lac/solver_richardson.h>
-#include <lac/vector_memory.h>
-#include <lac/precondition.h>
-#include <grid/tria.h>
-#include <grid/tria_iterator.h>
-#include <dofs/dof_constraints.h>
-#include <dofs/dof_accessor.h>
-#include <grid/grid_generator.h>
-#include <fe/fe_q.h>
-#include <fe/mapping_q1.h>
-#include <fe/fe_values.h>
-#include <dofs/dof_tools.h>
-
-#include <multigrid/mg_dof_handler.h>
-#include <multigrid/mg_dof_accessor.h>
-#include <multigrid/mg_dof_tools.h>
-#include <multigrid/mg_base.h>
-#include <multigrid/mg_smoother.h>
-#include <multigrid/multigrid.h>
-
-MappingQ1<2> mapping;
-
-//#include "helmholtz.h"
-
-#include <fstream>
-
-
-template<int dim>
-class RHSFunction : public Function<dim>
-{
-  public:
-    virtual double value (const Point<dim>&,
-                         const unsigned int) const;
-};
-
-
-class MGSmootherLAC : public MGSmootherBase
-{
-  private:
-    SmartPointer<MGLevelObject<SparseMatrix<double> > >matrices;
-  public:
-    MGSmootherLAC(MGLevelObject<SparseMatrix<double> >&);
-    
-    virtual void smooth (const unsigned int level,
-                        Vector<double> &u,
-                        const Vector<double> &rhs) const;
-    
-};
-
-int main()
-{
-  try 
-    {
-  std::ofstream logfile("mg.output");
-  deallog.attach(logfile);
-  deallog.depth_console(0);
-
-  Helmholtz equation;
-  RHSFunction<2> rhs;
-  QGauss5<2> quadrature;
-  
-  FE_Q<2> fe1(1);
-  FE_Q<2> fe2(2);
-  FE_Q<2> fe3(3);
-  FE_Q<2> fe4(4);
-  for (unsigned int degree=1;degree<=3;degree++)
-    {
-      Triangulation<2> tr;
-      MGDoFHandler<2> mgdof(tr);
-      DoFHandler<2>& dof(mgdof);
-  
-      GridGenerator::hyper_cube(tr,-1.,1.);
-  
-      FiniteElement<2>* fe;
-      switch(degree)
-       {
-         case 1: fe = &fe1; deallog.push("Q1"); break;
-         case 2: fe = &fe2; deallog.push("Q2"); break;
-         case 3: fe = &fe3; deallog.push("Q3"); break;
-         case 4: fe = &fe4; deallog.push("Q4"); break;
-       }
-      
-      for (unsigned int step=0;step < 3; ++step)
-       {
-         tr.refine_global(1);
-         dof.distribute_dofs(*fe);
-
-         ConstraintMatrix hanging_nodes;
-         DoFTools::make_hanging_node_constraints(dof, hanging_nodes);
-         hanging_nodes.close();
-
-         const unsigned int size = dof.n_dofs();
-         deallog << "DoFs " << size << std::endl;
-         deallog << "Levels: " << tr.n_levels() << std::endl;
-         
-         SparsityPattern structure(size, dof.max_couplings_between_dofs());
-         DoFTools::make_sparsity_pattern(dof, structure);
-         structure.compress();
-         
-         SparseMatrix<double> A(structure);
-         Vector<double> f(size);
-         
-         equation.build_all(A,f,dof, quadrature, rhs);
-         
-         Vector<double> u;
-         u.reinit(f);
-         PrimitiveVectorMemory<> mem;
-         SolverControl control(1000, 1.e-10,false , true);
-         SolverCG<>    solver(control, mem);
-         
-         u = 0.;
-         MGLevelObject<SparsityPattern> mgstruct(0, tr.n_levels()-1);
-         MGLevelObject<SparseMatrix<double> > mgA(0,tr.n_levels()-1);
-         for (unsigned int i=0;i<tr.n_levels();++i)
-           {
-             mgstruct[i].reinit(mgdof.n_dofs(i), mgdof.n_dofs(i),
-                                mgdof.max_couplings_between_dofs());
-             MGTools::make_sparsity_pattern(mgdof, mgstruct[i], i);
-             mgstruct[i].compress();
-             
-             mgA[i].reinit(mgstruct[i]);
-           }
-         equation.build_mgmatrix(mgA, mgdof, quadrature);
-
-         SolverControl cgcontrol(1000,1.e-12, false, false);
-         PrimitiveVectorMemory<> cgmem;
-         SolverCG<> cgsolver(cgcontrol, cgmem);
-         PreconditionIdentity cgprec;
-         MGCoarseGridLACIteration<SolverCG<>, SparseMatrix<double>, PreconditionIdentity>
-           coarse(cgsolver, mgA[0], cgprec);
-         
-         MGSmootherLAC smoother(mgA);
-         MGTransferPrebuilt transfer;
-         transfer.build_matrices(mgdof);
-
-
-Multigrid<2> multigrid(mgdof, hanging_nodes, mgstruct, mgA, transfer);
-         PreconditionMG<Multigrid<2> >
-           mgprecondition(multigrid, smoother, smoother, coarse);
-         
-         solver.solve(A, u, f, mgprecondition);
-       }
-      deallog.pop();
-    }
-
-    }
-  catch (const SolverControl::NoConvergence& e)
-    {
-      deallog << "Step " << e.last_step
-             << " Residual " << e.last_residual << std::endl;
-    }
-}
-
-template<int dim>
-double
-RHSFunction<dim>::value (const Point<dim>&,
-                        const unsigned int) const
-{
-  return 1.;
-}
-
-MGSmootherLAC::MGSmootherLAC(MGLevelObject<SparseMatrix<double> >& matrix)
-               :
-               matrices(&matrix)
-{}
-
-void
-MGSmootherLAC::smooth (const unsigned int level,
-                      Vector<double> &u,
-                      const Vector<double> &rhs) const
-{
-  SolverControl control(2,1.e-300,false,false);
-  PrimitiveVectorMemory<> mem;
-  SolverRichardson<> rich(control, mem);
-  PreconditionSSOR<> prec;
-  prec.initialize((*matrices)[level], 1.);
-
-  try
-    {
-      rich.solve((*matrices)[level], u, rhs, prec);
-    }
-  catch (const SolverControl::NoConvergence&)
-    {
-    }
-}
diff --git a/tests/multigrid/mg1.checked b/tests/multigrid/mg1.checked
deleted file mode 100644 (file)
index 6079191..0000000
+++ /dev/null
@@ -1,37 +0,0 @@
-
-DEAL:Q1::DoFs 9
-DEAL:Q1::Levels: 2
-DEAL:Q1:cg::Starting 
-DEAL:Q1:cg::Convergence step 3 
-DEAL:Q1::DoFs 25
-DEAL:Q1::Levels: 3
-DEAL:Q1:cg::Starting 
-DEAL:Q1:cg::Convergence step 4 
-DEAL:Q1::DoFs 81
-DEAL:Q1::Levels: 4
-DEAL:Q1:cg::Starting 
-DEAL:Q1:cg::Convergence step 5 
-DEAL:Q2::DoFs 25
-DEAL:Q2::Levels: 2
-DEAL:Q2:cg::Starting 
-DEAL:Q2:cg::Convergence step 5 
-DEAL:Q2::DoFs 81
-DEAL:Q2::Levels: 3
-DEAL:Q2:cg::Starting 
-DEAL:Q2:cg::Convergence step 6 
-DEAL:Q2::DoFs 289
-DEAL:Q2::Levels: 4
-DEAL:Q2:cg::Starting 
-DEAL:Q2:cg::Convergence step 6 
-DEAL:Q3::DoFs 49
-DEAL:Q3::Levels: 2
-DEAL:Q3:cg::Starting 
-DEAL:Q3:cg::Convergence step 6 
-DEAL:Q3::DoFs 169
-DEAL:Q3::Levels: 3
-DEAL:Q3:cg::Starting 
-DEAL:Q3:cg::Convergence step 6 
-DEAL:Q3::DoFs 625
-DEAL:Q3::Levels: 4
-DEAL:Q3:cg::Starting 
-DEAL:Q3:cg::Convergence step 6 
index fd7a3ad31203e86b1f13b52938b427efbd31fe54..70d213e08a639b32672c3c7fb90f0c5685a2c5dd 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2000, 2001, 2002, 2003, 2004 by the deal.II authors
+//    Copyright (C) 2000 - 2005 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -86,6 +86,7 @@ void check_select(const FiniteElement<dim>& fe,
   tr.refine_global(2);
   
   MGDoFHandler<dim> mgdof(tr);
+  DoFHandler<dim>& dof=mgdof;
   mgdof.distribute_dofs(fe);
   DoFRenumbering::component_wise(mgdof, target_component);
   vector<unsigned int> ndofs(fe.n_components());
@@ -109,7 +110,7 @@ void check_select(const FiniteElement<dim>& fe,
   
   
   MGTransferSelect<double> transfer;
-  transfer.build_matrices(mgdof, selected, mg_selected,
+  transfer.build_matrices(dof, mgdof, selected, mg_selected,
                          target_component, mg_target_component);
 
   Vector<double> u2(mg_ndofs[2][mg_selected]);
@@ -221,7 +222,7 @@ int main()
   std::ofstream logfile("transfer.output");
   logfile.precision(3);
   deallog.attach(logfile);
-  deallog.depth_console(10);
+  deallog.depth_console(0);
   
 //  check_simple (FE_DGP<2>(0));
 //  check_simple (FE_DGP<2>(1));
@@ -230,8 +231,6 @@ int main()
 //  check_simple (FE_Q<2>(1));
 //  check_simple (FE_Q<2>(2));
 
-  std::cerr << "Select" << std::endl;
-  
   std::vector<unsigned int> v1(4);
   std::vector<unsigned int> v2(4);
   std::vector<unsigned int> v3(4);

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.