From 4e7001dacd569a58c940580083a7ead7c3027b17 Mon Sep 17 00:00:00 2001 From: Guido Kanschat Date: Mon, 30 May 2005 09:08:21 +0000 Subject: [PATCH] add tests for new nodal RT element git-svn-id: https://svn.dealii.org/trunk@10775 0785d39b-7218-0410-832d-ea1e28bc413d --- tests/fe/Makefile | 3 +- tests/fe/fe_data_test.cc | 21 +++++- tests/fe/rtdiff.cc | 150 +++++++++++++++++++++++++++++++++++++++ tests/fe/rtn_1.cc | 4 +- 4 files changed, 174 insertions(+), 4 deletions(-) create mode 100644 tests/fe/rtdiff.cc diff --git a/tests/fe/Makefile b/tests/fe/Makefile index 150a22d78c..02b1a13a1e 100644 --- a/tests/fe/Makefile +++ b/tests/fe/Makefile @@ -49,6 +49,7 @@ rt_4.exe : rt_4.g.$(OBJEXT) $(libraries) rt_5.exe : rt_5.g.$(OBJEXT) $(libraries) rt_6.exe : rt_6.g.$(OBJEXT) $(libraries) rtn_1.exe : rtn_1.g.$(OBJEXT) $(libraries) +rtdiff.exe : rtdiff.g.$(OBJEXT) $(libraries) dgp_monomial_1.exe : dgp_monomial_1.g.$(OBJEXT) $(libraries) dgp_monomial_2.exe : dgp_monomial_2.g.$(OBJEXT) $(libraries) system_1.exe : system_1.g.$(OBJEXT) $(libraries) @@ -61,7 +62,7 @@ tests = fe_data_test traits fe_tools fe_tools_test mapping \ transfer internals \ dgq_1 \ q_1 q_2 q_3 q_4 \ - rt_1 rt_2 rt_3 rt_4 rt_5 rt_6 \ + rt_1 rt_2 rt_3 rt_4 rt_5 rt_6 rtdiff rtn_1 \ dgp_monomial_1 dgp_monomial_2 \ system_1 \ non_primitive_1 non_primitive_2 \ diff --git a/tests/fe/fe_data_test.cc b/tests/fe/fe_data_test.cc index add70c3e63..072a957536 100644 --- a/tests/fe/fe_data_test.cc +++ b/tests/fe/fe_data_test.cc @@ -59,12 +59,29 @@ void test_fe_datas() // Check Raviart-Thomas in 2d only if (dim==2) { - fe_datas.push_back(new FE_RaviartThomas(0)); + FE_RaviartThomas* rt0 = new FE_RaviartThomas(0); + FE_RaviartThomas* rt1 = new FE_RaviartThomas(1); + fe_datas.push_back(rt0); deallog << (*fe_datas.rbegin())->get_name() << std::endl; - fe_datas.push_back(new FE_RaviartThomas(1)); + fe_datas.push_back(rt1); deallog << (*fe_datas.rbegin())->get_name() << std::endl; fe_datas.push_back(new FE_RaviartThomas(2)); deallog << (*fe_datas.rbegin())->get_name() << std::endl; + fe_datas.push_back(new FESystem(*rt1, 1, + FE_DGQ (1), 1)); + deallog << (*fe_datas.rbegin())->get_name() << std::endl; + } + if (dim>1) + { + FE_RaviartThomasNodal* rt0 = new FE_RaviartThomasNodal(0); + FE_RaviartThomasNodal* rt1 = new FE_RaviartThomasNodal(1); + fe_datas.push_back(rt0); + deallog << (*fe_datas.rbegin())->get_name() << std::endl; + fe_datas.push_back(rt1); + deallog << (*fe_datas.rbegin())->get_name() << std::endl; + fe_datas.push_back(new FESystem(*rt1, 1, + FE_DGQ (1), 1)); + deallog << (*fe_datas.rbegin())->get_name() << std::endl; } diff --git a/tests/fe/rtdiff.cc b/tests/fe/rtdiff.cc new file mode 100644 index 0000000000..d8bd721ae1 --- /dev/null +++ b/tests/fe/rtdiff.cc @@ -0,0 +1,150 @@ +//---------------------------------------------------------------------- +// rt_1.cc,v 1.3 2003/06/09 16:00:38 wolf Exp +// Version: +// +// 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. +// +//---------------------------------------------------------------------- + +// FE_RaviartThomasNodal and FE_RaviartThomas + +// compare the shape funcions and shape values after converting to the +// same basis. + +// Summary: the different Raviart-Thomas implementations use the same +// polynomial spaces, but different basis functions. Here, we convert +// between the bases and test if the resulting functions are the same +// point-wise. + +#include "../tests.h" +#include +#include +#include +#include +#include +#include + +#include + +// This function copied from FERaviartThomasNodal. nodes is the +// element having the support points and the value of other in these +// points is computed. +template +void +initialize_node_matrix (const FiniteElement& other, + const FiniteElement& nodes, + FullMatrix& N) +{ + const unsigned int n_dofs = other.dofs_per_cell; + Assert (n_dofs == nodes.dofs_per_cell, + ExcDimensionMismatch(n_dofs, nodes.dofs_per_cell)); + + N.reinit(n_dofs, n_dofs); + + const std::vector >& unit_support_points = nodes.get_unit_support_points(); + + // The curent node functional index + unsigned int current = 0; + + // For each face and all quadrature + // points on it, the node value is + // the normal component of the + // shape function, possibly + // pointing in negative direction. + for (unsigned int face = 0; face::faces_per_cell; ++face) + for (unsigned int k=0;k::unit_normal_direction[face]); + ++current; + } + // Interior degrees of freedom in each direction + const unsigned int n_cell = (n_dofs - current) / dim; + + for (unsigned int d=0;d +void +compare_shapes (const FiniteElement& other, + const FiniteElement& nodes, + FullMatrix& M) +{ + QGauss quadrature(other.degree+1); + Table<3,double> other_values(quadrature.n_quadrature_points, other.dofs_per_cell, dim); + Table<3,double> nodes_values(quadrature.n_quadrature_points, other.dofs_per_cell, dim); + Table<3,Tensor<1,dim> > other_grads(quadrature.n_quadrature_points, other.dofs_per_cell, dim); + Table<3,Tensor<1,dim> > nodes_grads(quadrature.n_quadrature_points, other.dofs_per_cell, dim); + for (unsigned int k=0;k grad = other_grads[k][i][d]; + for (unsigned int j=0;j 1.e-12) + deallog << "Error value\t" << k << '\t' << i << '\t' << d << '\t' << value + << std::endl; + if (grad.norm() > 1.e-12) + deallog << "Error grad\t" << k << '\t' << i << '\t' << d << '\t' << grad + << '\t' << other_grads[k][i][d] + << std::endl; + } + deallog << std::endl; + } +} + + +template +void +test (unsigned int degree) +{ + FE_RaviartThomas rt1(degree); + FE_RaviartThomasNodal rtn1(degree); + FullMatrix N; + initialize_node_matrix(rt1, rtn1, N); + compare_shapes(rt1, rtn1, N); +} + +int main() +{ + std::ofstream logfile ("rtdiff.output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + test<2>(0); + test<2>(1); + test<2>(2); +} diff --git a/tests/fe/rtn_1.cc b/tests/fe/rtn_1.cc index 75f7ce463a..2fbb1fba7a 100644 --- a/tests/fe/rtn_1.cc +++ b/tests/fe/rtn_1.cc @@ -138,11 +138,13 @@ main() { std::ofstream logfile ("rtn_1.output"); deallog.attach(logfile); - deallog.depth_console(10); + deallog.depth_console(0); deallog.threshold_double(1.e-10); FE_RaviartThomasNodal<2> e20(0); check_support_points(e20); + FE_RaviartThomasNodal<3> e30(0); + check_support_points(e30); check_face_support_points(e20); FE_RaviartThomasNodal<2> e21(1); check_support_points(e21); -- 2.39.5