From: Guido Kanschat Date: Mon, 12 Oct 2009 17:27:12 +0000 (+0000) Subject: add first assembler test X-Git-Tag: v8.0.0~6910 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=d94c7bbba8af4b9b4ce7a41d040b77fd7ae816fe;p=dealii.git add first assembler test git-svn-id: https://svn.dealii.org/trunk@19826 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/tests/deal.II/Makefile b/tests/deal.II/Makefile index 4e19b4e270..6d2d3b7ac2 100644 --- a/tests/deal.II/Makefile +++ b/tests/deal.II/Makefile @@ -75,6 +75,7 @@ tests_x = block_info block_matrices \ iterators* \ mesh_smoothing_* \ create_* \ + mesh_worker_* \ line_coarsening_3d \ no_flux_* \ fe_values_view_* \ diff --git a/tests/deal.II/mesh_worker_01.cc b/tests/deal.II/mesh_worker_01.cc new file mode 100644 index 0000000000..a8101c9423 --- /dev/null +++ b/tests/deal.II/mesh_worker_01.cc @@ -0,0 +1,225 @@ +//---------------------------------------------------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2000, 2001, 2003, 2004, 2007, 2008, 2009 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 whether the various assember classes put the right data in the +// right place. + +#include "../tests.h" +#include +#include +#include + +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include + +using namespace dealii; + + +// Define a class that fills all available entries in the info objects +// with recognizable numbers. +template +class Local : public Subscriptor +{ + public: + typedef typename MeshWorker::IntegrationWorker::CellInfo CellInfo; + typedef typename MeshWorker::IntegrationWorker::FaceInfo FaceInfo; + + void cell(CellInfo& info) const; + void bdry(FaceInfo& info) const; + void face(FaceInfo& info1, FaceInfo& info2) const; + + bool cells; + bool faces; +}; + +// Fill local structures with the following format: +// 1. Residuals: CC.DDD +// 2. Matrices: CC.DDDEEE + +// CC : 2 digits number of cell +// +// DDD: degree of freedom of test function in cell, 1 digit block, 2 +// digits number in block +// EEE: degree of freedom of trial function in cell, same format as +// DDD + +template +void +Local::cell(CellInfo& info) const +{ + if (!cells) return; + const unsigned int cell = info.cell->user_index(); + + deallog << "Cell " << std::setw(2) << cell << ':'; + for (unsigned int i=0;i& M1 = info.M1[k].matrix; + for (unsigned int i=0;i +void +Local::bdry(FaceInfo& info) const +{ + const unsigned int cell = info.cell->user_index(); + deallog << "Bdry " << std::setw(2) << cell; + deallog << std::endl; +} + + +template +void +Local::face(FaceInfo& info1, FaceInfo& info2) const +{ + const unsigned int cell1 = info1.cell->user_index(); + const unsigned int cell2 = info2.cell->user_index(); + deallog << "Face " << cell1 << '|' << cell2; + deallog << std::endl; +} + + +template +void +assemble(const DoFHandler& dof_handler, WORKER& worker) +{ + const FiniteElement& fe = dof_handler.get_fe(); + MappingQ1 mapping; + + MeshWorker::IntegrationInfoBox info_box(dof_handler); + info_box.initialize(worker, fe, mapping); + + MeshWorker::integration_loop(dof_handler.begin_active(), dof_handler.end(), info_box, worker); +} + + +template +void +test_simple(MGDoFHandler& mgdofs) +{ + SparsityPattern pattern; + SparseMatrix matrix; + Vector v; + + const DoFHandler& dofs = mgdofs; + const FiniteElement& fe = dofs.get_fe(); + pattern.reinit (dofs.n_dofs(), dofs.n_dofs(), + (GeometryInfo::faces_per_cell + *GeometryInfo::max_children_per_face+1)*fe.dofs_per_cell); + DoFTools::make_flux_sparsity_pattern (dofs, pattern); + pattern.compress(); + matrix.reinit (pattern); + v.reinit (dofs.n_dofs()); + + Local local; + local.cells = true; + local.faces = false; + + MeshWorker::AssemblingIntegrator, Vector >, Local > + integrator(local); + integrator.initialize_gauss_quadrature(1, 1, 1); + integrator.initialize(matrix, v); + integrator.boundary_fluxes = local.faces; + integrator.interior_fluxes = local.faces; + + assemble(dofs, integrator); + for (unsigned int i=0;i +void +test(const FiniteElement& fe) +{ + Triangulation tr; + GridGenerator::hyper_cube(tr); + tr.refine_global(1); + tr.begin(1)->set_refine_flag(); + tr.execute_coarsening_and_refinement(); + tr.begin(2)->set_refine_flag(); + tr.execute_coarsening_and_refinement(); +// tr.refine_global(1); + deallog << "Triangulation levels"; + for (unsigned int l=0;l::cell_iterator cell = tr.begin(); + cell != tr.end(); ++cell, ++cn) + cell->set_user_index(cn); + + MGDoFHandler dofs(tr); + dofs.distribute_dofs(fe); + deallog << "DoFHandler " << dofs.n_dofs() << " levels"; + for (unsigned int l=0;l > > fe2; + fe2.push_back(boost::shared_ptr >(new FE_DGP<2>(1))); + fe2.push_back(boost::shared_ptr >(new FE_Q<2>(1))); + + for (unsigned int i=0;i