]> https://gitweb.dealii.org/ - dealii.git/commitdiff
A new test that checks that the order of various interpolations commute.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 16 Jan 2007 21:17:37 +0000 (21:17 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 16 Jan 2007 21:17:37 +0000 (21:17 +0000)
git-svn-id: https://svn.dealii.org/trunk@14314 0785d39b-7218-0410-832d-ea1e28bc413d

tests/fe/Makefile
tests/fe/injection_common.h [new file with mode: 0644]
tests/fe/injection_q.cc [new file with mode: 0644]

index 8c60260ac127eded46603c2ef1bffe2959aec9c7..a09a13accbd98d69dd2a024d9d17167cfc5fb12d 100644 (file)
@@ -1,6 +1,6 @@
 ############################################################
 # Makefile,v 1.15 2002/06/13 12:51:13 hartmann Exp
-# Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
+# Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007 by the deal.II authors
 ############################################################
 
 ############################################################
@@ -38,7 +38,8 @@ tests_x=fe_data_test traits fe_tools fe_tools_* mapping \
         non_primitive_* \
         nedelec* \
        up_and_down \
-       check_derivatives
+       check_derivatives \
+       injection_*
 
 # from above list of regular expressions, generate the real set of
 # tests
diff --git a/tests/fe/injection_common.h b/tests/fe/injection_common.h
new file mode 100644 (file)
index 0000000..cc6d8ea
--- /dev/null
@@ -0,0 +1,129 @@
+//----------------------------  injection_common.cc  ---------------------------
+//    $Id: injection_common.cc 12732 2006-03-28 23:15:45Z wolf $
+//    Version: $Name$ 
+//
+//    Copyright (C) 2006, 2007 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.
+//
+//----------------------------  injection_common.cc  ---------------------------
+
+
+// common framework to check the following: if we interpolate from one finite
+// element on a cell to a richer finite element on a finer cell, then it
+// shouldn't matter whether we go to the richer FE first and then to the finer
+// cells, or the other way around. test this
+
+#include "../tests.h"
+#include <base/logstream.h>
+#include <lac/full_matrix.h>
+
+#include <fe/fe_abf.h>
+#include <fe/fe_dgp.h>
+#include <fe/fe_dgp_monomial.h>
+#include <fe/fe_dgp_nonparametric.h>
+#include <fe/fe_dgq.h>
+#include <fe/fe_nedelec.h>
+#include <fe/fe_q.h>
+#include <fe/fe_q_hierarchical.h>
+#include <fe/fe_raviart_thomas.h>
+#include <fe/fe_system.h>
+
+#include <fstream>
+#include <vector>
+
+
+template <int dim>
+void test ();
+
+
+
+
+template <int dim>
+void do_check (const FiniteElement<dim> &coarse_fe,
+              const FiniteElement<dim> &fine_fe)
+{
+  FullMatrix<double> injection_1 (fine_fe.dofs_per_cell,
+                                 coarse_fe.dofs_per_cell);
+  FullMatrix<double> injection_2 (fine_fe.dofs_per_cell,
+                                 coarse_fe.dofs_per_cell);
+
+  for (unsigned int child_1=0; child_1<GeometryInfo<dim>::children_per_cell;
+       ++child_1)
+    for (unsigned int child_2=0; child_2<GeometryInfo<dim>::children_per_cell;
+        ++child_2)
+      {
+       injection_1 = 0;
+       injection_2 = 0;
+       
+                                        // check 1: first to finer fe, then
+                                        // to finer cells
+       {
+         FullMatrix<double> tmp1 (fine_fe.dofs_per_cell,
+                                  coarse_fe.dofs_per_cell);
+         FullMatrix<double> tmp2 (fine_fe.dofs_per_cell,
+                                  coarse_fe.dofs_per_cell);
+         fine_fe.get_interpolation_matrix (coarse_fe,
+                                           tmp1);
+         fine_fe.get_prolongation_matrix (child_1)
+           .mmult (tmp2, tmp1);
+         fine_fe.get_prolongation_matrix (child_2)
+           .mmult (injection_1, tmp2);
+       }
+
+                                        // check 2: first to finer cells,
+                                        // then to finer fe
+       {
+         FullMatrix<double> tmp1 (coarse_fe.dofs_per_cell,
+                                  coarse_fe.dofs_per_cell);
+         FullMatrix<double> tmp2 (fine_fe.dofs_per_cell,
+                                  coarse_fe.dofs_per_cell);
+
+         coarse_fe.get_prolongation_matrix (child_2)
+           .mmult (tmp1, coarse_fe.get_prolongation_matrix (child_1));
+
+         fine_fe.get_interpolation_matrix (coarse_fe,
+                                           tmp2);
+         tmp2.mmult (injection_2, tmp1);
+       }
+
+                                        // print one of the matrices. to
+                                        // reduce output, do so only for the
+                                        // some of the matrices
+       if (child_1 == ((child_2+1) % GeometryInfo<dim>::children_per_cell))
+         for (unsigned int i=0; i<fine_fe.dofs_per_cell; ++i)
+           for (unsigned int j=0; j<coarse_fe.dofs_per_cell; ++j)
+             deallog << i << ' ' << j << ' ' << injection_1(i,j)
+                     << std::endl;
+
+                                        // make sure that the two matrices
+                                        // are pretty much equal
+       for (unsigned int i=0; i<fine_fe.dofs_per_cell; ++i)
+         for (unsigned int j=0; j<coarse_fe.dofs_per_cell; ++j)
+           injection_2(i,j) -= injection_1(i,j);
+       Assert (injection_2.frobenius_norm() <=
+               1e-12 * injection_1.frobenius_norm(),
+               ExcInternalError());
+      }
+}
+
+
+
+
+int main ()
+{
+  std::ofstream logfile(logname);
+  logfile.precision (3);
+  
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  test<1>();
+  test<2>();
+  test<3>();
+}
+
diff --git a/tests/fe/injection_q.cc b/tests/fe/injection_q.cc
new file mode 100644 (file)
index 0000000..d063dc2
--- /dev/null
@@ -0,0 +1,29 @@
+//----------------------------  injection_q.cc  ---------------------------
+//    $Id: injection_q.cc 12732 2006-03-28 23:15:45Z wolf $
+//    Version: $Name$ 
+//
+//    Copyright (C) 2006, 2007 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.
+//
+//----------------------------  injection_q.cc  ---------------------------
+
+
+// check that computation of hp constraints works for Q elements correctly
+
+char logname[] = "injection_q/output";
+
+
+#include "injection_common.h"
+
+
+template <int dim>
+void test ()
+{
+  for (unsigned int i=1; i<4; ++i)
+    for (unsigned int j=i; j<4; ++j)
+      do_check (FE_Q<dim>(i), FE_Q<dim>(j));
+}

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.