From 43c3992a8d29956181f196c8d0892d54339df8c0 Mon Sep 17 00:00:00 2001 From: guido Date: Tue, 4 Jul 2006 05:49:15 +0000 Subject: [PATCH] working on mg boundary git-svn-id: https://svn.dealii.org/trunk@13329 0785d39b-7218-0410-832d-ea1e28bc413d --- .../deal.II/source/multigrid/mg_dof_tools.cc | 204 ++++++++++++++++++ tests/bits/Makefile | 1 + tests/bits/dof_tools_01c.cc | 29 +-- tests/dofs/Makefile | 59 +++++ tests/dofs/row_length_01.cc | 95 ++++++++ tests/multigrid/Makefile | 5 +- tests/multigrid/boundary_01.cc | 115 ++++++++++ tests/multigrid/boundary_01/.cvsignore | 1 + 8 files changed, 495 insertions(+), 14 deletions(-) create mode 100644 tests/dofs/Makefile create mode 100644 tests/dofs/row_length_01.cc create mode 100644 tests/multigrid/boundary_01.cc create mode 100644 tests/multigrid/boundary_01/.cvsignore diff --git a/deal.II/deal.II/source/multigrid/mg_dof_tools.cc b/deal.II/deal.II/source/multigrid/mg_dof_tools.cc index 9574b494bc..124e7b3450 100644 --- a/deal.II/deal.II/source/multigrid/mg_dof_tools.cc +++ b/deal.II/deal.II/source/multigrid/mg_dof_tools.cc @@ -1017,7 +1017,205 @@ MGTools::count_dofs_per_component ( } +#if deal_II_dimension == 1 + +template +void +MGTools::make_boundary_list( + const MGDoFHandler&, + const typename FunctionMap::type&, + std::vector >&, + const std::vector&) +{ + Assert(false, ExcNotImplemented()); +} + +#else + +template +void +MGTools::make_boundary_list( + const MGDoFHandler& dof, + const typename FunctionMap::type& function_map, + std::vector >& boundary_indices, + const std::vector& component_mask_) +{ + // if for whatever reason we were + // passed an empty map, return + // immediately + if (function_map.size() == 0) + return; + + const unsigned int n_components = DoFTools::n_components(dof); + const bool fe_is_system = (n_components != 1); + +// for (typename FunctionMap::type::const_iterator i=function_map.begin(); +// i!=function_map.end(); ++i) +// Assert (n_components == i->second->n_components, +// ExcInvalidFE()); + + // set the component mask to either + // the original value or a vector + // of @p{true}s + const std::vector component_mask ((component_mask_.size() == 0) ? + std::vector (n_components, true) : + component_mask_); +// Assert (std::count(component_mask.begin(), component_mask.end(), true) > 0, +// ExcComponentMismatch()); + + // field to store the indices + std::vector face_dofs; + face_dofs.reserve (DoFTools::max_dofs_per_face(dof)); + std::fill (face_dofs.begin (), face_dofs.end (), DoFHandler::invalid_dof_index); + + typename MGDoFHandler::cell_iterator cell = dof.begin(), + endc = dof.end(); + for (; cell!=endc; ++cell) + for (unsigned int face_no = 0; face_no < GeometryInfo::faces_per_cell; + ++face_no) + { + const FiniteElement &fe = cell->get_fe(); + const unsigned int level = cell->level(); + + // we can presently deal only with + // primitive elements for boundary + // values. this does not preclude + // us using non-primitive elements + // in components that we aren't + // interested in, however. make + // sure that all shape functions + // that are non-zero for the + // components we are interested in, + // are in fact primitive + for (unsigned int i=0; iget_fe().dofs_per_cell; ++i) + { + const std::vector &nonzero_component_array + = cell->get_fe().get_nonzero_components (i); + for (unsigned int c=0; cget_fe().is_primitive (i), + ExcMessage ("This function can only deal with requested boundary " + "values that correspond to primitive (scalar) base " + "elements")); + } + + typename MGDoFHandler::face_iterator face = cell->face(face_no); + const unsigned char boundary_component = face->boundary_indicator(); + if (function_map.find(boundary_component) != function_map.end()) + // face is of the right component + { + // get indices, physical location and + // boundary values of dofs on this + // face + face_dofs.resize (fe.dofs_per_face); + face->get_mg_dof_indices (face_dofs); + if (fe_is_system) + { + // enter those dofs + // into the list that + // match the + // component + // signature. avoid + // the usual + // complication that + // we can't just use + // *_system_to_component_index + // for non-primitive + // FEs + for (unsigned int i=0; i void MGTools::reinit_vector (const MGDoFHandler &mg_dof, @@ -1218,3 +1416,9 @@ template void MGTools::count_dofs_per_component ( const MGDoFHandler&, std::vector >&, std::vector); + +template void MGTools::make_boundary_list( + const MGDoFHandler&, + const FunctionMap::type&, + std::vector >&, + const std::vector&); diff --git a/tests/bits/Makefile b/tests/bits/Makefile index e0522f2222..52042bb456 100644 --- a/tests/bits/Makefile +++ b/tests/bits/Makefile @@ -32,6 +32,7 @@ tests_x = geometry_info_* \ data_out_rotation_0* \ data_out_stack_0* \ dof_tools_[0-9]* \ + dof_tools_row_[0-9] \ fe_tools_[0-9]* \ fe_tools_cpfqpm* \ anna_? \ diff --git a/tests/bits/dof_tools_01c.cc b/tests/bits/dof_tools_01c.cc index 7a402f7ccd..a6e3d1d604 100644 --- a/tests/bits/dof_tools_01c.cc +++ b/tests/bits/dof_tools_01c.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2003, 2004 by the deal.II authors +// Copyright (C) 2003, 2004, 2006 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -38,17 +38,22 @@ check_this (const DoFHandler &dof_handler) // create sparsity pattern const unsigned int n_components = dof_handler.get_fe().n_components(); - BlockSparsityPattern sp (n_components, - n_components); - std::vector dofs_per_component(n_components); - DoFTools::count_dofs_per_component (dof_handler, - dofs_per_component); - for (unsigned int i=0; i dofs_per_component(n_blocks); + DoFTools::count_dofs_per_block (dof_handler, + dofs_per_component); + BlockIndices block_indices(dofs_per_component); + + std::vector > + row_length(n_blocks, std::vector(dof_handler.n_dofs())); + Table<2,DoFTools::Coupling> couplings(n_components, n_components); + couplings.fill(DoFTools::always); + DoFTools::compute_row_length_vector(dof_handler, row_length, + couplings, couplings); + sp.reinit(block_indices, block_indices, row_length); DoFTools::make_sparsity_pattern (dof_handler, sp); sp.compress (); diff --git a/tests/dofs/Makefile b/tests/dofs/Makefile new file mode 100644 index 0000000000..e169d15021 --- /dev/null +++ b/tests/dofs/Makefile @@ -0,0 +1,59 @@ +############################################################ +# 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 +############################################################ + +############################################################ +# Include general settings for including DEAL libraries +############################################################ + +include ../Makefile.paths +include $D/common/Make.global_options +debug-mode = on + +libraries = $(lib-deal2-1d.g) \ + $(lib-deal2-2d.g) \ + $(lib-deal2-3d.g) \ + $(lib-lac.g) \ + $(lib-base.g) + +default: run-tests + +############################################################ + +tests_x = row_length_* + +# tests for the hp branch: +# fe_collection_* + +ifeq ($(USE_CONTRIB_PETSC),yes) + tests_x += petsc_* +endif + +ifeq ($(USE_CONTRIB_METIS),yes) + tests_x += metis_* +endif + +ifeq ($(USE_CONTRIB_HSL),yes) + tests_x += hsl_* +endif + +ifeq ($(USE_CONTRIB_UMFPACK),yes) + tests_x += umfpack_* +endif + + +# from above list of regular expressions, generate the real set of +# tests +expand = $(shell echo $(addsuffix .cc,$(1)) \ + | $(PERL) -pi -e 's/\.cc//g;') +tests = $(call expand,$(tests_x)) + +############################################################ + + +include ../Makefile.rules +include Makefile.depend +include Makefile.tests + +.PHONY: default diff --git a/tests/dofs/row_length_01.cc b/tests/dofs/row_length_01.cc new file mode 100644 index 0000000000..9774b8e8df --- /dev/null +++ b/tests/dofs/row_length_01.cc @@ -0,0 +1,95 @@ +//---------------------------------------------------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2003, 2004, 2005, 2006 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 "../lib/dof_tools_frame.h" +#include "../lib/test_grids.h" + +#include +#include + +template +void +check_this (const DoFHandler &dof) +{ + const unsigned int n = dof.n_dofs(); + deallog << "dofs: " << n + << "\tcell: " << dof.get_fe().dofs_per_cell + << std::endl; + + std::vector row_length(n); + DoFTools::compute_row_length_vector(dof, row_length); + SparsityPattern sparsity(n, n, row_length); + DoFTools::make_sparsity_pattern(dof, sparsity); + sparsity.compress(); + unsigned int sum = std::accumulate(row_length.begin(), row_length.end(), 0U); + deallog << std::endl << "Entries estimated/actual " << sum + << '/' << sparsity.n_nonzero_elements() << std::endl; +// sparsity.clear(); + output_vector(row_length); + + DoFTools::compute_row_length_vector(dof, row_length, DoFTools::always); + output_vector(row_length); +} + + +template +void check() +{ + Triangulation tr; + TestGrids::hypercube(tr); + deallog << "cube" << dim << std::endl; + check_grid(tr); + tr.refine_global(1); + deallog << "refined cube" << dim << std::endl; + check_grid(tr); +} + + +int +main() +{ + try + { + std::ofstream logfile("row_length_01/output"); + logfile.precision (2); + deallog.attach(logfile); +// deallog.depth_console(0); + deallog.threshold_double(1.e-10); + check<2>(); + return 0; + } + catch (std::exception &exc) + { + std::cerr << std::endl << std::endl + << "----------------------------------------------------" + << std::endl; + std::cerr << "Exception on processing: " << std::endl + << exc.what() << std::endl + << "Aborting!" << std::endl + << "----------------------------------------------------" + << std::endl; + return 1; + } + catch (...) + { + std::cerr << std::endl << std::endl + << "----------------------------------------------------" + << std::endl; + std::cerr << "Unknown exception!" << std::endl + << "Aborting!" << std::endl + << "----------------------------------------------------" + << std::endl; + return 1; + }; +} + diff --git a/tests/multigrid/Makefile b/tests/multigrid/Makefile index 06559e9d0a..19074c3432 100644 --- a/tests/multigrid/Makefile +++ b/tests/multigrid/Makefile @@ -21,8 +21,9 @@ default: run-tests ############################################################ -tests_x = cycles transfer transfer_select transfer_block \ - smoother_block count_* renumbering_* +tests_x = cycles count_* boundary_* renumbering_* \ + transfer transfer_select transfer_block \ + smoother_block # from above list of regular expressions, generate the real set of # tests diff --git a/tests/multigrid/boundary_01.cc b/tests/multigrid/boundary_01.cc new file mode 100644 index 0000000000..50693e9b38 --- /dev/null +++ b/tests/multigrid/boundary_01.cc @@ -0,0 +1,115 @@ +//---------------------------------------------------------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2006 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. +// +//---------------------------------------------------------------------------- + +// check MGTools::make_boundary_list + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include + +using namespace std; + +void log_vector (const std::vector >& count) +{ + for (unsigned int l=0;l +void face_levels(const MGDoFHandler& dof) +{ + typename MGDoFHandler::face_iterator face; + const typename MGDoFHandler::face_iterator end = dof.end_face(); + + for (face = dof.begin_face(); face != end; ++face) + { + deallog << "Face " << face->level(); + for (unsigned int i=0;i::vertices_per_face;++i) + deallog << " v" << face->vertex(i); + deallog << " dofs "; + deallog << std::endl; + } +} + + +template +void check_fe(FiniteElement& fe) +{ + deallog << fe.get_name() << std::endl; + + Triangulation tr; + GridGenerator::hyper_cube(tr); + tr.refine_global(2); + ZeroFunction zero; + typename FunctionMap::type fmap; + fmap.insert(std::make_pair(0, &zero)); + + MGDoFHandler mgdof(tr); + mgdof.distribute_dofs(fe); + face_levels(mgdof); + + std::vector > boundary_indices(tr.n_levels()); + MGTools::make_boundary_list(mgdof, fmap, boundary_indices); + log_vector(boundary_indices); +} + + +template +void check() +{ + FE_Q q1(1); + FE_Q q2(2); +// FE_DGQ dq1(1); + + FESystem s1(q1, 2, q2,1); + + check_fe(q1); + check_fe(q2); + check_fe(s1); +} + +int main() +{ + std::ofstream logfile("boundary_01/output"); + logfile.precision(3); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + check<2> (); + check<3> (); +} diff --git a/tests/multigrid/boundary_01/.cvsignore b/tests/multigrid/boundary_01/.cvsignore new file mode 100644 index 0000000000..d196dfb299 --- /dev/null +++ b/tests/multigrid/boundary_01/.cvsignore @@ -0,0 +1 @@ +obj.* exe OK output -- 2.39.5