]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Test new functionality also for MappingQ2.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 12 Aug 2011 23:13:58 +0000 (23:13 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 12 Aug 2011 23:13:58 +0000 (23:13 +0000)
git-svn-id: https://svn.dealii.org/trunk@24058 0785d39b-7218-0410-832d-ea1e28bc413d

tests/fe/fe_face_values_1d_mapping_q2.cc [new file with mode: 0644]
tests/fe/fe_face_values_1d_mapping_q2/cmp/generic [new file with mode: 0644]

diff --git a/tests/fe/fe_face_values_1d_mapping_q2.cc b/tests/fe/fe_face_values_1d_mapping_q2.cc
new file mode 100644 (file)
index 0000000..44adb63
--- /dev/null
@@ -0,0 +1,188 @@
+//-------------------------------------------------------
+//    rt_11.cc,v 1.3 2003/06/09 16:00:38 wolf Exp
+//    Version: 
+//
+//    Copyright (C) 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.
+//
+//-------------------------------------------------------
+
+// like fe_face_values_1d but testing that one can use MappingQ(2)
+
+#include "../tests.h"
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/logstream.h>
+#include <deal.II/lac/vector.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/tria_boundary_lib.h>
+#include <deal.II/fe/mapping_cartesian.h>
+#include <deal.II/fe/mapping_q1.h>
+#include <deal.II/fe/mapping_q.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_values.h>
+#include <deal.II/fe/fe.h>
+#include <vector>
+#include <fstream>
+#include <iomanip>
+#include <string>
+#include <sstream>
+
+#define PRECISION 4
+
+
+
+template<int dim>
+inline void
+plot_faces(Mapping<dim> &mapping,
+          FiniteElement<dim> &fe,
+          typename DoFHandler<dim>::cell_iterator &cell)
+{
+                                  // create a QGauss<0>(4), which should
+                                  // still only have 1 quadrature point
+  QGauss<dim-1> q(4);
+  Assert (q.size() == 1, ExcInternalError());
+
+  FEFaceValues<dim> fe_values(mapping, fe, q,
+                             UpdateFlags(update_q_points
+                                         | update_JxW_values
+                                         | update_values
+                                         | update_gradients
+                                         | update_hessians
+                                         | update_normal_vectors));
+
+  for (unsigned int face_nr=0;
+       face_nr < GeometryInfo<dim>::faces_per_cell;
+       ++ face_nr)
+    {
+      deallog << "Face=" << face_nr << std::endl;
+      fe_values.reinit(cell, face_nr);
+
+                                      // print some data on the location
+      deallog << "x=" << fe_values.quadrature_point(0) << std::endl;
+      deallog << "n=" << fe_values.normal_vector(0) << std::endl;
+      deallog << "JxW=" << fe_values.JxW(0) << std::endl;
+
+                                      // now print some data on the shape
+                                      // functions
+      for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
+       deallog << "shape_function " << i << ":" << std::endl
+               << "  phi=" << fe_values.shape_value(i,0) << std::endl
+               << "  grad phi=" << fe_values.shape_grad(i,0) << std::endl
+               << "  grad^2 phi=" << fe_values.shape_hessian(i,0) << std::endl;
+
+      deallog << std::endl;
+    }
+}
+
+
+
+template<int dim>
+inline void
+plot_subfaces(Mapping<dim> &mapping,
+             FiniteElement<dim> &fe,
+             typename DoFHandler<dim>::cell_iterator &cell)
+{
+
+                                  // create a QGauss<0>(4), which should
+                                  // still only have 1 quadrature point
+  QGauss<dim-1> q(4);
+  Assert (q.size() == 1, ExcInternalError());
+
+  FESubfaceValues<dim> fe_values(mapping, fe, q,
+                                UpdateFlags(update_q_points
+                                            | update_JxW_values
+                                            | update_values
+                                            | update_gradients
+                                            | update_hessians
+                                            | update_normal_vectors));
+  for (unsigned int face_nr=0;
+       face_nr < GeometryInfo<dim>::faces_per_cell;
+       ++ face_nr)
+    for (unsigned int sub_nr=0;
+        sub_nr < GeometryInfo<dim>::max_children_per_face;
+        ++ sub_nr)
+      {
+       deallog << "Face=" << face_nr
+               << ", subface=" << sub_nr
+               << std::endl;
+       
+       fe_values.reinit(cell, face_nr, sub_nr);
+       
+                                        // print some data on the location
+       deallog << "x=" << fe_values.quadrature_point(0) << std::endl;
+       deallog << "n=" << fe_values.normal_vector(0) << std::endl;
+       deallog << "JxW=" << fe_values.JxW(0) << std::endl;
+
+                                        // now print some data on the shape
+                                        // functions
+       for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
+         deallog << "shape_function " << i << ":" << std::endl
+                 << "  phi=" << fe_values.shape_value(i,0) << std::endl
+                 << "  grad phi=" << fe_values.shape_grad(i,0) << std::endl
+                 << "  grad^2 phi=" << fe_values.shape_hessian(i,0) << std::endl;
+
+       deallog << std::endl;
+      }
+}
+
+  
+template<int dim>
+void mapping_test()
+{
+  deallog << "dim=" << dim << std::endl;
+  
+  std::vector<Mapping<dim> *> mapping_ptr;
+  std::vector<std::string> mapping_strings;
+
+  MappingQ<dim> mapping(2);
+  std::string mapping_name = "MappingQ(2)";
+
+  Triangulation<dim> tria;
+  GridGenerator::hyper_cube (tria);
+
+  FE_Q<dim> fe_q4(2);
+  DoFHandler<dim> dof(tria);
+  dof.distribute_dofs(fe_q4);      
+
+  typename DoFHandler<dim>::cell_iterator cell = dof.begin_active();
+  {
+    std::ostringstream ost;
+    ost << "MappingFace" << dim << "d-"
+       << mapping_name;
+    deallog << ost.str() << std::endl;     
+    plot_faces(mapping, fe_q4, cell);
+  }
+
+  {
+    std::ostringstream ost;
+    ost << "MappingSubface" << dim << "d-"
+       << mapping_name;
+    deallog << ost.str() << std::endl;     
+    plot_subfaces(mapping, fe_q4, cell);
+  }
+}
+
+
+
+int main()
+{
+  std::ofstream logfile ("fe_face_values_1d_mapping_q2/output");
+  deallog << std::setprecision(PRECISION);
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+  
+                                  // ----------------------- 
+                                  // Tests for dim=1
+                                  // -----------------------
+  mapping_test<1>();
+}
+
diff --git a/tests/fe/fe_face_values_1d_mapping_q2/cmp/generic b/tests/fe/fe_face_values_1d_mapping_q2/cmp/generic
new file mode 100644 (file)
index 0000000..4bb2441
--- /dev/null
@@ -0,0 +1,72 @@
+
+DEAL::dim=1
+DEAL::MappingFace1d-MappingQ(2)
+DEAL::Face=0
+DEAL::x=0.000
+DEAL::n=-1.000
+DEAL::JxW=1.000
+DEAL::shape_function 0:
+DEAL::  phi=1.000
+DEAL::  grad phi=-3.000
+DEAL::  grad^2 phi=4.000
+DEAL::shape_function 1:
+DEAL::  phi=0
+DEAL::  grad phi=-1.000
+DEAL::  grad^2 phi=4.000
+DEAL::shape_function 2:
+DEAL::  phi=0
+DEAL::  grad phi=4.000
+DEAL::  grad^2 phi=-8.000
+DEAL::
+DEAL::Face=1
+DEAL::x=1.000
+DEAL::n=1.000
+DEAL::JxW=1.000
+DEAL::shape_function 0:
+DEAL::  phi=0
+DEAL::  grad phi=1.000
+DEAL::  grad^2 phi=4.000
+DEAL::shape_function 1:
+DEAL::  phi=1.000
+DEAL::  grad phi=3.000
+DEAL::  grad^2 phi=4.000
+DEAL::shape_function 2:
+DEAL::  phi=0
+DEAL::  grad phi=-4.000
+DEAL::  grad^2 phi=-8.000
+DEAL::
+DEAL::MappingSubface1d-MappingQ(2)
+DEAL::Face=0, subface=0
+DEAL::x=0.000
+DEAL::n=-1.000
+DEAL::JxW=1.000
+DEAL::shape_function 0:
+DEAL::  phi=1.000
+DEAL::  grad phi=-3.000
+DEAL::  grad^2 phi=4.000
+DEAL::shape_function 1:
+DEAL::  phi=0
+DEAL::  grad phi=-1.000
+DEAL::  grad^2 phi=4.000
+DEAL::shape_function 2:
+DEAL::  phi=0
+DEAL::  grad phi=4.000
+DEAL::  grad^2 phi=-8.000
+DEAL::
+DEAL::Face=1, subface=0
+DEAL::x=1.000
+DEAL::n=1.000
+DEAL::JxW=1.000
+DEAL::shape_function 0:
+DEAL::  phi=0
+DEAL::  grad phi=1.000
+DEAL::  grad^2 phi=4.000
+DEAL::shape_function 1:
+DEAL::  phi=1.000
+DEAL::  grad phi=3.000
+DEAL::  grad^2 phi=4.000
+DEAL::shape_function 2:
+DEAL::  phi=0
+DEAL::  grad phi=-4.000
+DEAL::  grad^2 phi=-8.000
+DEAL::

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.