From: Martin Kronbichler Date: Wed, 29 Apr 2009 16:31:59 +0000 (+0000) Subject: Implement a function that can distribute vectors with inhomogeneous constraints,... X-Git-Tag: v8.0.0~7725 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=16d4c4279777e7d47cd828751d925a62fefb3e35;p=dealii.git Implement a function that can distribute vectors with inhomogeneous constraints, but operating only on vectors. git-svn-id: https://svn.dealii.org/trunk@18790 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/lac/include/lac/constraint_matrix.h b/deal.II/lac/include/lac/constraint_matrix.h index b26c2b2616..0d0027a4e2 100644 --- a/deal.II/lac/include/lac/constraint_matrix.h +++ b/deal.II/lac/include/lac/constraint_matrix.h @@ -928,12 +928,7 @@ class ConstraintMatrix : public Subscriptor * local_dof_indices have the same * number of elements, this function is * happy with whatever it is - * given. Note that this function will - * apply all constraints as if they - * were homogeneous. For correctly - * setting inhomogeneous constraints, - * use the function with both matrix - * and vector arguments. + * given. * * In contrast to the similar function * in the DoFAccessor class, this @@ -956,32 +951,47 @@ class ConstraintMatrix : public Subscriptor * vectors and matrices are fully * assembled. * - * In order to do its work properly, - * this function has to know which - * degrees of freedom are fixed, for - * example boundary values. For this, - * the third argument is a map between - * the numbers of the DoFs that are - * fixed and the values they are fixed - * to. One can pass an empty map in for - * this argument, but note that you - * will then have to fix these nodes - * later on again, for example by using - * MatrixTools::apply_boundary_values - * to the resulting matrix. However, - * since the present function was - * written for the express purpose of - * not having to use tools that later - * modify the matrix, it is advisable - * to have the list of fixed nodes - * available when calling the present - * function. + * Note that, unless an optional + * FullMatrix object is provided, this + * function will apply all constraints + * as if they were homogeneous and + * throw an exception in case it + * encounters inhomogeneous + * constraints. For correctly setting + * inhomogeneous constraints, you + * should provide an additional matrix + * argument or use one of the functions + * with both matrix and vector + * arguments. + * + * The optional argument + * local_matrix is intended to + * be used in case one wants to apply + * inhomogeneous constraints on the + * vector only. Such a situation could + * be where one wants to assemble of a + * right hand side vector on a problem + * with inhomogeneous constraints, but + * the global matrix has been assembled + * previously. A typical example of + * this is a time stepping algorithm + * where the stiffness matrix is + * assembled once, and the right hand + * side updated every time step. Note + * that, however, the entries in the + * columns of the local matrix have to + * be exactly the same as those that + * have been written into the global + * matrix. Otherwise, this function + * will not be able to correctly handle + * inhomogeneities. */ template void distribute_local_to_global (const Vector &local_vector, const std::vector &local_dof_indices, - VectorType &global_vector) const; + VectorType &global_vector, + const FullMatrix &local_matrix = FullMatrix()) const; /** * This function takes a matrix of diff --git a/deal.II/lac/include/lac/constraint_matrix.templates.h b/deal.II/lac/include/lac/constraint_matrix.templates.h index e5b6084994..10e10a3b34 100644 --- a/deal.II/lac/include/lac/constraint_matrix.templates.h +++ b/deal.II/lac/include/lac/constraint_matrix.templates.h @@ -727,35 +727,39 @@ void ConstraintMatrix:: distribute_local_to_global (const Vector &local_vector, const std::vector &local_dof_indices, - VectorType &global_vector) const + VectorType &global_vector, + const FullMatrix &local_matrix) const { Assert (local_vector.size() == local_dof_indices.size(), ExcDimensionMismatch(local_vector.size(), local_dof_indices.size())); Assert (sorted == true, ExcMatrixNotClosed()); + const bool use_matrix = local_matrix.m() != 0 ? true : false; + if (local_matrix.m() != 0) + { + Assert (local_matrix.m() == local_dof_indices.size(), + ExcDimensionMismatch(local_matrix.m(), local_dof_indices.size())); + Assert (local_matrix.n() == local_dof_indices.size(), + ExcDimensionMismatch(local_matrix.n(), local_dof_indices.size())); + } const unsigned int n_local_dofs = local_vector.size(); - // have a special case where there are no - // constraints at all, since then we can be - // a lot faster Threads::ThreadMutex::ScopedLock lock(mutex); for (unsigned int i=0; i &local_vector, lines.end(), index_comparison); - // if the line is not - // constrained, then simply - // copy the data. otherwise - // distribute it, but make - // sure we don't touch the - // entries of fixed dofs - // - // there is one critical - // point: sometimes a dof - // may be both constrained - // and fixed, for example - // hanging nodes in 3d at - // the boundary. in that - // case, we don't quite - // know what to do -- - // handle the constraint or - // the fixed - // value. however, this - // isn't so hard if all the - // nodes that this node is - // constrained to are also - // fixed nodes, in which - // case we could do both - // but opt to copy the - // element. however, we - // have to check that all - // the nodes to which it is - // constrained are also - // fixed + // check whether we've really found the + // right constraint. Assert (position->line == local_dof_indices[i], ExcInternalError()); + + if (use_matrix) + { + const double val = position->inhomogeneity; + if (val != 0) + for (unsigned int j=0; j::const_iterator + position_j = std::lower_bound (lines.begin(), + lines.end(), + index_comparison_j); + for (unsigned int q=0; qentries.size(); ++q) + { + Assert (is_constrained(position_j->entries[q].first) == false, + ExcMessage ("Tried to distribute to a fixed dof.")); + global_vector(position_j->entries[j].first) + -= val * position_j->entries[q].second * matrix_entry; + } + } + } + else + // in case the constraint is + // inhomogeneous and we have no matrix + // available, this function is not + // appropriate. Throw an exception. + Assert (position->inhomogeneity == 0., + ExcMessage ("Inhomogeneous constraint cannot be condensed " + "without any matrix specified.")); + + // now distribute the constraint, + // but make sure we don't touch + // the entries of fixed dofs for (unsigned int j=0; jentries.size(); ++j) - global_vector(position->entries[j].first) - += local_vector(i) * position->entries[j].second; + { + Assert (is_constrained(position->entries[j].first) == false, + ExcMessage ("Tried to distribute to a fixed dof.")); + global_vector(position->entries[j].first) + += local_vector(i) * position->entries[j].second; + } } } diff --git a/deal.II/lac/source/constraint_matrix.cc b/deal.II/lac/source/constraint_matrix.cc index 7d165e7dc5..d5d26297f8 100644 --- a/deal.II/lac/source/constraint_matrix.cc +++ b/deal.II/lac/source/constraint_matrix.cc @@ -2004,7 +2004,8 @@ ConstraintMatrix::memory_consumption () const template void ConstraintMatrix:: \ distribute_local_to_global (const Vector &, \ const std::vector &, \ - VectorType &) const; \ + VectorType &, \ + const FullMatrix &) const; \ template void ConstraintMatrix::distribute(const VectorType &condensed,\ VectorType &uncondensed) const;\ template void ConstraintMatrix::distribute(VectorType &vec) const