From 371789051fa38590e2b5972d54e8e6169a13c867 Mon Sep 17 00:00:00 2001 From: wolf Date: Fri, 14 Apr 2000 17:09:32 +0000 Subject: [PATCH] Move the generation of intergrid constraints to the library. git-svn-id: https://svn.dealii.org/trunk@2728 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/dofs/dof_tools.h | 40 ++- deal.II/deal.II/source/dofs/dof_tools.cc | 437 +++++++++++++++++++++++ 2 files changed, 466 insertions(+), 11 deletions(-) diff --git a/deal.II/deal.II/include/dofs/dof_tools.h b/deal.II/deal.II/include/dofs/dof_tools.h index 5bcb95a0a4..985cf0b36d 100644 --- a/deal.II/deal.II/include/dofs/dof_tools.h +++ b/deal.II/deal.II/include/dofs/dof_tools.h @@ -248,8 +248,9 @@ class DoFTools * considered. */ template - static void make_flux_sparsity_pattern (const DoFHandler &dof_handler, - SparsityPattern &sparsity_pattern); + static void + make_flux_sparsity_pattern (const DoFHandler &dof_handler, + SparsityPattern &sparsity_pattern); /** * Make up the constraints which @@ -353,22 +354,35 @@ class DoFTools * the number of components in the * finite element used by #dof#. */ - template - static void extract_dofs(const DoFHandler &dof_handler, - const vector &select, - vector &selected_dofs); + template + static void + extract_dofs(const DoFHandler &dof_handler, + const vector &select, + vector &selected_dofs); /** * Do the same thing as #extract_dofs# * for one level of a multi-grid DoF * numbering. */ - template - static void extract_level_dofs(const unsigned int level, - const MGDoFHandler &dof, - const vector &select, - vector &selected_dofs); + template + static void + extract_level_dofs(const unsigned int level, + const MGDoFHandler &dof, + const vector &select, + vector &selected_dofs); + + template + static void + compute_intergrid_constraints (const DoFHandler &coarse_grid, + const unsigned int coarse_component, + const DoFHandler &fine_grid, + const unsigned int fine_component, + const InterGridMap &coarse_to_fine_grid_map, + ConstraintMatrix &constraints); + + /** * Exception */ @@ -385,6 +399,10 @@ class DoFTools << "is invalid with respect to the number " << "of components in the finite element " << "(" << arg2 << ")"); + /** + * Exception + */ + DeclException0 (ExcFiniteElementsDontMatch); }; diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index e7a841a06d..1a7cceb263 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -14,6 +14,7 @@ #include #include +#include #include #include #include @@ -709,6 +710,433 @@ DoFTools::extract_level_dofs(const unsigned int level, } + + +template +void +DoFTools::compute_intergrid_constraints (const DoFHandler &coarse_grid, + const unsigned int coarse_component, + const DoFHandler &fine_grid, + const unsigned int fine_component, + const InterGridMap &coarse_to_fine_grid_map, + ConstraintMatrix &constraints) +{ + // aliases to the finite elements + // used by the dof handlers: + const FiniteElement &coarse_fe = coarse_grid.get_fe(), + &fine_fe = fine_grid.get_fe(); + + // global numbers of dofs + const unsigned int n_coarse_dofs = coarse_grid.n_dofs(), + n_fine_dofs = fine_grid.n_dofs(); + + // local numbers of dofs + const unsigned int coarse_dofs_per_cell = coarse_fe.dofs_per_cell, + fine_dofs_per_cell = fine_fe.dofs_per_cell; + + // alias the number of dofs per + // cell belonging to the + // coarse_component which is to be + // the restriction of the fine + // grid: + const unsigned int coarse_dofs_per_cell_component + = coarse_fe.base_element(coarse_fe.component_to_base(coarse_component)).dofs_per_cell; + + + // check whether component numbers + // are valid + Assert (coarse_component < coarse_fe.n_components(), + ExcInvalidComponent (coarse_component, coarse_fe.n_components())); + Assert (fine_component < fine_fe.n_components(), + ExcInvalidComponent (fine_component, fine_fe.n_components())); + // check whether respective finite + // elements are equal + Assert (coarse_fe.base_element (coarse_fe.component_to_base(coarse_component)) + == + fine_fe.base_element (fine_fe.component_to_base(fine_component)), + ExcFiniteElementsDontMatch()); + + +/* + * From here on: the term `parameter' refers to the selected component + * on the coarse grid and its analogon on the fine grid. The naming of + * variables containing this term is due to the fact that + * `selected_component' is longer, but also due to the fact that the + * code of this function was initially written for a program where the + * component which we wanted to match between grids was actually the + * `parameter' variable. + * + * Likewise, the terms `parameter grid' and `state grid' refer to the + * coarse and fine grids, respectively. + * + * Changing the names of variables would in principle be a good idea, + * but would not make things simpler and would be another source of + * errors. If anyone feels like doing so: patches would be welcome! + */ + + + + // set up vectors of cell-local + // data; each vector represents one + // degree of freedom of the + // coarse-grid variable in the + // fine-grid element + vector > parameter_dofs (coarse_dofs_per_cell_component, + Vector(fine_dofs_per_cell)); + // for each coarse dof: find its + // position within the fine element + // and set this value to one in the + // respective vector (all other values + // are zero by construction) + for (unsigned int local_coarse_dof=0; + local_coarse_dof dof_is_interesting (fine_grid.n_dofs(), false); + vector local_dof_indices (fine_fe.dofs_per_cell); + + for (DoFHandler::active_cell_iterator cell=fine_grid.begin_active(); + cell!=fine_grid.end(); ++cell) + { + cell->get_dof_indices (local_dof_indices); + for (unsigned int i=0; i global_parameter_representation (n_fine_dofs); + + // store the weights with which a dof + // on the parameter grid contributes + // to a dof on the fine grid. see the + // long doc below for more info + // + // allocate as many rows as there are + // parameter dofs on the coarse grid + // and as many columns as there are + // parameter dofs on the fine grid. + // + // weight_mapping is used to map the + // global (fine grid) parameter dof + // indices to the columns + FullMatrix weights (n_coarse_dofs, + n_parameters_on_fine_grid); + // this is this mapping. there is one + // entry for each dof on the fine grid; + // if it is a parameter dof, then its + // value is the column in weights for + // that parameter dof, if it is any + // other dof, then its value is -1, + // indicating an error + vector weight_mapping (n_fine_dofs, -1); + + // set up this mapping + if (true) + { + // list of local dof numbers + // belonging to the parameter + // part of the FE on the fine grid + vector parameter_dofs; + for (unsigned int i=0; i local_dof_indices(fine_fe.dofs_per_cell); + unsigned int next_free_index=0; + for (typename DoFHandler::active_cell_iterator cell=fine_grid.begin_active(); + cell != fine_grid.end(); ++cell) + { + cell->get_dof_indices (local_dof_indices); + for (unsigned int i=0; i::active_cell_iterator cell; + vector parameter_dof_indices (coarse_fe.dofs_per_cell); + + for (cell=coarse_grid.begin_active(); cell!=coarse_grid.end(); ++cell) + { + // get the global indices of the + // parameter dofs on this parameter + // grid cell + cell->get_dof_indices (parameter_dof_indices); + + // loop over all dofs on this + // cell and check whether they + // are interesting for us + for (unsigned int local_parameter_dof=0; + local_parameter_dof + set_dof_values_by_interpolation (parameter_dofs[local_parameter_dof], + global_parameter_representation); + // now that we've got the global + // representation of each parameter + // dof, we've only got to clobber the + // non-zero entries in that vector and + // store the result + // + // what we have learned: if entry #i# + // of the global vector holds the value + // #v[i]#, then this is the weight with + // which the present dof contributes + // to #i#. there may be several such + // #i#s and their weights' sum should + // be one. Then, #v[i]# should be + // equal to #\sum_j w_{ij} p[j]# with + // #p[j]# be the values of the degrees + // of freedom on the coarse grid. we + // can thus compute constraints which + // link the degrees of freedom #v[i]# + // on the fine grid to those on the + // coarse grid, #p[j]#. Now to use + // these as real constraints, rather + // than as additional equations, we + // have to identify representants + // among the #i# for each #j#. this will + // be done by simply taking the first + // #i# for which #w_{ij}==1#. + for (unsigned int i=0; i1) && (sum==0)), ExcInternalError()); + }; +#endif + + + // now we know that the weights in each + // row constitute a constraint. enter this + // into the constraints object + // + // first task: for each parameter dof on + // the parameter grid, find a representant + // on the fine, global grid. this is + // possible since we use conforming + // finite element. we take this representant + // to be the first element in this row with + // weight identical to one + vector representants(weights.m(), -1); + for (unsigned int parameter_dof=0; parameter_dof(column)) + break; + Assert (global_dof < weight_mapping.size(), ExcInternalError()); + + // now enter the representants global + // index into our list + representants[parameter_dof] = global_dof; + }; + + for (unsigned int global_dof=0; global_dof(global_dof))) + // dof unconstrained + continue; + + // otherwise enter all constraints + constraints.add_line (global_dof); + for (unsigned int row=first_used_row; row 1 template void @@ -755,3 +1183,12 @@ template void DoFTools::extract_level_dofs(unsigned int level, const vector& local_select, vector& selected_dofs); + +template +void +DoFTools::compute_intergrid_constraints (const DoFHandler &, + const unsigned int , + const DoFHandler &, + const unsigned int , + const InterGridMap &, + ConstraintMatrix &); -- 2.39.5