]> https://gitweb.dealii.org/ - dealii.git/commitdiff
tests for RTN and RT
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Wed, 19 Oct 2005 04:56:45 +0000 (04:56 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Wed, 19 Oct 2005 04:56:45 +0000 (04:56 +0000)
git-svn-id: https://svn.dealii.org/trunk@11622 0785d39b-7218-0410-832d-ea1e28bc413d

tests/fe/Makefile
tests/fe/interpolate_rt.cc [new file with mode: 0644]
tests/fe/interpolate_rtn.cc

index 1ebe5c3f2fd70508d7eca8abcf614aa139d7bc39..5c6c3b1d2832a6125bb80e55fc450c75106d4f53 100644 (file)
@@ -30,7 +30,7 @@ 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_rtn \
+       rt_1 rt_2 rt_3 rt_4 rt_5 rt_6 rtdiff rtn_1 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/interpolate_rt.cc b/tests/fe/interpolate_rt.cc
new file mode 100644 (file)
index 0000000..4d61217
--- /dev/null
@@ -0,0 +1,73 @@
+//----------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$ 
+//
+//    Copyright (C) 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.
+//
+//----------------------------------------------------------------------
+
+#include "interpolate_common.h"
+#include <base/logstream.h>
+#include <base/quadrature_lib.h>
+
+#include <fe/fe_raviart_thomas.h>
+
+#include <fstream>
+
+// FE_RaviartThomas<dim>::interpolate(...)
+
+template <int dim>
+void check1(const Function<dim>& f,
+           const unsigned int degree)
+{
+  FE_RaviartThomas<dim> fe(degree);
+  deallog << fe.get_name() << ' ';
+  deallog << fe.get_generalized_support_points().size() << ' ';
+  
+  std::vector<double> dofs(fe.dofs_per_cell);
+  
+  std::vector<std::vector<double> >
+    values(dim, std::vector<double>(fe.get_generalized_support_points().size()));
+  std::vector<Vector<double> >
+    vectors(fe.get_generalized_support_points().size(),
+           Vector<double>(dim));
+  f.vector_value_list(fe.get_generalized_support_points(), vectors);
+
+  for (unsigned int c=0;c<values.size();++c)
+    for (unsigned int k=0;k<values[c].size();++k)
+      values[c][k] = vectors[k](c);
+  
+  fe.interpolate(dofs, values);
+  deallog << " vector " << vector_difference(fe,dofs,f,0);
+  
+  fe.interpolate(dofs, vectors, 0);
+  deallog << " Vector " << vector_difference(fe,dofs,f,0) << std::endl;
+}
+
+int main()
+{
+  std::ofstream logfile ("interpolate_rt.output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-15);
+
+//   Q1WedgeFunction<1,1,2> w1;
+//   check1(w1,1,2);
+//   check1(w1,2,2);
+//   check1(w1,3,2);
+  Q1WedgeFunction<2,1,2> w21;
+  check1(w21,1);
+  check1(w21,2);
+  Q1WedgeFunction<2,2,2> w22;
+  check1(w22,2);
+  Q1WedgeFunction<2,3,2> w23;
+  check1(w23,3);
+//  Q1WedgeFunction<3,1,3> w3;
+//  check1(w3,1);
+//  check1(w3,2);
+}
index f22e6e52652757bb2bb617937b2807adfe6dad27..f18bf86caa36d48840142c4f87682e52d25d733c 100644 (file)
@@ -54,18 +54,20 @@ int main()
   std::ofstream logfile ("interpolate_rtn.output");
   deallog.attach(logfile);
   deallog.depth_console(0);
-  deallog.threshold_double(1.e-15);
+  deallog.threshold_double(1.e-13);
 
-//   Q1WedgeFunction<1,1,2> w1;
-//   check1(w1,1,2);
-//   check1(w1,2,2);
-//   check1(w1,3,2);
+//    Q1WedgeFunction<1,1,2> w1;
+//    check1(w1,1,2);
+//    check1(w1,2,2);
+//    check1(w1,3,2);
   Q1WedgeFunction<2,1,2> w2;
   check1(w2,1);
   check1(w2,2);
   Q1WedgeFunction<2,2,2> w22;
   check1(w22,2);
-//  Q1WedgeFunction<3,1,3> w3;
-//  check1(w3,1);
-//  check1(w3,2);
+  Q1WedgeFunction<2,3,2> w23;
+  check1(w23,3);
+  Q1WedgeFunction<3,1,3> w3;
+  check1(w3,1);
+  check1(w3,2);
 }

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.