*/
void make_sparsity_pattern (SparseMatrixStruct &) const;
+ /**
+ * This function does mistly the same as
+ * the above one, but it is specialized for
+ * mixed finite elements and allows to
+ * specify which variables couple in which
+ * equation. For example, if wanted to solve
+ * the Stokes equations,
+ * \begin{verbatim}
+ * -\Delta \vec u + \nabla p = 0,
+ * \div u = 0
+ * \end{verbatim}
+ * in, say, two space dimensions, using
+ * nonstabilized Q2/Q2/Q1 mixed elements
+ * (i.e. using the #FESystem# class), then
+ * you don't want all degrees of freedom
+ * to couple in each equation. You rather
+ * may want to give the following pattern
+ * of couplings:
+ * \begin{verbatim}
+ * 1 0 1
+ * 0 1 1
+ * 1 1 0
+ * \end{verbatim}
+ * where "1" indicates that two variables
+ * (i.e. components of the #FESystem#)
+ * couple in the respective equation, and
+ * a "0" means no coupling, in which case
+ * it is not necessary to allocate space
+ * in the matrix structure. Obviously, the
+ * mask refers to components of the
+ * composed #FESystem#, rather than to the
+ * degrees of freedom contained in there.
+ *
+ * This function is designed to accept
+ * a mask, like the one shown above,
+ * through the #mask# parameter, which
+ * contains boolean values. It builds
+ * the matrix structure just like the
+ * previous function, but does not create
+ * elements if not specified by the mask.
+ */
+ void make_sparsity_pattern (const vector<vector<bool> > &mask,
+ SparseMatrixStruct &) const;
+
/**
* Write the sparsity structure of the
* matrix composed of the basis functions
int,
<< "The given list of new dof indices is not consecutive: "
<< "the index " << arg1 << " does not exist.");
-
+ /**
+ * Exception
+ */
+ DeclException2 (ExcInvalidMaskDimension,
+ int, int,
+ << "The dimension of the mask " << arg1 << " does not match"
+ << " the number of components in the finite element object.");
+
protected:
/**
};
+
+template <int dim>
+void
+DoFHandler<dim>::make_sparsity_pattern (const vector<vector<bool> > &mask,
+ SparseMatrixStruct &sparsity) const
+{
+ const unsigned int total_dofs = selected_fe->total_dofs;
+
+ Assert (selected_fe != 0, ExcNoFESelected());
+ Assert (sparsity.n_rows() == n_dofs(),
+ ExcDifferentDimensions (sparsity.n_rows(), n_dofs()));
+ Assert (sparsity.n_cols() == n_dofs(),
+ ExcDifferentDimensions (sparsity.n_cols(), n_dofs()));
+ Assert (mask.size() == selected_fe->n_components,
+ ExcInvalidMaskDimension(mask.size(), selected_fe->n_components));
+ for (unsigned int i=0; i<mask.size(); ++i)
+ Assert (mask[i].size() == selected_fe->n_components,
+ ExcInvalidMaskDimension(mask[i].size(), selected_fe->n_components));
+
+ // first build a mask for each dof,
+ // not like the one given which represents
+ // components
+ vector<vector<bool> > dof_mask(total_dofs,
+ vector<bool>(total_dofs, false));
+ for (unsigned int i=0; i<total_dofs; ++i)
+ for (unsigned int j=0; j<total_dofs; ++j)
+ dof_mask[i][j] = mask
+ [selected_fe->system_to_component_index(i).first]
+ [selected_fe->system_to_component_index(j).first];
+
+
+ vector<int> dofs_on_this_cell(total_dofs);
+ active_cell_iterator cell = begin_active(),
+ endc = end();
+ for (; cell!=endc; ++cell)
+ {
+ cell->get_dof_indices (dofs_on_this_cell);
+ // make sparsity pattern for this cell
+ for (unsigned int i=0; i<total_dofs; ++i)
+ for (unsigned int j=0; j<total_dofs; ++j)
+ if (dof_mask[i][j] == true)
+ sparsity.add (dofs_on_this_cell[i],
+ dofs_on_this_cell[j]);
+ };
+};
+
+
+
+
#if deal_II_dimension == 1
template <>