* All member functions are static, so there is no need to create an
* object of class #MGDoFTools#.
*
- * @author Wolfgang Bangerth and others, 1999
+ * @author Wolfgang Bangerth, Guido Kanschat, 1999, 2000
*/
class MGDoFTools
{
SparsityPattern &sparsity,
const unsigned int level);
+ /**
+ * Create sparsity pattern for
+ * the fluxes at refinement
+ * edges.
+ * @see{make_flux_sparsity_pattern}
+ * @see{DoFTools}
+ */
+ template <int dim>
+ static void
+ make_flux_sparsity_pattern_edge (const MGDoFHandler<dim> &dof_handler,
+ SparsityPattern &sparsity,
+ const unsigned int level);
+
/**
* This function does the same as
* the other with the same name,
for (unsigned int i=0; i<total_dofs; ++i)
{
for (unsigned int j=0; j<total_dofs; ++j)
- if (flux_dof_mask[i][j])
- {
+ {
+ if (flux_dof_mask[i][j])
sparsity.add (dofs_on_this_cell[i],
dofs_on_other_cell[j]);
+ if (flux_dof_mask[j][i])
sparsity.add (dofs_on_other_cell[i],
dofs_on_this_cell[j]);
- }
+ }
}
sub_neighbor->face(neighbor_face)->set_user_flag ();
}
for (unsigned int i=0; i<total_dofs; ++i)
{
for (unsigned int j=0; j<total_dofs; ++j)
- if (flux_dof_mask[i][j])
- {
+ {
+ if (flux_dof_mask[i][j])
sparsity.add (dofs_on_this_cell[i],
dofs_on_other_cell[j]);
+ if (flux_dof_mask[j][i])
sparsity.add (dofs_on_other_cell[i],
dofs_on_this_cell[j]);
}
}
}
+
+
template<int dim>
void
MGDoFTools::make_flux_sparsity_pattern (const MGDoFHandler<dim> &dof,
+template<int dim>
+void
+MGDoFTools::make_flux_sparsity_pattern_edge (const MGDoFHandler<dim> &dof,
+ SparsityPattern &sparsity,
+ const unsigned int level)
+{
+ Assert ((level>=1) && (level<dof.get_tria().n_levels()),
+ ExcIndexRange(level, 1, dof.get_tria().n_levels()));
+
+ const unsigned int fine_dofs = dof.n_dofs(level);
+ const unsigned int coarse_dofs = dof.n_dofs(level-1);
+
+ Assert (sparsity.n_rows() == coarse_dofs,
+ ExcDimensionMismatch (sparsity.n_rows(), coarse_dofs));
+ Assert (sparsity.n_cols() == fine_dofs,
+ ExcDimensionMismatch (sparsity.n_cols(), fine_dofs));
+
+ const unsigned int dofs_per_cell = dof.get_fe().dofs_per_cell;
+ vector<unsigned int> dofs_on_this_cell(dofs_per_cell);
+ vector<unsigned int> dofs_on_other_cell(dofs_per_cell);
+ MGDoFHandler<dim>::cell_iterator cell = dof.begin(level),
+ endc = dof.end(level);
+ for (; cell!=endc; ++cell)
+ {
+ cell->get_mg_dof_indices (dofs_on_this_cell);
+ // Loop over all interior neighbors
+ for (unsigned int face = 0;
+ face < GeometryInfo<dim>::faces_per_cell;
+ ++face)
+ {
+ if ( (! cell->at_boundary(face)) &&
+ (static_cast<unsigned int>(cell->neighbor_level(face)) != level) )
+ {
+ MGDoFHandler<dim>::cell_iterator neighbor = cell->neighbor(face);
+ neighbor->get_mg_dof_indices (dofs_on_other_cell);
+
+ for (unsigned int i=0; i<dofs_per_cell; ++i)
+ {
+ for (unsigned int j=0; j<dofs_per_cell; ++j)
+ {
+ sparsity.add (dofs_on_this_cell[i],
+ dofs_on_other_cell[j]);
+ sparsity.add (dofs_on_this_cell[j],
+ dofs_on_other_cell[i]);
+ }
+ }
+ }
+ }
+ }
+}
+
+
+
template<int dim>
void
MGDoFTools::make_flux_sparsity_pattern (const MGDoFHandler<dim> &dof,
SparsityPattern &,
const unsigned int);
+template void MGDoFTools::make_flux_sparsity_pattern_edge (const MGDoFHandler<deal_II_dimension> &,
+ SparsityPattern &,
+ const unsigned int);
+
#if deal_II_dimension > 1
template void MGDoFTools::make_flux_sparsity_pattern (const MGDoFHandler<deal_II_dimension> &,
SparsityPattern &,