]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
add a directory with useful functions for tests
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 9 Jun 2006 07:33:39 +0000 (07:33 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 9 Jun 2006 07:33:39 +0000 (07:33 +0000)
git-svn-id: https://svn.dealii.org/trunk@13212 0785d39b-7218-0410-832d-ea1e28bc413d

tests/lib/dof_tools_frame.h [new file with mode: 0644]
tests/lib/test_grids.h [new file with mode: 0644]

diff --git a/tests/lib/dof_tools_frame.h b/tests/lib/dof_tools_frame.h
new file mode 100644 (file)
index 0000000..859389e
--- /dev/null
@@ -0,0 +1,127 @@
+//----------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$ 
+//
+//    Copyright (C) 2003, 2004, 2005, 2006 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.
+//
+//----------------------------------------------------------------------
+
+
+#include "../tests.h"
+#include <base/logstream.h>
+#include <grid/tria.h>
+#include <grid/tria_iterator.h>
+#include <grid/grid_generator.h>
+#include <dofs/dof_accessor.h>
+#include <dofs/dof_handler.h>
+#include <dofs/dof_tools.h>
+#include <fe/fe_q.h>
+#include <fe/fe_dgq.h>
+#include <fe/fe_dgp.h>
+#include <fe/fe_nedelec.h>
+#include <fe/fe_raviart_thomas.h>
+#include <fe/fe_system.h>
+
+#include <fstream>
+#include <iostream>
+#include <iomanip>
+#include <string>
+
+
+// forward declaration of the function that must be provided in the
+// .cc files
+template <int dim>
+void
+check_this (const DoFHandler<dim> &dof_handler);
+
+
+void
+output_vector (std::vector<bool> &v)
+{
+for (unsigned int i=0; i<v.size(); ++i)
+    deallog << (v[i] ? '1' : '0');
+  deallog << std::endl;
+}
+
+
+void
+output_vector (std::vector<unsigned int> &v)
+{
+for (unsigned int i=0; i<v.size(); ++i)
+  deallog << ' ' << v[i];
+  deallog << std::endl;
+}
+
+
+
+template <int dim>
+void
+check (const Triangulation<dim>& tria,
+       const FiniteElement<dim> &fe)
+{
+  deallog << fe.get_name() << std::endl;
+  
+  DoFHandler<dim> dof_handler (tria);
+  dof_handler.distribute_dofs (fe);
+
+                                   // call main function in .cc files
+  check_this (dof_handler);
+}
+
+
+template <int dim>
+void check_grid(const Triangulation<dim>& tr)
+{
+  FE_Q<dim> q1(1);
+  check(tr, q1);
+  FE_Q<dim> q2(2);
+  check(tr, q2);
+  FE_Q<dim> q3(3);
+  check(tr, q3);
+  
+  FE_DGQ<dim> dgq0(0);
+  check(tr, dgq0);
+  FE_DGQ<dim> dgq1(1);
+  check(tr, dgq1);
+  FE_DGQ<dim> dgq2(2);
+  check(tr, dgq2);
+
+  FE_DGP<dim> dgp1(1);
+  check(tr, dgp1);
+  FE_DGP<dim> dgp2(2);
+  check(tr, dgp2);
+  
+  FE_Nedelec<dim> nedelec1(1);
+  check(tr, nedelec1);
+
+  FE_RaviartThomas<dim> rt0(0);
+  check(tr, rt0);
+  FE_RaviartThomas<dim> rt1(1);
+  check(tr, rt1);
+  FE_RaviartThomas<dim> rt2(2);
+  check(tr, rt2);
+
+  FESystem<dim> s1(q1, 3);
+  check(tr, s1);
+  FESystem<dim> s2(dgq1, 2, q1, 1);
+  check(tr, s2);
+  FESystem<dim> s3(q1, 2, dgq0, 3);
+  check(tr, s3);
+  FESystem<dim> s4(q1, 3, dgq0, 2, dgp1, 1);
+  check(tr, s4);
+
+  FESystem<dim> s10(rt1, 1, dgq1, 1);
+  check(tr, s10);
+  FESystem<dim> s11(rt0, 2, rt1, 1);
+  check(tr, s11);
+
+  FESystem<dim> ss1(s1, 2, s3, 1);
+  check(tr, ss1);
+}
+
+
diff --git a/tests/lib/test_grids.h b/tests/lib/test_grids.h
new file mode 100644 (file)
index 0000000..dcd1564
--- /dev/null
@@ -0,0 +1,116 @@
+//----------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$ 
+//
+//    Copyright (C) 2006 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.
+//
+//----------------------------------------------------------------------
+
+#include <grid/tria.h>
+
+/**
+ * A set of test meshes for the deal.II test suite.
+ *
+ * These meshes exhibit certain key features for writing tests. If you
+ * want to test certain properties of algorithms, the following table
+ * might be of help.
+ *
+ * <table border=1>
+ * <tr><th>Mesh</th><th>Feature tested</th></tr>
+ * <tr><td>#hypercube(tr)</td><td>works at all on a single
+ * cell</td></tr>
+ * <tr><td>#hypercube(tr,2)</td><td>works on uniform meshes</td></tr>
+ * <tr><td>#hypercube(tr,3,true)</td><td>works with local
+ * refinement</td></tr>
+ * <tr><td>#star_shaped(tr,1)</td><td>method is robust if more than
+ * usual cells meet in one vertex</td></tr>
+ * <tr><td>#star_shaped(tr,2,true)</td><td>method is robust if more than
+ * usual cells meet in one vertex and local refinement exceeds one
+ * level</td></tr>
+ * </table>
+ *
+ * @author Guido Kanschat, 2006
+ */
+namespace TestGrids
+{
+                                  /**
+                                   * Generate grids based on
+                                   * hypercube. These meshes have a
+                                   * regular geometry and topology.
+                                   *
+                                   * @param <tt>refinement</tt>
+                                   * denotes the number of refinement
+                                   * steps of the root cell.
+                                   *
+                                   * @param if <tt>local</tt> is
+                                   * <tt>true</tt>, refine only the
+                                   * cell containing the corner with
+                                   * only negative coordinates.
+                                   */
+  template <int dim>
+  void hypercube(Triangulation<dim>& tr,
+                unsigned int refinement = 0,
+                bool local = false);
+                                  /**
+                                   * Create a star-shaped mesh,
+                                   * having more than the average
+                                   * <tt>2<sup>dim</sup></tt> cells
+                                   * in the central vertex.
+                                   *
+                                   * @param <tt>refinement</tt>
+                                   * denotes the number of refinement
+                                   * steps of the root mesh.
+                                   *
+                                   * @param if <tt>local</tt> is
+                                   * <tt>true</tt>, refine only one
+                                   * of the coarse cells.
+                                   */
+  template <int dim>
+  void star_shaped(Triangulation<dim>& tr,
+                  unsigned int refinement = 0,
+                  bool local = false);
+                                  /**
+                                   * Local refinement of every other
+                                   * cell in a checkerboard fashion.
+                                   */
+  template <int dim>
+  void checkers(Triangulation<dim>& tr);
+                                  /**
+                                   * Islands of local refinement
+                                   */
+  template <int dim>
+  void islands(Triangulation<dim>& tr);
+                                  /**
+                                   * Local refinement with an
+                                   * unrefined hole.
+                                   */
+  template <int dim>
+  void laguna(Triangulation<dim>& tr);
+
+
+  template <int dim>
+  void hypercube(Triangulation<dim>& tr,
+                unsigned int refinement,
+                bool local)
+  {
+    GridGenerator::hyper_cube(tr);
+    if (!local)
+      {
+       if (refinement > 0)
+         tr.refine_global(refinement);
+      }
+    else
+      {
+       for (unsigned int i=0;i<refinement;++i)
+         {
+           tr.begin_active()->set_refine_flag();
+           tr.execute_coarsening_and_refinement();
+         }
+      }
+  }
+}

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.