* 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
* 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
+ * <tt>local_matrix</tt> 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 <typename VectorType>
void
distribute_local_to_global (const Vector<double> &local_vector,
const std::vector<unsigned int> &local_dof_indices,
- VectorType &global_vector) const;
+ VectorType &global_vector,
+ const FullMatrix<double> &local_matrix = FullMatrix<double>()) const;
/**
* This function takes a matrix of
ConstraintMatrix::
distribute_local_to_global (const Vector<double> &local_vector,
const std::vector<unsigned int> &local_dof_indices,
- VectorType &global_vector) const
+ VectorType &global_vector,
+ const FullMatrix<double> &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<n_local_dofs; ++i)
{
// let's see if we can use the bool
// vector that tells about whether a
- // certain constraint exists.
- if (constraint_line_exists.size() <= local_dof_indices[i])
- {
- global_vector(local_dof_indices[i]) += local_vector(i);
- continue;
- }
- if (constraint_line_exists[local_dof_indices[i]] == false)
+ // certain constraint exists. then, we
+ // simply copy over the data.
+ if (is_constrained (local_dof_indices[i]) == false)
{
global_vector(local_dof_indices[i]) += local_vector(i);
continue;
}
- // first figure out whether this
- // dof is constrained
+
+ // if the dof is constrained, we
+ // have to find it in the list of
+ // all constraints.
ConstraintLine index_comparison;
index_comparison.line = local_dof_indices[i];
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<n_local_dofs; ++j)
+ if (is_constrained (local_dof_indices[j]) == false)
+ {
+ global_vector(local_dof_indices[j]) -= val * local_matrix(j,i);
+ }
+ else
+ {
+ const double matrix_entry = local_matrix(j,i);
+ if (matrix_entry == 0)
+ continue;
+
+ ConstraintLine index_comparison_j;
+ index_comparison_j.line = local_dof_indices[j];
+
+ const std::vector<ConstraintLine>::const_iterator
+ position_j = std::lower_bound (lines.begin(),
+ lines.end(),
+ index_comparison_j);
+ for (unsigned int q=0; q<position_j->entries.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; j<position->entries.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;
+ }
}
}