From: kronbichler Date: Wed, 22 Oct 2008 07:48:42 +0000 (+0000) Subject: Added the possibility to omit constrained dofs. Apply this to step-31 and step-32. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=4502f734844e397cfb3a433ad853c96da41ed0ad;p=dealii-svn.git Added the possibility to omit constrained dofs. Apply this to step-31 and step-32. git-svn-id: https://svn.dealii.org/trunk@17302 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/dofs/dof_constraints.h b/deal.II/deal.II/include/dofs/dof_constraints.h index d04118ec41..f3be1b5ec8 100644 --- a/deal.II/deal.II/include/dofs/dof_constraints.h +++ b/deal.II/deal.II/include/dofs/dof_constraints.h @@ -17,6 +17,7 @@ #include #include #include +#include #include #include @@ -24,6 +25,7 @@ DEAL_II_NAMESPACE_OPEN +template class Table; template class Vector; template class FullMatrix; class SparsityPattern; @@ -991,23 +993,52 @@ class ConstraintMatrix : public Subscriptor * second input argument is not necessary * here. * - * The last argument to this function, - * keep_constrained_entries determines - * whether the function shall allocate - * entries in the sparsity pattern at all - * for entries that will later be set to - * zero upon condensation of the - * matrix. These entries are necessary if - * the matrix is built unconstrained, and - * only later condensed. They are not - * necessary if the matrix is built using - * the distribute_local_to_global() + * The third argument to this + * function, + * keep_constrained_entries + * determines whether the + * function shall allocate + * entries in the sparsity + * pattern at all for entries + * that will later be set to zero + * upon condensation of the + * matrix. These entries are + * necessary if the matrix is + * built unconstrained, and only + * later condensed. They are not + * necessary if the matrix is + * built using the + * distribute_local_to_global() * function of this class which - * distributes entries right away when - * copying a local matrix into a global - * object. The default of this argument - * is true, meaning to allocate the few - * entries that may later be set to zero. + * distributes entries right away + * when copying a local matrix + * into a global object. The + * default of this argument is + * true, meaning to allocate the + * few entries that may later be + * set to zero. + * + * By default, the function adds + * entries for all pairs of + * indices given in the first + * argument to the sparsity + * pattern (unless + * keep_constrained_entries is + * false). However, sometimes one + * would like to only add a + * subset of all of these + * pairs. In that case, the last + * argument can be used which + * specifies a boolean mask which + * of the pairs of indices should + * be considered. If the mask is + * false for a pair of indices, + * then no entry will be added to + * the sparsity pattern for this + * pair, irrespective of whether + * one or both of the indices + * correspond to constrained + * degrees of freedom. * * This function is not typically called * from user code, but is used in the @@ -1019,7 +1050,8 @@ class ConstraintMatrix : public Subscriptor void add_entries_local_to_global (const std::vector &local_dof_indices, SparsityType &sparsity_pattern, - const bool keep_constrained_entries = true) const; + const bool keep_constrained_entries = true, + const Table<2,bool> &dof_mask = default_empty_table) const; /** * Delete hanging nodes in a @@ -1297,6 +1329,14 @@ class ConstraintMatrix : public Subscriptor * weight. */ static bool check_zero_weight (const std::pair &p); + + //public: + /** + * Dummy table that serves as + * default argument for function + * add_entries_local_to_global(). + */ + static const Table<2,bool> default_empty_table; }; diff --git a/deal.II/deal.II/include/dofs/dof_constraints.templates.h b/deal.II/deal.II/include/dofs/dof_constraints.templates.h index 3a89ede482..c9b41d4762 100644 --- a/deal.II/deal.II/include/dofs/dof_constraints.templates.h +++ b/deal.II/deal.II/include/dofs/dof_constraints.templates.h @@ -824,18 +824,50 @@ void ConstraintMatrix:: add_entries_local_to_global (const std::vector &local_dof_indices, SparsityType &sparsity_pattern, - const bool keep_constrained_entries) const + const bool keep_constrained_entries, + const Table<2,bool> &dof_mask) const { // similar to the function for distributing // matrix entries; see there for comments. const unsigned int n_local_dofs = local_dof_indices.size(); - + bool dof_mask_is_active = false; + if (dof_mask.n_rows() == n_local_dofs) + { + dof_mask_is_active = true; + Assert (dof_mask.n_cols() == n_local_dofs, + ExcDimensionMismatch(dof_mask.n_cols(), n_local_dofs)); + } + if (lines.size() == 0) { + bool add_these_indices; for (unsigned int i=0; i &local_dof_indices, } + bool add_these_indices; for (unsigned int i=0; ientries.size(); ++q) - sparsity_pattern.add (position_i->entries[q].first, + + + if ((is_constrained_i == false) && + (is_constrained_j == false) && + (keep_constrained_entries == false)) + { + sparsity_pattern.add (local_dof_indices[i], local_dof_indices[j]); - } - else if ((is_constrained_i == false) && - (is_constrained_j == true)) - { - for (unsigned int q=0; qentries.size(); ++q) - sparsity_pattern.add (local_dof_indices[i], - position_j->entries[q].first); - } - else if ((is_constrained_i == true) && - (is_constrained_j == true)) - { - for (unsigned int p=0; pentries.size(); ++p) - for (unsigned int q=0; qentries.size(); ++q) - sparsity_pattern.add (position_i->entries[p].first, + } + else if ((is_constrained_i == true) && + (is_constrained_j == false)) + { + for (unsigned int q=0; qentries.size(); ++q) + sparsity_pattern.add (position_i->entries[q].first, + local_dof_indices[j]); + } + else if ((is_constrained_i == false) && + (is_constrained_j == true)) + { + for (unsigned int q=0; qentries.size(); ++q) + sparsity_pattern.add (local_dof_indices[i], position_j->entries[q].first); - - if (i == j) - sparsity_pattern.add (local_dof_indices[i], - local_dof_indices[i]); - } - else + } + else if ((is_constrained_i == true) && + (is_constrained_j == true)) + { + for (unsigned int p=0; pentries.size(); ++p) + for (unsigned int q=0; qentries.size(); ++q) + sparsity_pattern.add (position_i->entries[p].first, + position_j->entries[q].first); + + if (i == j) + sparsity_pattern.add (local_dof_indices[i], + local_dof_indices[i]); + } + else // the only case that can // happen here is one that we // have already taken care of - Assert ((is_constrained_i == false) && - (is_constrained_j == false) && - (keep_constrained_entries == true), - ExcInternalError()); - } + Assert ((is_constrained_i == false) && + (is_constrained_j == false) && + (keep_constrained_entries == true), + ExcInternalError()); + } + } } } } diff --git a/deal.II/deal.II/include/dofs/dof_tools.h b/deal.II/deal.II/include/dofs/dof_tools.h index e2a3962092..3e9a12263b 100644 --- a/deal.II/deal.II/include/dofs/dof_tools.h +++ b/deal.II/deal.II/include/dofs/dof_tools.h @@ -168,7 +168,7 @@ template class Mapping; * * * @ingroup dofs - * @author Wolfgang Bangerth, Guido Kanschat and others, 1998 - 2005 + * @author Wolfgang Bangerth, Guido Kanschat and others, 1998 - 2008 */ class DoFTools { @@ -375,6 +375,7 @@ class DoFTools * BlockSparsityPattern, * BlockCompressedSparsityPattern, * BlockCompressedSetSparsityPattern, + * BlockCompressedSimpleSparsityPattern, * or any other class that * satisfies similar * requirements. It is assumed @@ -403,18 +404,38 @@ class DoFTools * taken into account at the time of * creating the sparsity pattern. For * this, pass the ConstraintMatrix object - * as the last argument to the current + * as the third argument to the current * function. No call to * ConstraintMatrix::condense() is then * necessary. This process is explained * in @ref step_27 "step-27". + * + * In case the constraints are + * already taken care of in this + * function, it is possible to + * neglect off-diagonal entries + * in the sparsity pattern. When + * using + * ConstraintMatrix::distribute_local_to_global + * during assembling, no entries + * will ever be written into + * these matrix position, so that + * one can save some computing + * time in matrix-vector products + * by not even creating these + * elements. In that case, the + * variable + * keep_constrained_dofs + * needs to be set to + * false. */ template static void - make_sparsity_pattern (const DH &dof, - SparsityPattern &sparsity_pattern, - const ConstraintMatrix &constraints = ConstraintMatrix()); + make_sparsity_pattern (const DH &dof, + SparsityPattern &sparsity_pattern, + const ConstraintMatrix &constraints = ConstraintMatrix(), + const bool keep_constrained_dofs = true); /** * Locate non-zero entries for @@ -499,13 +520,60 @@ class DoFTools * * Not implemented for * hp::DoFHandler. + * + * As mentioned before, the + * creation of the sparsity + * pattern is a purely local + * process and the sparsity + * pattern does not provide for + * entries introduced by the + * elimination of hanging + * nodes. They have to be taken + * care of by a call to + * ConstraintMatrix::condense() + * afterwards. + * + * Alternatively, the constraints + * on degrees of freedom can + * already be taken into account + * at the time of creating the + * sparsity pattern. For this, + * pass the ConstraintMatrix + * object as the third argument + * to the current function. No + * call to + * ConstraintMatrix::condense() + * is then necessary. This + * process is explained in @ref + * step_27 "step-27". + * + * In case the constraints are + * already taken care of in this + * function, it is possible to + * neglect off-diagonal entries + * in the sparsity pattern. When + * using + * ConstraintMatrix::distribute_local_to_global + * during assembling, no entries + * will ever be written into + * these matrix position, so that + * one can save some computing + * time in matrix-vector products + * by not even creating these + * elements. In that case, the + * variable + * keep_constrained_dofs + * needs to be set to + * false. */ template static void make_sparsity_pattern (const DH &dof, const Table<2, Coupling> &coupling, - SparsityPattern& sparsity_pattern); + SparsityPattern &sparsity_pattern, + const ConstraintMatrix &constraints = ConstraintMatrix(), + const bool keep_constrained_dofs = true); /** * @deprecated This is the old diff --git a/deal.II/deal.II/source/dofs/dof_constraints.cc b/deal.II/deal.II/source/dofs/dof_constraints.cc index 348d11b051..5dbd6097a9 100644 --- a/deal.II/deal.II/source/dofs/dof_constraints.cc +++ b/deal.II/deal.II/source/dofs/dof_constraints.cc @@ -52,6 +52,12 @@ DEAL_II_NAMESPACE_OPEN + + // Static member variable +const Table<2,bool> ConstraintMatrix::default_empty_table = Table<2,bool>(); + + + bool ConstraintMatrix::check_zero_weight (const std::pair &p) { @@ -1802,37 +1808,45 @@ MATRIX_FUNCTIONS(TrilinosWrappers::BlockSparseMatrix); template void ConstraintMatrix:: add_entries_local_to_global (const std::vector &, SparsityPattern &, - const bool) const; + const bool, + const Table<2,bool> &) const; template void ConstraintMatrix:: add_entries_local_to_global (const std::vector &, CompressedSparsityPattern &, - const bool) const; + const bool, + const Table<2,bool> &) const; template void ConstraintMatrix:: add_entries_local_to_global (const std::vector &, CompressedSetSparsityPattern &, - const bool) const; + const bool, + const Table<2,bool> &) const; template void ConstraintMatrix:: add_entries_local_to_global (const std::vector &, CompressedSimpleSparsityPattern &, - const bool) const; + const bool, + const Table<2,bool> &) const; template void ConstraintMatrix:: add_entries_local_to_global (const std::vector &, BlockSparsityPattern &, - const bool) const; + const bool, + const Table<2,bool> &) const; template void ConstraintMatrix:: add_entries_local_to_global (const std::vector &, BlockCompressedSparsityPattern &, - const bool) const; + const bool, + const Table<2,bool> &) const; template void ConstraintMatrix:: add_entries_local_to_global (const std::vector &, BlockCompressedSetSparsityPattern &, - const bool) const; + const bool, + const Table<2,bool> &) const; template void ConstraintMatrix:: add_entries_local_to_global (const std::vector &, BlockCompressedSimpleSparsityPattern &, - const bool) const; + const bool, + const Table<2,bool> &) const; DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index 980f15b1c9..710cf47970 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -48,9 +48,10 @@ DEAL_II_NAMESPACE_OPEN template void -DoFTools::make_sparsity_pattern (const DH &dof, - SparsityPattern &sparsity, - const ConstraintMatrix &constraints) +DoFTools::make_sparsity_pattern (const DH &dof, + SparsityPattern &sparsity, + const ConstraintMatrix &constraints, + const bool keep_constrained_dofs) { const unsigned int n_dofs = dof.n_dofs(); @@ -74,7 +75,8 @@ DoFTools::make_sparsity_pattern (const DH &dof, // given, then the following call acts // as if simply no constraints existed constraints.add_entries_local_to_global (dofs_on_this_cell, - sparsity); + sparsity, + keep_constrained_dofs); } } @@ -85,7 +87,9 @@ void DoFTools::make_sparsity_pattern ( const DH &dof, const Table<2,Coupling> &couplings, - SparsityPattern& sparsity) + SparsityPattern &sparsity, + const ConstraintMatrix &constraints, + const bool keep_constrained_dofs) { const unsigned int n_dofs = dof.n_dofs(); @@ -107,13 +111,12 @@ DoFTools::make_sparsity_pattern ( // respect to non-primitive shape // functions, which takes some additional // thought - std::vector > > dof_mask(fe_collection.size()); + std::vector > dof_mask(fe_collection.size()); for (unsigned int f=0; f(dofs_per_cell, false)); + dof_mask[f].reinit (dofs_per_cell, dofs_per_cell); for (unsigned int i=0; iactive_fe_index(); - + const unsigned int dofs_per_cell = fe_collection[fe_index].dofs_per_cell; - + dofs_on_this_cell.resize (dofs_per_cell); cell->get_dof_indices (dofs_on_this_cell); - - // make sparsity pattern for this cell - for (unsigned int i=0; i, SparsityPattern> (const DoFHandler &dof, SparsityPattern &sparsity, - const ConstraintMatrix &); + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, CompressedSparsityPattern> (const DoFHandler &dof, CompressedSparsityPattern &sparsity, - const ConstraintMatrix &); + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, CompressedSetSparsityPattern> (const DoFHandler &dof, CompressedSetSparsityPattern &sparsity, - const ConstraintMatrix &); + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, CompressedSimpleSparsityPattern> (const DoFHandler &dof, CompressedSimpleSparsityPattern &sparsity, - const ConstraintMatrix &); + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, BlockSparsityPattern> (const DoFHandler &dof, BlockSparsityPattern &sparsity, - const ConstraintMatrix &); + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, BlockCompressedSparsityPattern> (const DoFHandler &dof, BlockCompressedSparsityPattern &sparsity, - const ConstraintMatrix &); + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, BlockCompressedSetSparsityPattern> (const DoFHandler &dof, BlockCompressedSetSparsityPattern &sparsity, - const ConstraintMatrix &); + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, BlockCompressedSimpleSparsityPattern> (const DoFHandler &dof, BlockCompressedSimpleSparsityPattern &sparsity, - const ConstraintMatrix &); + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, SparsityPattern> (const hp::DoFHandler &dof, SparsityPattern &sparsity, - const ConstraintMatrix &); + const ConstraintMatrix &, + const bool); template void @@ -5170,21 +5185,24 @@ DoFTools::make_sparsity_pattern, CompressedSparsityPattern> (const hp::DoFHandler &dof, CompressedSparsityPattern &sparsity, - const ConstraintMatrix &); + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, CompressedSetSparsityPattern> (const hp::DoFHandler &dof, CompressedSetSparsityPattern &sparsity, - const ConstraintMatrix &); + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, CompressedSimpleSparsityPattern> (const hp::DoFHandler &dof, CompressedSimpleSparsityPattern &sparsity, - const ConstraintMatrix &); + const ConstraintMatrix &, + const bool); template void @@ -5192,25 +5210,29 @@ DoFTools::make_sparsity_pattern, BlockSparsityPattern> (const hp::DoFHandler &dof, BlockSparsityPattern &sparsity, - const ConstraintMatrix &); + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, BlockCompressedSparsityPattern> (const hp::DoFHandler &dof, BlockCompressedSparsityPattern &sparsity, - const ConstraintMatrix &); + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, BlockCompressedSetSparsityPattern> (const hp::DoFHandler &dof, BlockCompressedSetSparsityPattern &sparsity, - const ConstraintMatrix &); + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, BlockCompressedSimpleSparsityPattern> (const hp::DoFHandler &dof, BlockCompressedSimpleSparsityPattern &sparsity, - const ConstraintMatrix &); + const ConstraintMatrix &, + const bool); template void @@ -5218,98 +5240,130 @@ DoFTools::make_sparsity_pattern, SparsityPattern> (const DoFHandler&, const Table<2,Coupling>&, - SparsityPattern&); + SparsityPattern&, + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, CompressedSparsityPattern> (const DoFHandler&, const Table<2,Coupling>&, - CompressedSparsityPattern&); + CompressedSparsityPattern&, + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, CompressedSetSparsityPattern> (const DoFHandler&, const Table<2,Coupling>&, - CompressedSetSparsityPattern&); + CompressedSetSparsityPattern&, + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, CompressedSimpleSparsityPattern> (const DoFHandler&, const Table<2,Coupling>&, - CompressedSimpleSparsityPattern&); + CompressedSimpleSparsityPattern&, + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, BlockSparsityPattern> (const DoFHandler&, const Table<2,Coupling>&, - BlockSparsityPattern&); + BlockSparsityPattern&, + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, BlockCompressedSparsityPattern> (const DoFHandler&, const Table<2,Coupling>&, - BlockCompressedSparsityPattern&); + BlockCompressedSparsityPattern&, + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, BlockCompressedSetSparsityPattern> (const DoFHandler&, const Table<2,Coupling>&, - BlockCompressedSetSparsityPattern&); + BlockCompressedSetSparsityPattern&, + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, BlockCompressedSimpleSparsityPattern> (const DoFHandler&, const Table<2,Coupling>&, - BlockCompressedSimpleSparsityPattern&); + BlockCompressedSimpleSparsityPattern&, + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, SparsityPattern> (const hp::DoFHandler&, const Table<2,Coupling>&, - SparsityPattern&); + SparsityPattern&, + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, CompressedSparsityPattern> (const hp::DoFHandler&, const Table<2,Coupling>&, - CompressedSparsityPattern&); + CompressedSparsityPattern&, + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, CompressedSetSparsityPattern> (const hp::DoFHandler&, const Table<2,Coupling>&, - CompressedSetSparsityPattern&); + CompressedSetSparsityPattern&, + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, CompressedSimpleSparsityPattern> (const hp::DoFHandler&, const Table<2,Coupling>&, - CompressedSimpleSparsityPattern&); + CompressedSimpleSparsityPattern&, + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, BlockSparsityPattern> (const hp::DoFHandler&, const Table<2,Coupling>&, - BlockSparsityPattern&); + BlockSparsityPattern&, + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, BlockCompressedSparsityPattern> (const hp::DoFHandler&, const Table<2,Coupling>&, - BlockCompressedSparsityPattern&); + BlockCompressedSparsityPattern&, + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, BlockCompressedSetSparsityPattern> (const hp::DoFHandler&, const Table<2,Coupling>&, - BlockCompressedSetSparsityPattern&); + BlockCompressedSetSparsityPattern&, + const ConstraintMatrix &, + const bool); template void DoFTools::make_sparsity_pattern, BlockCompressedSimpleSparsityPattern> (const hp::DoFHandler&, const Table<2,Coupling>&, - BlockCompressedSimpleSparsityPattern&); + BlockCompressedSimpleSparsityPattern&, + const ConstraintMatrix &, + const bool); template void diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index 94005cd216..6e67d2d137 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -353,6 +353,13 @@ inconvenience this causes.

deal.II

    +
  1. +

    + New: When calling function DoFTools::make_sparsity_pattern with a ConstraintMatrix, it is now possible to set a bool argument keep_constrained_dofs. When this flag is set to false, constrained rows and columns will not be part of the sparsity pattern, which increases the performance of matrix operations and decrease memory consumption in case there are many constraints. +
    + (Martin Kronbichler 2008/10/21) +

    +
  2. New: There is now a second DoFTools::count_dofs_with_subdomain_association function that diff --git a/deal.II/examples/step-31/step-31.cc b/deal.II/examples/step-31/step-31.cc index f2e8620555..2cdf2f6e24 100644 --- a/deal.II/examples/step-31/step-31.cc +++ b/deal.II/examples/step-31/step-31.cc @@ -1,5 +1,5 @@ /* $Id$ */ -/* Author: Martin Kronbichler, University of Uppsala, +/* Author: Martin Kronbichler, Uppsala University, Wolfgang Bangerth, Texas A&M University 2007, 2008 */ /* $Id$ */ @@ -55,9 +55,10 @@ // interfaces to the matrix and vector // classes based on Trilinos as well as // Trilinos preconditioners: -#include #include #include +#include +#include #include // Finally, here are two C++ headers that @@ -1053,13 +1054,36 @@ void BoussinesqFlowProblem::setup_dofs () // components at the boundary // again. // - // Then, constraints are applied to - // the temporary sparsity patterns, - // which are finally copied into an - // object of type SparsityPattern - // and used to initialize the - // nonzero pattern of the Trilinos - // matrix objects we use. + // When generating the sparsity + // pattern, we directly apply the + // constraints from hanging nodes + // and no-flux boundary + // conditions. This approach was + // already used in step-27, but is + // different from the one in early + // tutorial programs. The reason + // for doing so is that later + // during assembly we are going to + // distribute the constraints + // immediately when transferring + // local to global + // dofs. Consequently, there will + // be no data written at positions + // of constrained degrees of + // freedom, so we can let the + // DoFTools::make_sparsity_pattern + // function omit these entries by + // setting the last boolean flag to + // false. Once the + // sparsity pattern is ready, we + // can use it to initialize the + // Trilinos matrices. Note that the + // Trilinos matrices store the + // sparsity pattern internally, so + // there is no need to keep the + // sparsity pattern around after + // the initialization of the + // matrix. stokes_block_sizes.resize (2); stokes_block_sizes[0] = n_u; stokes_block_sizes[1] = n_p; @@ -1084,8 +1108,8 @@ void BoussinesqFlowProblem::setup_dofs () else coupling[c][d] = DoFTools::none; - DoFTools::make_sparsity_pattern (stokes_dof_handler, coupling, csp); - stokes_constraints.condense (csp); + DoFTools::make_sparsity_pattern (stokes_dof_handler, coupling, csp, + stokes_constraints, false); BlockSparsityPattern stokes_sparsity_pattern; stokes_sparsity_pattern.copy_from (csp); @@ -1116,8 +1140,8 @@ void BoussinesqFlowProblem::setup_dofs () else coupling[c][d] = DoFTools::none; - DoFTools::make_sparsity_pattern (stokes_dof_handler, coupling, csp); - stokes_constraints.condense (csp); + DoFTools::make_sparsity_pattern (stokes_dof_handler, coupling, csp, + stokes_constraints, false); BlockSparsityPattern stokes_preconditioner_sparsity_pattern; stokes_preconditioner_sparsity_pattern.copy_from (csp); @@ -1126,14 +1150,20 @@ void BoussinesqFlowProblem::setup_dofs () stokes_preconditioner_matrix.collect_sizes(); } + // The generation of the + // temperature matrix follows the + // generation of the Stokes matrix + // – except that it is much + // easier since we do not need to + // take care of any blocks. { temperature_mass_matrix.clear (); temperature_stiffness_matrix.clear (); temperature_matrix.clear (); CompressedSetSparsityPattern csp (n_T, n_T); - DoFTools::make_sparsity_pattern (temperature_dof_handler, csp); - temperature_constraints.condense (csp); + DoFTools::make_sparsity_pattern (temperature_dof_handler, csp, + temperature_constraints, false); SparsityPattern temperature_sparsity_pattern; temperature_sparsity_pattern.copy_from (csp); @@ -1440,36 +1470,33 @@ BoussinesqFlowProblem::build_stokes_preconditioner () // work by doing the full assembly // only when it is needed. // - // Regarding the technical details - // of implementation, not much has + // Regarding the technical details of + // implementation, not much has // changed from step-22. We reset - // matrix and vector, create - // a quadrature formula on the - // cells and one on cell faces - // (for implementing Neumann - // boundary conditions). Then, - // we create a respective - // FEValues object for both the - // cell and the face integration. - // For the the update flags of - // the first, we perform the + // matrix and vector, create a + // quadrature formula on the cells + // and one on cell faces (for + // implementing Neumann boundary + // conditions). Then, we create a + // respective FEValues object for + // both the cell and the face + // integration. For the the update + // flags of the first, we perform the // calculations of basis function - // derivatives only in - // case of a full assembly, since - // they are not needed otherwise, - // which makes the call of - // the FEValues::reinit function - // further down in the program - // more efficient. + // derivatives only in case of a full + // assembly, since they are not + // needed otherwise, which makes the + // call of the FEValues::reinit + // function further down in the + // program more efficient. // - // The declarations proceed - // with some shortcuts for - // array sizes, the creation of - // the local matrix and right - // hand side as well as the - // vector for the indices of - // the local dofs compared to - // the global system. + // The declarations proceed with some + // shortcuts for array sizes, the + // creation of the local matrix and + // right hand side as well as the + // vector for the indices of the + // local dofs compared to the global + // system. template void BoussinesqFlowProblem::assemble_stokes_system () { @@ -1520,20 +1547,20 @@ void BoussinesqFlowProblem::assemble_stokes_system () // term in the momentum equation. // // The set of vectors we create - // next hold the evaluations of - // the basis functions that will - // be used for creating the + // next hold the evaluations of the + // basis functions that will be + // used for creating the // matrices. This gives faster // access to that data, which - // increases the performance - // of the assembly. See step-22 - // for details. + // increases the performance of the + // assembly. See step-22 for + // details. // - // The last two declarations - // are used to extract the - // individual blocks (velocity, - // pressure, temperature) from - // the total FE system. + // The last two declarations are + // used to extract the individual + // blocks (velocity, pressure, + // temperature) from the total FE + // system. std::vector boundary_values (n_face_q_points); std::vector old_temperature_values(n_q_points); diff --git a/deal.II/examples/step-32/step-32.cc b/deal.II/examples/step-32/step-32.cc index 668a1df565..b0dd2d33ec 100644 --- a/deal.II/examples/step-32/step-32.cc +++ b/deal.II/examples/step-32/step-32.cc @@ -1,5 +1,5 @@ /* $Id$ */ -/* Author: Martin Kronbichler, University of Uppsala, +/* Author: Martin Kronbichler, Uppsala University, Wolfgang Bangerth, Texas A&M University 2007, 2008 */ /* */ /* Copyright (C) 2008 by the deal.II authors */ @@ -659,8 +659,8 @@ void BoussinesqFlowProblem::setup_dofs () else coupling[c][d] = DoFTools::none; - DoFTools::make_sparsity_pattern (stokes_dof_handler, coupling, csp); - stokes_constraints.condense (csp); + DoFTools::make_sparsity_pattern (stokes_dof_handler, coupling, csp, + stokes_constraints, false); stokes_sparsity_pattern.copy_from (csp); @@ -691,7 +691,8 @@ void BoussinesqFlowProblem::setup_dofs () else coupling[c][d] = DoFTools::none; - DoFTools::make_sparsity_pattern (stokes_dof_handler, coupling, csp); + DoFTools::make_sparsity_pattern (stokes_dof_handler, coupling, csp, + stokes_constraints, false); stokes_constraints.condense (csp); stokes_preconditioner_sparsity_pattern.copy_from (csp); @@ -716,7 +717,8 @@ void BoussinesqFlowProblem::setup_dofs () SparsityPattern temperature_sparsity_pattern; CompressedSetSparsityPattern csp (n_T, n_T); - DoFTools::make_sparsity_pattern (temperature_dof_handler, csp); + DoFTools::make_sparsity_pattern (temperature_dof_handler, csp, + temperature_constraints, false); temperature_constraints.condense (csp); temperature_sparsity_pattern.copy_from (csp);