]> https://gitweb.dealii.org/ - dealii.git/commitdiff
test Assembler::Functional
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Tue, 15 Feb 2011 17:58:08 +0000 (17:58 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Tue, 15 Feb 2011 17:58:08 +0000 (17:58 +0000)
git-svn-id: https://svn.dealii.org/trunk@23367 0785d39b-7218-0410-832d-ea1e28bc413d

tests/integrators/empty_info.h [new file with mode: 0644]
tests/integrators/functional_01.cc [new file with mode: 0644]
tests/integrators/functional_01/cmp/generic [new file with mode: 0644]

diff --git a/tests/integrators/empty_info.h b/tests/integrators/empty_info.h
new file mode 100644 (file)
index 0000000..9fcd35e
--- /dev/null
@@ -0,0 +1,46 @@
+//----------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2000, 2001, 2003, 2004, 2007, 2008, 2009, 2010, 2011 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.
+//
+//----------------------------------------------------------------------
+
+// Provide scratch data objects with no contents for the MeshWorker::loop()
+
+#include <numerics/mesh_worker_info.h>
+
+class EmptyInfo
+{
+  public:
+    template <class DOFINFO>
+    void reinit(const DOFINFO&)
+      {}
+};
+
+
+class EmptyInfoBox
+{
+   public:
+    typedef EmptyInfo CellInfo;
+    template <int dim, class DOFINFO>
+    void post_cell(const MeshWorker::DoFInfoBox<dim, DOFINFO>&)
+      {}
+    
+    template <int dim, class DOFINFO>
+    void post_faces(const MeshWorker::DoFInfoBox<dim, DOFINFO>&)
+      {}
+    
+    EmptyInfo cell;
+    EmptyInfo boundary;
+    EmptyInfo face;
+    EmptyInfo subface;
+    EmptyInfo neighbor;
+};
+
+
diff --git a/tests/integrators/functional_01.cc b/tests/integrators/functional_01.cc
new file mode 100644 (file)
index 0000000..0de3225
--- /dev/null
@@ -0,0 +1,169 @@
+//----------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2000, 2001, 2003, 2004, 2007, 2008, 2009, 2010, 2011 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 whether Assember::Functional adds up correctly.
+
+#include "../tests.h"
+#include "empty_info.h"
+#include <numerics/mesh_worker.h>
+#include <numerics/mesh_worker_assembler.h>
+#include <numerics/mesh_worker_loop.h>
+
+#include <base/std_cxx1x/function.h>
+#include <base/logstream.h>
+#include <grid/grid_generator.h>
+#include <dofs/dof_tools.h>
+#include <fe/fe_dgp.h>
+
+#include <fstream>
+#include <iomanip>
+
+using namespace dealii;
+
+
+// Define a class that fills all available entries in the info objects
+// with recognizable numbers.
+
+// Fills the following functionals:
+// 0: number of cells
+// 1: number of interior faces
+// 2: number of boundary faces
+
+const unsigned int n_functionals = 3;
+
+template <int dim>
+class Local : public Subscriptor
+{
+  public:
+    typedef EmptyInfo CellInfo;
+    
+    void cell(MeshWorker::DoFInfo<dim>& dinfo, CellInfo& info) const;
+    void bdry(MeshWorker::DoFInfo<dim>& dinfo, CellInfo& info) const;
+    void face(MeshWorker::DoFInfo<dim>& dinfo1, MeshWorker::DoFInfo<dim>& dinfo2,
+             CellInfo& info1, CellInfo& info2) const;
+};
+
+
+template <int dim>
+void
+Local<dim>::cell(MeshWorker::DoFInfo<dim>& info, CellInfo&) const
+{
+  info.value(0) = 1.;
+}
+
+
+template <int dim>
+void
+Local<dim>::bdry(MeshWorker::DoFInfo<dim>&  info, CellInfo&) const
+{
+  info.value(2) = 1.;
+}
+
+
+template <int dim>
+void
+Local<dim>::face(MeshWorker::DoFInfo<dim>&  info1, MeshWorker::DoFInfo<dim>& info2,
+                CellInfo&, CellInfo&) const
+{
+  info1.value(1) = 1./2.;
+  info2.value(1) = 1./2.;
+}
+
+
+template <int dim>
+void
+test_mesh(MGDoFHandler<dim>& mgdofs)
+{
+  const DoFHandler<dim>& dofs = mgdofs;
+  
+  Local<dim> local;
+  EmptyInfoBox info_box;
+  MeshWorker::DoFInfo<dim> dof_info(dofs);
+  
+  MeshWorker::Assembler::Functional<double> assembler;
+  assembler.initialize(n_functionals);
+  
+  MeshWorker::loop<dim, dim, MeshWorker::DoFInfo<dim>, EmptyInfoBox>
+    (dofs.begin_active(), dofs.end(),
+     dof_info, info_box,
+     std_cxx1x::bind (&Local<dim>::cell, local, _1, _2),
+     std_cxx1x::bind (&Local<dim>::bdry, local, _1, _2),
+     std_cxx1x::bind (&Local<dim>::face, local, _1, _2, _3, _4),
+     assembler, true);
+
+  deallog << "  Results";
+  for (unsigned int i=0;i<n_functionals;++i)
+    deallog << '\t' << assembler(i);
+  deallog << std::endl;
+
+  assembler.initialize(n_functionals);
+  MeshWorker::DoFInfo<dim> mg_dof_info(mgdofs);
+  MeshWorker::loop<dim, dim, MeshWorker::DoFInfo<dim>, EmptyInfoBox>
+    (mgdofs.begin(), mgdofs.end(),
+     mg_dof_info, info_box,
+     std_cxx1x::bind (&Local<dim>::cell, local, _1, _2),
+     std_cxx1x::bind (&Local<dim>::bdry, local, _1, _2),
+     std_cxx1x::bind (&Local<dim>::face, local, _1, _2, _3, _4),
+     assembler, true);
+
+  deallog << "MGResults";
+  for (unsigned int i=0;i<n_functionals;++i)
+    deallog << '\t' << assembler(i);
+  deallog << std::endl;
+}
+
+
+template<int dim>
+void
+test(const FiniteElement<dim>& fe)
+{
+  Triangulation<dim> tr;
+  MGDoFHandler<dim> dofs(tr);
+  GridGenerator::hyper_cube(tr);
+  tr.refine_global(1);
+  deallog.push("1");
+  dofs.distribute_dofs(fe);
+  test_mesh(dofs);
+  deallog.pop();
+  tr.begin(1)->set_refine_flag();
+  tr.execute_coarsening_and_refinement();
+  deallog.push("2");
+  dofs.distribute_dofs(fe);
+  test_mesh(dofs);
+  deallog.pop();
+  tr.begin(2)->set_refine_flag();
+  tr.execute_coarsening_and_refinement();
+  deallog.push("3");
+  dofs.distribute_dofs(fe);
+  test_mesh(dofs);
+  deallog.pop();
+}
+
+
+int main ()
+{
+  const std::string logname = JobIdentifier::base_name(__FILE__) + std::string("/output");
+  std::ofstream logfile(logname.c_str());
+  deallog.attach(logfile);
+//  deallog.depth_console (0);
+
+  FE_DGP<2> el2(0);
+  FE_DGP<3> el3(0);
+
+  deallog.push("2D");  
+  test(el2);
+  deallog.pop();
+  deallog.push("3D");
+  test(el3);
+  deallog.pop();
+}
diff --git a/tests/integrators/functional_01/cmp/generic b/tests/integrators/functional_01/cmp/generic
new file mode 100644 (file)
index 0000000..0c72628
--- /dev/null
@@ -0,0 +1,13 @@
+
+DEAL:2D:1::  Results   4.00000 4.00000 8.00000
+DEAL:2D:1::MGResults   5.00000 4.00000 12.0000
+DEAL:2D:2::  Results   7.00000 10.0000 10.0000
+DEAL:2D:2::MGResults   9.00000 12.0000 16.0000
+DEAL:2D:3::  Results   10.0000 16.0000 12.0000
+DEAL:2D:3::MGResults   13.0000 20.0000 20.0000
+DEAL:3D:1::  Results   8.00000 12.0000 24.0000
+DEAL:3D:1::MGResults   9.00000 12.0000 30.0000
+DEAL:3D:2::  Results   15.0000 33.0000 33.0000
+DEAL:3D:2::MGResults   17.0000 36.0000 42.0000
+DEAL:3D:3::  Results   22.0000 54.0000 42.0000
+DEAL:3D:3::MGResults   25.0000 60.0000 54.0000

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.