50 "The Triangulation has not been classified. You need to call the "
51 "reclassify()-function before using this function.");
55 "The incoming cell does not belong to the triangulation passed to "
63 template <
typename VectorType>
67 const auto min_max_element =
68 std::minmax_element(local_levelset_values.begin(),
69 local_levelset_values.end());
71 if (*min_max_element.second < 0)
73 if (0 < *min_max_element.first)
85 template <
int dim,
typename VectorType>
107 &cell)
const override;
116 const unsigned int face_index,
134 template <
int dim,
typename VectorType>
144 template <
int dim,
typename VectorType>
153 template <
int dim,
typename VectorType>
157 const unsigned int face_index,
160 const auto cell_with_dofs = cell->as_dof_handler_iterator(*
dof_handler);
162 const unsigned int n_dofs_per_face =
164 std::vector<types::global_dof_index> dof_indices(n_dofs_per_face);
165 cell_with_dofs->face(face_index)->get_dof_indices(dof_indices);
167 local_levelset_values.
reinit(dof_indices.size());
169 for (
unsigned int i = 0; i < dof_indices.size(); i++)
170 local_levelset_values[i] =
177 template <
int dim,
typename VectorType>
182 const auto cell_with_dofs = cell->as_dof_handler_iterator(*
dof_handler);
184 return cell_with_dofs->active_fe_index();
216 &cell)
const override;
225 const unsigned int face_index,
257 element.get_unit_face_support_points()),
267 const unsigned int face_index,
274 const std::vector<Point<dim>> &points =
277 for (
unsigned int i = 0; i < points.size(); i++)
278 local_levelset_values[i] =
level_set->value(points[i]);
305 template <
typename VectorType>
307 const VectorType &level_set)
310 std::make_unique<
internal::MeshClassifierImplementation::
311 DiscreteLevelSetDescription<dim, VectorType>>(
315#ifdef DEAL_II_WITH_LAPACK
318 for (
unsigned int i = 0; i < fe_collection.
size(); i++)
323 Assert(fe_collection[i].has_face_support_points(),
325 "The elements in the FECollection of the incoming DoFHandler "
326 "must have face support points."));
341 std::make_unique<
internal::MeshClassifierImplementation::
342 AnalyticLevelSetDescription<dim>>(level_set,
364 for (
const auto &cell :
triangulation->active_cell_iterators())
365 if (!cell->is_artificial())
373 for (
unsigned int f = 1; f < GeometryInfo<dim>::faces_per_cell; ++f)
380 if (face_location != face0_location)
393 const unsigned int face_index)
405 const unsigned int n_local_dofs =
411 local_levelset_values);
421 const bool is_linear = fe_q_iso_q1 !=
nullptr ||
422 (fe_poly !=
nullptr && fe_poly->
get_degree() == 1);
428 local_levelset_values);
432 local_levelset_values);
435 local_levelset_values);
459 const unsigned int face_index)
const
484 for (
unsigned int i = 0; i < fe_collection.
size(); i++)
494 for (
unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; f++)
500 *fe_q, face_interpolation_matrix, f);
510#include "mesh_classifier.inst"
TriaIterator< TriaAccessor< dim - 1, dim, spacedim > > face(const unsigned int i) const
unsigned int active_cell_index() const
const hp::FECollection< dim, spacedim > & get_fe_collection() const
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix, const unsigned int face_no=0) const override
unsigned int get_degree() const
const unsigned int dofs_per_face
unsigned int n_components() const
const unsigned int n_components
LocationToLevelSet location_to_level_set(const typename Triangulation< dim >::cell_iterator &cell) const
std::vector< LocationToLevelSet > face_locations
LocationToLevelSet determine_face_location_to_levelset(const typename Triangulation< dim >::active_cell_iterator &cell, const unsigned int face_index)
const std::unique_ptr< internal::MeshClassifierImplementation::LevelSetDescription< dim > > level_set_description
const SmartPointer< const Triangulation< dim > > triangulation
std::vector< std::array< LAPACKFullMatrix< double >, GeometryInfo< dim >::faces_per_cell > > lagrange_to_bernstein_face
std::vector< LocationToLevelSet > cell_locations
MeshClassifier(const DoFHandler< dim > &level_set_dof_handler, const VectorType &level_set)
AnalyticLevelSetDescription(const Function< dim > &level_set, const FiniteElement< dim > &element)
unsigned int active_fe_index(const typename Triangulation< dim >::active_cell_iterator &cell) const override
FEFaceValues< dim > fe_face_values
const hp::FECollection< dim > & get_fe_collection() const override
const hp::FECollection< dim > fe_collection
const SmartPointer< const Function< dim > > level_set
void get_local_level_set_values(const typename Triangulation< dim >::active_cell_iterator &cell, const unsigned int face_index, Vector< double > &local_levelset_values) override
void get_local_level_set_values(const typename Triangulation< dim >::active_cell_iterator &cell, const unsigned int face_index, Vector< double > &local_levelset_values) override
const SmartPointer< const DoFHandler< dim > > dof_handler
unsigned int active_fe_index(const typename Triangulation< dim >::active_cell_iterator &cell) const override
const hp::FECollection< dim > & get_fe_collection() const override
const SmartPointer< const VectorType > level_set
DiscreteLevelSetDescription(const DoFHandler< dim > &dof_handler, const VectorType &level_set)
const Triangulation< dim, spacedim > & get_triangulation() const
virtual size_type size() const override
virtual void reinit(const size_type N, const bool omit_zeroing_entries=false)
unsigned int size() const
unsigned int n_components() const
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
#define DeclExceptionMsg(Exception, defaulttext)
static ::ExceptionBase & ExcNeedsLAPACK()
static ::ExceptionBase & ExcReclassifyNotCalled()
static ::ExceptionBase & ExcMessage(std::string arg1)
static ::ExceptionBase & ExcTriangulationMismatch()
#define AssertThrow(cond, exc)
TriaActiveIterator< CellAccessor< dim, spacedim > > active_cell_iterator
TriaIterator< CellAccessor< dim, spacedim > > cell_iterator
@ update_quadrature_points
Transformed quadrature points.
LocationToLevelSet location_from_dof_signs(const VectorType &local_levelset_values)
static constexpr unsigned int faces_per_cell
static VectorType::value_type get(const VectorType &V, const types::global_dof_index i)