From: Guido Kanschat Date: Thu, 14 Mar 2013 22:24:11 +0000 (+0000) Subject: local blocks are in DoFInfo now only X-Git-Tag: v8.0.0~985 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=d27b9b449c64b6982b7e8a9707e4e63b185d429e;p=dealii.git local blocks are in DoFInfo now only git-svn-id: https://svn.dealii.org/trunk@28905 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/tests/integrators/assembler_simple_matrix_01.cc b/tests/integrators/assembler_simple_matrix_01.cc index 50e180a961..43de184d5c 100644 --- a/tests/integrators/assembler_simple_matrix_01.cc +++ b/tests/integrators/assembler_simple_matrix_01.cc @@ -1,7 +1,7 @@ //---------------------------------------------------------------------- // $Id$ // -// Copyright (C) 2012 by the deal.II authors +// Copyright (C) 2012, 2013 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -12,30 +12,35 @@ /** * @file Test initialization of Assembler::MatrixSimple and - * LocalResults + * DoFInfo */ #include "../tests.h" #include #include #include +#include +#include +#include #include +#include #include +#include +#include + using namespace dealii; -template -void test(MeshWorker::Assembler::MatrixSimple& ass) -{ - MeshWorker::LocalResults info1; - ass.initialize_info(info1, false); +template +void test(DOFINFO& info, MeshWorker::Assembler::MatrixSimple& ass) +{ + ass.initialize_info(info, false); deallog << "No faces" << std::endl; - info1.print_debug(deallog); + info.print_debug(deallog); - MeshWorker::LocalResults info2; - ass.initialize_info(info2, true); + ass.initialize_info(info, true); deallog << "With faces" << std::endl; - info2.print_debug(deallog); + info.print_debug(deallog); } int main() @@ -45,23 +50,32 @@ int main() deallog.attach(logfile); deallog.depth_console (0); + Triangulation<2,2> tr; + GridGenerator::hyper_cube(tr); + FE_DGP<2,2> fe1(1); + FE_DGP<2,2> fe2(2); + FE_DGP<2,2> fe3(3); + FE_DGP<2,2> fe5(5); + FESystem<2,2> fes1(fe3,1,fe5,1,fe1,1); + FESystem<2,2> fes2(fe3,1,fe5,1,fe1,1, fe2,1); + + DoFHandler<2,2> dof1(tr); + dof1.distribute_dofs(fes1); + DoFHandler<2,2> dof2(tr); + dof2.distribute_dofs(fes2); + dof1.initialize_local_block_info(); + dof2.initialize_local_block_info(); + MeshWorker::DoFInfo<2,2,double> info1(dof1); + MeshWorker::DoFInfo<2,2,double> info1b(dof1.block_info()); + MeshWorker::DoFInfo<2,2,double> info2b(dof2.block_info()); + MeshWorker::Assembler::MatrixSimple > ass1; deallog.push("Single block"); - test(ass1); - BlockIndices ind; - ind.push_back(3); - ind.push_back(5); - ind.push_back(1); - ass1.initialize_local_blocks(ind); + test(info1, ass1); deallog.pop(); deallog.push("Multiple blocks"); - test(ass1); - ind.push_back(2); - deallog.pop(); - deallog.push("Same blocks"); - test(ass1); - ass1.initialize_local_blocks(ind); + test(info1b, ass1); deallog.pop(); deallog.push("More blocks"); - test(ass1); + test(info2b, ass1); } diff --git a/tests/integrators/assembler_simple_matrix_01/cmp/generic b/tests/integrators/assembler_simple_matrix_01/cmp/generic index 72033b0d24..6852e6de90 100644 --- a/tests/integrators/assembler_simple_matrix_01/cmp/generic +++ b/tests/integrators/assembler_simple_matrix_01/cmp/generic @@ -35,32 +35,6 @@ DEAL:Multiple blocks:: 1,2 0x0 face 1,2 0x0 DEAL:Multiple blocks:: 2,0 0x0 face 2,0 0x0 DEAL:Multiple blocks:: 2,1 0x0 face 2,1 0x0 DEAL:Multiple blocks:: 2,2 0x0 face 2,2 0x0 -DEAL:Same blocks::No faces -DEAL:Same blocks::J: 0 -DEAL:Same blocks::R: 0 -DEAL:Same blocks::M: 9 face 0 -DEAL:Same blocks:: 0,0 0x0 -DEAL:Same blocks:: 0,1 0x0 -DEAL:Same blocks:: 0,2 0x0 -DEAL:Same blocks:: 1,0 0x0 -DEAL:Same blocks:: 1,1 0x0 -DEAL:Same blocks:: 1,2 0x0 -DEAL:Same blocks:: 2,0 0x0 -DEAL:Same blocks:: 2,1 0x0 -DEAL:Same blocks:: 2,2 0x0 -DEAL:Same blocks::With faces -DEAL:Same blocks::J: 0 -DEAL:Same blocks::R: 0 -DEAL:Same blocks::M: 9 face 9 -DEAL:Same blocks:: 0,0 0x0 face 0,0 0x0 -DEAL:Same blocks:: 0,1 0x0 face 0,1 0x0 -DEAL:Same blocks:: 0,2 0x0 face 0,2 0x0 -DEAL:Same blocks:: 1,0 0x0 face 1,0 0x0 -DEAL:Same blocks:: 1,1 0x0 face 1,1 0x0 -DEAL:Same blocks:: 1,2 0x0 face 1,2 0x0 -DEAL:Same blocks:: 2,0 0x0 face 2,0 0x0 -DEAL:Same blocks:: 2,1 0x0 face 2,1 0x0 -DEAL:Same blocks:: 2,2 0x0 face 2,2 0x0 DEAL:More blocks::No faces DEAL:More blocks::J: 0 DEAL:More blocks::R: 0 diff --git a/tests/integrators/assembler_simple_mgmatrix_01.cc b/tests/integrators/assembler_simple_mgmatrix_01.cc index 9a62ca644a..d20bd8b562 100644 --- a/tests/integrators/assembler_simple_mgmatrix_01.cc +++ b/tests/integrators/assembler_simple_mgmatrix_01.cc @@ -1,7 +1,7 @@ //---------------------------------------------------------------------- // $Id$ // -// Copyright (C) 2012 by the deal.II authors +// Copyright (C) 2012, 2013 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -12,30 +12,35 @@ /** * @file Test initialization of Assembler::MatrixSimple and - * LocalResults + * DoFInfo */ #include "../tests.h" #include #include #include +#include +#include +#include #include +#include #include +#include +#include + using namespace dealii; -template -void test(MeshWorker::Assembler::MGMatrixSimple& ass) -{ - MeshWorker::LocalResults info1; - ass.initialize_info(info1, false); +template +void test(DOFINFO& info, MeshWorker::Assembler::MGMatrixSimple& ass) +{ + ass.initialize_info(info, false); deallog << "No faces" << std::endl; - info1.print_debug(deallog); + info.print_debug(deallog); - MeshWorker::LocalResults info2; - ass.initialize_info(info2, true); + ass.initialize_info(info, true); deallog << "With faces" << std::endl; - info2.print_debug(deallog); + info.print_debug(deallog); } int main() @@ -45,23 +50,32 @@ int main() deallog.attach(logfile); deallog.depth_console (0); + Triangulation<2,2> tr; + GridGenerator::hyper_cube(tr); + FE_DGP<2,2> fe1(1); + FE_DGP<2,2> fe2(2); + FE_DGP<2,2> fe3(3); + FE_DGP<2,2> fe5(5); + FESystem<2,2> fes1(fe3,1,fe5,1,fe1,1); + FESystem<2,2> fes2(fe3,1,fe5,1,fe1,1, fe2,1); + + DoFHandler<2,2> dof1(tr); + dof1.distribute_dofs(fes1); + DoFHandler<2,2> dof2(tr); + dof2.distribute_dofs(fes2); + dof1.initialize_local_block_info(); + dof2.initialize_local_block_info(); + MeshWorker::DoFInfo<2,2,double> info1(dof1); + MeshWorker::DoFInfo<2,2,double> info1b(dof1.block_info()); + MeshWorker::DoFInfo<2,2,double> info2b(dof2.block_info()); + MeshWorker::Assembler::MGMatrixSimple > ass1; deallog.push("Single block"); - test(ass1); - BlockIndices ind; - ind.push_back(3); - ind.push_back(5); - ind.push_back(1); - ass1.initialize_local_blocks(ind); + test(info1, ass1); deallog.pop(); deallog.push("Multiple blocks"); - test(ass1); - ind.push_back(2); - deallog.pop(); - deallog.push("Same blocks"); - test(ass1); - ass1.initialize_local_blocks(ind); + test(info1b, ass1); deallog.pop(); deallog.push("More blocks"); - test(ass1); + test(info2b, ass1); } diff --git a/tests/integrators/assembler_simple_mgmatrix_01/cmp/generic b/tests/integrators/assembler_simple_mgmatrix_01/cmp/generic index 72033b0d24..6852e6de90 100644 --- a/tests/integrators/assembler_simple_mgmatrix_01/cmp/generic +++ b/tests/integrators/assembler_simple_mgmatrix_01/cmp/generic @@ -35,32 +35,6 @@ DEAL:Multiple blocks:: 1,2 0x0 face 1,2 0x0 DEAL:Multiple blocks:: 2,0 0x0 face 2,0 0x0 DEAL:Multiple blocks:: 2,1 0x0 face 2,1 0x0 DEAL:Multiple blocks:: 2,2 0x0 face 2,2 0x0 -DEAL:Same blocks::No faces -DEAL:Same blocks::J: 0 -DEAL:Same blocks::R: 0 -DEAL:Same blocks::M: 9 face 0 -DEAL:Same blocks:: 0,0 0x0 -DEAL:Same blocks:: 0,1 0x0 -DEAL:Same blocks:: 0,2 0x0 -DEAL:Same blocks:: 1,0 0x0 -DEAL:Same blocks:: 1,1 0x0 -DEAL:Same blocks:: 1,2 0x0 -DEAL:Same blocks:: 2,0 0x0 -DEAL:Same blocks:: 2,1 0x0 -DEAL:Same blocks:: 2,2 0x0 -DEAL:Same blocks::With faces -DEAL:Same blocks::J: 0 -DEAL:Same blocks::R: 0 -DEAL:Same blocks::M: 9 face 9 -DEAL:Same blocks:: 0,0 0x0 face 0,0 0x0 -DEAL:Same blocks:: 0,1 0x0 face 0,1 0x0 -DEAL:Same blocks:: 0,2 0x0 face 0,2 0x0 -DEAL:Same blocks:: 1,0 0x0 face 1,0 0x0 -DEAL:Same blocks:: 1,1 0x0 face 1,1 0x0 -DEAL:Same blocks:: 1,2 0x0 face 1,2 0x0 -DEAL:Same blocks:: 2,0 0x0 face 2,0 0x0 -DEAL:Same blocks:: 2,1 0x0 face 2,1 0x0 -DEAL:Same blocks:: 2,2 0x0 face 2,2 0x0 DEAL:More blocks::No faces DEAL:More blocks::J: 0 DEAL:More blocks::R: 0 diff --git a/tests/integrators/elasticity_01.cc b/tests/integrators/elasticity_01.cc new file mode 100644 index 0000000000..b3234e4e29 --- /dev/null +++ b/tests/integrators/elasticity_01.cc @@ -0,0 +1,242 @@ +//-------------------------------------------------------------------- +// $Id$ +// +// Copyright (C) 2012, 2013 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. +// +//-------------------------------------------------------------------- + +// Test the functions in integrators/elasticity.h +// Output matrices and assert consistency of residuals + +#include "../tests.h" +#include "../lib/test_grids.h" + +#include +#include + +#include +#include +#include +#include + +using namespace LocalIntegrators::Laplace; + +template +void test_cell(const FEValuesBase& fev) +{ + const unsigned int n = fev.dofs_per_cell; + FullMatrix M(n,n); + cell_matrix(M,fev); + M.print(deallog,8); + + Vector u(n), v(n), w(n); + std::vector > > + ugrad(dim,std::vector >(fev.n_quadrature_points)); + + std::vector indices(n); + for (unsigned int i=0;i +void test_boundary(const FEValuesBase& fev) +{ + const unsigned int n = fev.dofs_per_cell; + unsigned int d=fev.get_fe().n_components(); + FullMatrix M(n,n); + nitsche_matrix(M, fev, 17); + M.print(deallog,8); + + Vector u(n), v(n), w(n); + std::vector > + uval (d,std::vector(fev.n_quadrature_points)), + null_val(d,std::vector(fev.n_quadrature_points, 0.)); + std::vector > > + ugrad (d,std::vector >(fev.n_quadrature_points)); + + std::vector indices(n); + for (unsigned int i=0;i +void test_face(const FEValuesBase& fev1, + const FEValuesBase& fev2) +{ + const unsigned int n1 = fev1.dofs_per_cell; + const unsigned int n2 = fev1.dofs_per_cell; + unsigned int d=fev1.get_fe().n_components(); + FullMatrix M11(n1,n1); + FullMatrix M12(n1,n2); + FullMatrix M21(n2,n1); + FullMatrix M22(n2,n2); + + ip_matrix(M11, M12, M21, M22, fev1, fev2, 17); + deallog << "M11" << std::endl; + M11.print(deallog,8); + deallog << "M22" << std::endl; + M22.print(deallog,8); + deallog << "M12" << std::endl; + M12.print(deallog,8); + deallog << "M21" << std::endl; + M21.print(deallog,8); + + Vector u1(n1), v1(n1), w1(n1); + Vector u2(n2), v2(n2), w2(n2); + std::vector > + u1val (d,std::vector(fev1.n_quadrature_points)), + u2val (d,std::vector(fev2.n_quadrature_points)); + std::vector > > + u1grad(d,std::vector >(fev1.n_quadrature_points)), + u2grad(d,std::vector >(fev2.n_quadrature_points)); + + std::vector indices1(n1), indices2(n2); + for (unsigned int i=0;i +void +test_fe(Triangulation& tr, FiniteElement& fe) +{ + deallog << fe.get_name() << std::endl << "cell matrix" << std::endl; + QGauss quadrature(fe.tensor_degree()+1); + FEValues fev(fe, quadrature, update_gradients); + + typename Triangulation::cell_iterator cell1 = tr.begin(1); + fev.reinit(cell1); + test_cell(fev); + + QGauss face_quadrature(fe.tensor_degree()+1); + FEFaceValues fef1(fe, face_quadrature, update_values | update_gradients | update_normal_vectors); + for (unsigned int i=0;i::faces_per_cell;++i) + { + deallog << "boundary_matrix " << i << std::endl; + fef1.reinit(cell1, i); + test_boundary(fef1); + } + + FEFaceValues fef2(fe, face_quadrature, update_values | update_gradients); + typename Triangulation::cell_iterator cell2 = cell1->neighbor(1); + + deallog << "face_matrix " << 0 << std::endl; + cell1 = tr.begin(1); + fef1.reinit(cell1, 1); + fef2.reinit(cell2, 0); + test_face(fef1, fef2); +} + + +template +void +test(Triangulation& tr) +{ + FE_DGQ q1(1); + FESystem fe1(q1,dim); + test_fe(tr, fe1); + + FE_DGQ q2(2); + FESystem fe2(q2,dim); + test_fe(tr, fe2); + + FE_Nedelec n1(1); + test_fe(tr, n1); + + FE_RaviartThomas rt1(1); + test_fe(tr, rt1); +} + + +int main() +{ + initlog(__FILE__); + deallog.threshold_double(1.e-10); + + Triangulation<2> tr2; + TestGrids::hypercube(tr2, 1); + test(tr2); + + Triangulation<3> tr3; + TestGrids::hypercube(tr3, 1); + test(tr3); +}