Table<2,Coupling> flux_couplings)
{
Assert(false, ExcNotImplemented());
-
- const FiniteElement<dim>& fe = dofs.get_fe();
- const unsigned int ncomp = fe.n_components;
-
- Assert (row_lengths.size() == dofs.n_dofs(),
- ExcDimensionMismatch(row_lengths.size(), dofs.n_dofs()));
-
- Assert (couplings.n_rows() == couplings.n_cols(),
- ExcDimensionMismatch(couplings.n_rows(), couplings.n_cols()));
- Assert (flux_couplings.n_rows() == flux_couplings.n_cols(),
- ExcDimensionMismatch(flux_couplings.n_rows(), flux_couplings.n_cols()));
-
- if (couplings.n_rows() > 1)
- Assert (couplings.n_rows() == ncomp,
- ExcDimensionMismatch(couplings.n_rows(), ncomp));
- if (flux_couplings.n_rows() > 1)
- Assert (flux_couplings.n_rows() == ncomp,
- ExcDimensionMismatch(flux_couplings.n_rows(), ncomp));
-
- // Function starts here by
- // resetting the counters.
- std::fill(row_lengths.begin(), row_lengths.end(), 0);
- // We need the user flags, so we
- // save them for later restoration
- std::vector<bool> old_flags;
- user_flags_triangulation.save_user_flags(old_flags);
- user_flags_triangulation.clear_user_flags();
-
- const typename DoFHandler<dim>::cell_iterator end = dofs.end();
- typename DoFHandler<dim>::active_cell_iterator cell;
- std::vector<unsigned int> cell_indices(fe.dofs_per_cell);
- std::vector<unsigned int> face_indices(fe.dofs_per_face);
-
- for (cell = dofs.begin(); cell != end; ++cell)
- {
- cell->get_dof_indices(cell_indices);
- // At this point, we have
- // counted all dofs
- // contributiong from cells
- // coupled topologically to the
- // adjacent cells, but we
- // subtracted some faces.
-
- // Now, let's go by the faces
- // and add the missing
- // contribution as well as the
- // flux contributions.
- for (unsigned int iface=0;iface<GeometryInfo<dim>::faces_per_cell;++iface)
- {
- if (cell->at_boundary(iface))
- continue;
- const bool neighbor_same_level = (cell->neighbor(iface)->level()
- == cell->level());
- const unsigned int
- nface = (neighbor_same_level)
- ? cell->neighbor_of_neighbor(iface)
- : cell->neighbor_of_coarser_neighbor(iface);
- typename DoFHandler<dim>::face_iterator face = cell->face(iface);
-
- // Flux couplings are
- // computed from both sides
- // for simplicity.
- for (unsigned int i=0;i<fe.dofs_per_face;++i)
- for (unsigned int j=0;j<fe.dofs_per_face;++j)
- if (flux_couplings(fe.system_to_component_index(i).first,
- fe.system_to_component_index(j).first) == always
- ||
- flux_couplings(fe.system_to_component_index(i).first,
- fe.system_to_component_index(j).first) == nonzero
- && fe.has_support_on_face(i, iface)
- && fe.has_support_on_face(j, nface))
- ;// ???
- }
- }
- user_flags_triangulation.load_user_flags(old_flags);
}