]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add two tests.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 19 Oct 2005 05:33:52 +0000 (05:33 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Wed, 19 Oct 2005 05:33:52 +0000 (05:33 +0000)
git-svn-id: https://svn.dealii.org/trunk@11624 0785d39b-7218-0410-832d-ea1e28bc413d

tests/fe/Makefile
tests/fe/rt_7.cc [new file with mode: 0644]
tests/fe/rtn_2.cc [new file with mode: 0644]

index 5c6c3b1d2832a6125bb80e55fc450c75106d4f53..5434b8becff1b09c7aa00a5a2d14b2ab38953a40 100644 (file)
@@ -30,7 +30,8 @@ tests = fe_data_test traits fe_tools fe_tools_test mapping \
         transfer internals \
        dgq_1 \
        q_1 q_2 q_3 q_4 interpolate_q1 \
-       rt_1 rt_2 rt_3 rt_4 rt_5 rt_6 rtdiff rtn_1 interpolate_rt interpolate_rtn \
+       rt_1 rt_2 rt_3 rt_4 rt_5 rt_6 rt_7 rtdiff rtn_1 rtn_2 \
+       interpolate_rt interpolate_rtn \
        dgp_monomial_1 dgp_monomial_2 \
        system_1 interpolate_system \
         non_primitive_1 non_primitive_2 \
diff --git a/tests/fe/rt_7.cc b/tests/fe/rt_7.cc
new file mode 100644 (file)
index 0000000..483eb74
--- /dev/null
@@ -0,0 +1,88 @@
+//----------------------------  rt_7.cc  ---------------------------
+//    rt_7.cc,v 1.3 2003/06/09 16:00:38 wolf Exp
+//    Version: 
+//
+//    Copyright (C) 2003, 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.
+//
+//----------------------------  rt_7.cc  ---------------------------
+
+// RT(2) had some problems with shape functions...
+
+#include "../tests.h"
+#include <base/quadrature_lib.h>
+#include <base/logstream.h>
+#include <lac/vector.h>
+#include <grid/tria.h>
+#include <grid/tria_iterator.h>
+#include <dofs/dof_accessor.h>
+#include <grid/grid_generator.h>
+#include <grid/grid_tools.h>
+#include <fe/fe_raviart_thomas.h>
+#include <fe/fe_values.h>
+
+#include <vector>
+#include <fstream>
+#include <string>
+
+#define PRECISION 2
+
+
+
+template<int dim>
+void
+plot_shape_functions(const unsigned int degree)
+{
+  FE_RaviartThomasNodal<dim> fe_rt(degree);
+  Triangulation<dim> tr;
+  GridGenerator::hyper_cube(tr, 0., 1.);
+
+  DoFHandler<dim> dof(tr);
+  typename DoFHandler<dim>::cell_iterator c = dof.begin();
+  dof.distribute_dofs(fe_rt);
+      
+  QTrapez<1> q_trapez;
+  const unsigned int div=10;
+  QIterated<dim> q(q_trapez, div);
+  FEValues<dim> fe(fe_rt, q, update_values|update_gradients|update_q_points);
+  fe.reinit(c);
+
+  Assert (fe.get_fe().n_components() == dim, ExcInternalError());
+  
+  for (unsigned int q_point=0; q_point<q.n_quadrature_points; ++q_point)
+    {
+      if (q_point % QIterated<1>(q_trapez,div).n_quadrature_points == 0)
+        deallog << std::endl;
+      
+      deallog << fe.quadrature_point(q_point) << " ";
+             
+      for (unsigned int i=0;i<fe_rt.dofs_per_cell;++i)
+        for (unsigned int c=0; c<fe.get_fe().n_components(); ++c)
+          deallog << " " << fe.shape_value_component(i,q_point,c);
+
+      deallog << std::endl;
+    }
+}
+
+
+int
+main()
+{
+  std::ofstream logfile ("rt_7.output");
+  logfile.precision (PRECISION);
+  logfile.setf(std::ios::fixed);  
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  plot_shape_functions<2>(2);
+  
+  return 0;
+}
+
+
+
diff --git a/tests/fe/rtn_2.cc b/tests/fe/rtn_2.cc
new file mode 100644 (file)
index 0000000..c9b7a18
--- /dev/null
@@ -0,0 +1,89 @@
+//----------------------------  rtn_2.cc  ---------------------------
+//    rtn_2.cc,v 1.3 2003/06/09 16:00:38 wolf Exp
+//    Version: 
+//
+//    Copyright (C) 2003, 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.
+//
+//----------------------------  rtn_2.cc  ---------------------------
+
+// rt(2) had some problems with shape functions (tested in rt_7.cc). because
+// it is so simple, just run the same test for the RT-Nodal element as well
+
+#include "../tests.h"
+#include <base/quadrature_lib.h>
+#include <base/logstream.h>
+#include <lac/vector.h>
+#include <grid/tria.h>
+#include <grid/tria_iterator.h>
+#include <dofs/dof_accessor.h>
+#include <grid/grid_generator.h>
+#include <grid/grid_tools.h>
+#include <fe/fe_raviart_thomas.h>
+#include <fe/fe_values.h>
+
+#include <vector>
+#include <fstream>
+#include <string>
+
+#define PRECISION 2
+
+
+
+template<int dim>
+void
+plot_shape_functions(const unsigned int degree)
+{
+  FE_RaviartThomasNodal<dim> fe_rt(degree);
+  Triangulation<dim> tr;
+  GridGenerator::hyper_cube(tr, 0., 1.);
+
+  DoFHandler<dim> dof(tr);
+  typename DoFHandler<dim>::cell_iterator c = dof.begin();
+  dof.distribute_dofs(fe_rt);
+      
+  QTrapez<1> q_trapez;
+  const unsigned int div=10;
+  QIterated<dim> q(q_trapez, div);
+  FEValues<dim> fe(fe_rt, q, update_values|update_gradients|update_q_points);
+  fe.reinit(c);
+
+  Assert (fe.get_fe().n_components() == dim, ExcInternalError());
+  
+  for (unsigned int q_point=0; q_point<q.n_quadrature_points; ++q_point)
+    {
+      if (q_point % QIterated<1>(q_trapez,div).n_quadrature_points == 0)
+        deallog << std::endl;
+      
+      deallog << fe.quadrature_point(q_point) << " ";
+             
+      for (unsigned int i=0;i<fe_rt.dofs_per_cell;++i)
+        for (unsigned int c=0; c<fe.get_fe().n_components(); ++c)
+          deallog << " " << fe.shape_value_component(i,q_point,c);
+
+      deallog << std::endl;
+    }
+}
+
+
+int
+main()
+{
+  std::ofstream logfile ("rtn_2.output");
+  logfile.precision (PRECISION);
+  logfile.setf(std::ios::fixed);  
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  plot_shape_functions<2>(2);
+  
+  return 0;
+}
+
+
+

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.