]> https://gitweb.dealii.org/ - dealii.git/commitdiff
new tests for FETools and FEValues
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Tue, 25 May 2004 12:41:18 +0000 (12:41 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Tue, 25 May 2004 12:41:18 +0000 (12:41 +0000)
git-svn-id: https://svn.dealii.org/trunk@9319 0785d39b-7218-0410-832d-ea1e28bc413d

tests/fe/Makefile
tests/fe/fe_tools.cc
tests/fe/function.cc [new file with mode: 0644]

index f9bb3af03146b58fbf6e6c3a46900cd53a352b21..8a70014ed71ffa4dbcbbd52670d251ae486812c0 100644 (file)
@@ -26,6 +26,7 @@ derivatives.exe         : derivatives.g.$(OBJEXT)         $(libraries)
 fe_data_test.exe        : fe_data_test.g.$(OBJEXT)        $(libraries)
 fe_tools.exe            : fe_tools.g.$(OBJEXT)            $(libraries)
 fe_tools_test.exe       : fe_tools_test.g.$(OBJEXT)       $(libraries)
+function.exe            : function.g.$(OBJEXT)            $(libraries)
 traits.exe              : traits.g.$(OBJEXT)              $(libraries)
 internals.exe           : internals.g.$(OBJEXT)           $(libraries)
 mapping.exe             : mapping.g.$(OBJEXT)             $(libraries)
@@ -56,7 +57,7 @@ transfer.exe          : transfer.g.$(OBJEXT)            $(libraries)
 up_and_down.exe         : up_and_down.g.$(OBJEXT)         $(libraries)
 
 tests = fe_data_test traits fe_tools fe_tools_test mapping \
-        mapping_c1 shapes derivatives numbering mapping_q1_eulerian \
+        mapping_c1 shapes derivatives function numbering mapping_q1_eulerian \
         transfer internals \
        dgq_1 \
        q_1 q_2 q_3 q_4 \
index ccc34c7b262ab4215d9193ca4168348485c93076..639dcba9bc4e4d62ec62cb7a4f9eb83ac8c21f65 100644 (file)
@@ -2,7 +2,7 @@
 //    fe_tools.cc,v 1.1 2003/11/28 15:03:26 guido Exp
 //    Version: 
 //
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003 by the deal.II authors
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -11,7 +11,7 @@
 //
 //--------------------------------------------------------------------------
 
-// Test the cell matrices generated in FETools.
+// Test the cell matrices generated in FETools and the local renumbering vector.
 
 #include "../tests.h"
 #include <iostream>
@@ -83,11 +83,60 @@ void test_projection (std::ostream& out)
 }
 
 
+template<int dim>
+void test_renumbering(const FiniteElement<dim>& fe)
+{
+  std::vector<unsigned int> v(fe.dofs_per_cell);
+  std::vector<std::vector<unsigned int> > start;
+  FETools::compute_component_wise(fe, v, start);
+
+  deallog << fe.get_name() << std::endl << '[';
+  for (unsigned int i=0;i<v.size();++i)
+    deallog << ' ' << v[i];
+  deallog << " ]" << std::endl;
+
+  for (unsigned int i=0;i<start.size();++i)
+    {
+      deallog << "Base " << i << ':';
+      for (unsigned int j=0;j<start[i].size();++j)
+       deallog << ' ' << start[i][j];
+      deallog << std::endl;
+    }
+  
+}
+
+
+template<int dim>
+void test_renumbering()
+{
+  deallog.push("Renumber");
+  
+  FE_Q<dim> q1(1);
+  FE_Q<dim> q3(3);
+//Todo:[GK] Test Raviart-Thomas and Nedelec?
+  
+  FESystem<dim> q1_3(q1,3);
+  test_renumbering(q1_3);
+  FESystem<dim> q3_3(q3,3);
+  test_renumbering(q3_3);
+  FESystem<dim> q3_2_q1_3(q3, 2, q1, 3);
+  test_renumbering(q3_2_q1_3);
+
+  deallog.pop();
+}
+
+
 int main()
 {
   std::ofstream logfile("fe_tools.output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
   
   test_projection<1>(logfile);
   test_projection<2>(logfile);
   test_projection<3>(logfile);
+  
+  test_renumbering<1>();
+  test_renumbering<2>();
+  test_renumbering<3>();
 }
diff --git a/tests/fe/function.cc b/tests/fe/function.cc
new file mode 100644 (file)
index 0000000..48d2768
--- /dev/null
@@ -0,0 +1,93 @@
+// shapes.cc,v 1.18 2003/04/09 15:49:55 wolf Exp
+// (c) Guido Kanschat
+
+// Test the different FEValuesBase::get_function_values
+
+#include <fstream>
+
+#include <base/quadrature_lib.h>
+#include <base/logstream.h>
+
+#include <lac/vector.h>
+
+#include <fe/fe_q.h>
+#include <fe/fe_system.h>
+#include <fe/fe_values.h>
+#include <fe/fe_tools.h>
+
+#include <dofs/dof_handler.h>
+#include <dofs/dof_renumbering.h>
+
+#include <grid/grid_generator.h>
+
+
+// Call this function with a system consisting of several copies of
+// the SAME element
+template<int dim>
+void vector_values(const FiniteElement<dim>& fe)
+{
+  Assert(fe.n_base_elements() == 1, ExcNotImplemented());
+  deallog.push(fe.get_name());
+  
+  QTrapez<dim> quadrature;
+  std::vector<unsigned int> renumbering(fe.dofs_per_cell);
+  std::vector<std::vector<unsigned int> > component_start;
+  FETools::compute_component_wise(fe, renumbering, component_start);
+
+  Triangulation<dim> tr;
+  GridGenerator::hyper_cube(tr);
+  tr.refine_global(1);
+
+  DoFHandler<dim> dof(tr);
+  dof.distribute_dofs(fe);
+  DoFRenumbering::component_wise(dof);
+  
+  Vector<float> v(dof.n_dofs());
+  for (unsigned int i=0;i<v.size();++i)
+    v(i) = i;
+  
+  FEValues<dim> feval(fe, quadrature, update_values);
+  std::vector<Vector<double> > local(quadrature.n_quadrature_points,
+                                    Vector<double>(fe.n_components()));
+  
+  typename DoFHandler<dim>::active_cell_iterator cell = dof.begin_active();
+  const typename DoFHandler<dim>::active_cell_iterator end = dof.end();
+
+  unsigned int cell_no = 0;
+  while (cell != end)
+    {
+      deallog << "Cell " << cell_no++ << std::endl;
+      feval.reinit(cell);
+      feval.get_function_values(v, local);
+      for (unsigned int c=0;c<fe.n_components();++c)
+       {
+         deallog << "Component " << c;
+         for (unsigned int k=0;k<quadrature.n_quadrature_points;++k)
+           deallog << '\t' << (int) local[k](c);
+         deallog << std::endl;
+       }
+      ++cell;
+    }
+  deallog.pop();
+}
+
+template<int dim>
+void test_vectors()
+{
+  FE_Q<dim> q1(1);
+
+  FESystem<dim> q1_3(q1,3);
+  vector_values(q1_3);
+}
+
+
+int main()
+{
+  std::ofstream logfile("function.output");
+  deallog.attach(logfile);
+//  deallog.depth_console(0);
+
+  test_vectors<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.