[P]arallel [Hi]gh-order [Li]brary for [P]DEs  Latest
Parallel High-Order Library for PDEs through hp-adaptive Discontinuous Galerkin methods
PHiLiP::GridRefinement::ElementIsotropic< dim, real > Class Template Reference

Isotropic element class. More...

#include <field.h>

Inheritance diagram for PHiLiP::GridRefinement::ElementIsotropic< dim, real >:
Collaboration diagram for PHiLiP::GridRefinement::ElementIsotropic< dim, real >:

Public Member Functions

real & scale () override
 Reference for element size. More...
 
void set_scale (const real val) override
 Set the scale for the element.
 
real get_scale () const override
 Get the scale for the element.
 
void set_anisotropic_ratio (const std::array< real, dim > &ratio) override
 Set the anisotropic ratio for each reference axis. More...
 
std::array< real, dim > get_anisotropic_ratio () override
 Get the anisotropic ratio of each reference axis as an array. More...
 
real get_anisotropic_ratio (const unsigned int j) override
 Get the anisotropic ratio corresponding to the \(j^{th}\) reference axis. More...
 
void set_unit_axis (const std::array< dealii::Tensor< 1, dim, real >, dim > &unit_axis) override
 Set the unit axes of the element. More...
 
std::array< dealii::Tensor< 1, dim, real >, dim > get_unit_axis () override
 Get unit axis directions as an array.
 
dealii::Tensor< 1, dim, real > get_unit_axis (const unsigned int j) override
 Get the \(j^{th}\) reference axis.
 
void set_axis (const std::array< dealii::Tensor< 1, dim, real >, dim > &axis) override
 Set the scaled local frame axes based on vector set (length and direction) More...
 
std::array< dealii::Tensor< 1, dim, real >, dim > get_axis () override
 Get the array of axes of the local frame field \((v_1,v_2,v_3)\).
 
dealii::Tensor< 1, dim, real > get_axis (const unsigned int j) override
 Get the vector corresponding to the \(j^{th}\) frame axis \(v_j\).
 
dealii::Tensor< 2, dim, real > get_metric () override
 Get metric matrix at point describing mapping from reference element. More...
 
dealii::Tensor< 2, dim, real > get_inverse_metric () override
 Get inverse metric matrix for the reference element. More...
 
void set_anisotropy (const std::array< real, dim > &derivative_value, const std::array< dealii::Tensor< 1, dim, real >, dim > &derivative_direction, const unsigned int order) override
 Set anisotropy from reconstructed directional derivatives. More...
 
void apply_anisotropic_limit (const real anisotropic_ratio_min, const real anisotropic_ratio_max)
 Limits the anisotropic ratios to a given bandwidth. More...
 
- Public Member Functions inherited from PHiLiP::GridRefinement::Element< dim, real >
virtual ~Element ()=default
 Destructor.
 
template<typename DoFCellAccessorType >
void set_cell (const DoFCellAccessorType &cell)
 Set the Element based on the input cell (from current mesh) More...
 

Protected Member Functions

void set_cell_internal (const typename Element< dim, real >::VertexList &vertices) override
 Set element to match geometry of input vertex set. More...
 
- Protected Member Functions inherited from PHiLiP::GridRefinement::Element< dim, real >
ChordList get_chord_list (const VertexList &vertices)
 Get the chord list from an input set of vertices.
 
virtual void set_cell_internal (const VertexList &vertices)=0
 Set element to match geometry of input vertex set. More...
 

Private Attributes

real m_scale
 element size
 

Friends

std::ostream & operator<< (std::ostream &os, const ElementIsotropic< dim, real > &element)
 Write properties of element to ostream. More...
 

Additional Inherited Members

- Public Types inherited from PHiLiP::GridRefinement::Element< dim, real >
using VertexList = std::array< dealii::Tensor< 1, dim, real >, dealii::GeometryInfo< dim >::vertices_per_cell >
 Type alias for array of vertex coordinates for element.
 
using ChordList = std::array< dealii::Tensor< 1, dim, real >, dim >
 Type alias for array of chord veectors (face center to face center) in Deal.II ordering.
 

Detailed Description

template<int dim, typename real>
class PHiLiP::GridRefinement::ElementIsotropic< dim, real >

Isotropic element class.

Specialization of the element type for the isotropic remeshing case. Contains only the scale of local element (size field) . A collection of elements is contained in the corresponding Field class. Note: virtual functions controlling axes or anisotropic ratio do nothing, assert(0).

Definition at line 161 of file field.h.

Member Function Documentation

◆ apply_anisotropic_limit()

template<int dim, typename real >
void PHiLiP::GridRefinement::ElementIsotropic< dim, real >::apply_anisotropic_limit ( const real  anisotropic_ratio_min,
const real  anisotropic_ratio_max 
)
virtual

Limits the anisotropic ratios to a given bandwidth.

First finds ratio above max, redistributes length change to maintain constant volume and scale. Then the process is repeated with a lower bound. Note: does nothing in the isotropic case.

Implements PHiLiP::GridRefinement::Element< dim, real >.

Definition at line 228 of file field.cpp.

◆ get_anisotropic_ratio() [1/2]

template<int dim, typename real >
std::array< real, dim > PHiLiP::GridRefinement::ElementIsotropic< dim, real >::get_anisotropic_ratio ( )
overridevirtual

Get the anisotropic ratio of each reference axis as an array.

Note: equals 1 for each axis in isotropic case.

Implements PHiLiP::GridRefinement::Element< dim, real >.

Definition at line 99 of file field.cpp.

◆ get_anisotropic_ratio() [2/2]

template<int dim, typename real >
real PHiLiP::GridRefinement::ElementIsotropic< dim, real >::get_anisotropic_ratio ( const unsigned int  j)
overridevirtual

Get the anisotropic ratio corresponding to the \(j^{th}\) reference axis.

Note: equals 1 for each axis in isotropic case.

Implements PHiLiP::GridRefinement::Element< dim, real >.

Definition at line 110 of file field.cpp.

◆ get_inverse_metric()

template<int dim, typename real >
dealii::Tensor< 2, dim, real > PHiLiP::GridRefinement::ElementIsotropic< dim, real >::get_inverse_metric ( )
overridevirtual

Get inverse metric matrix for the reference element.

In 2D orthogonal case, \(V = 1/h * \mathrm{diag}{\rho,1/\rho} * R(-\theta)\).

Implements PHiLiP::GridRefinement::Element< dim, real >.

Definition at line 181 of file field.cpp.

◆ get_metric()

template<int dim, typename real >
dealii::Tensor< 2, dim, real > PHiLiP::GridRefinement::ElementIsotropic< dim, real >::get_metric ( )
overridevirtual

Get metric matrix at point describing mapping from reference element.

In 2D orthorgonal case, \(V = [v,w] = h * R(\theta) * \mathrm{diag}{\rho,1/\rho}\). Under transformation, order of axes is maintained with \((i,j,k)\) vectors mapping to \((v_1,v_2,v_3)\)

Implements PHiLiP::GridRefinement::Element< dim, real >.

Definition at line 171 of file field.cpp.

◆ scale()

template<int dim, typename real >
real & PHiLiP::GridRefinement::ElementIsotropic< dim, real >::scale ( )
overridevirtual

Reference for element size.

Allows direct read/write of scale of mean element axis. Measure of element (length, area, volume) will be \(scale^{dim}\)

Implements PHiLiP::GridRefinement::Element< dim, real >.

Definition at line 73 of file field.cpp.

◆ set_anisotropic_ratio()

template<int dim, typename real >
void PHiLiP::GridRefinement::ElementIsotropic< dim, real >::set_anisotropic_ratio ( const std::array< real, dim > &  ratio)
overridevirtual

Set the anisotropic ratio for each reference axis.

Requires array of order matching axis definition. Each reference axis will have length \(l = \alpha * scale\). Note: does nothing in the isotropic case.

Implements PHiLiP::GridRefinement::Element< dim, real >.

Definition at line 92 of file field.cpp.

◆ set_anisotropy()

template<int dim, typename real >
void PHiLiP::GridRefinement::ElementIsotropic< dim, real >::set_anisotropy ( const std::array< real, dim > &  derivative_value,
const std::array< dealii::Tensor< 1, dim, real >, dim > &  derivative_direction,
const unsigned int  order 
)
overridevirtual

Set anisotropy from reconstructed directional derivatives.

Based on Dolejsi's method for simplices, uses values obtained from reconstructed polynomial on neighbouring cell patch.

Implements PHiLiP::GridRefinement::Element< dim, real >.

Definition at line 218 of file field.cpp.

◆ set_axis()

template<int dim, typename real >
void PHiLiP::GridRefinement::ElementIsotropic< dim, real >::set_axis ( const std::array< dealii::Tensor< 1, dim, real >, dim > &  axis)
overridevirtual

Set the scaled local frame axes based on vector set (length and direction)

Describe the axes of reference parrelogram or parrelopiped at point. Note: does nothing in the isotropic case.

Implements PHiLiP::GridRefinement::Element< dim, real >.

Definition at line 146 of file field.cpp.

◆ set_cell_internal()

template<int dim, typename real >
void PHiLiP::GridRefinement::ElementIsotropic< dim, real >::set_cell_internal ( const typename Element< dim, real >::VertexList vertices)
overrideprotected

Set element to match geometry of input vertex set.

Vertices describe the tensor product element in Deal.II ordering. Internal function used in handling of set_cell().

Definition at line 201 of file field.cpp.

◆ set_unit_axis()

template<int dim, typename real >
void PHiLiP::GridRefinement::ElementIsotropic< dim, real >::set_unit_axis ( const std::array< dealii::Tensor< 1, dim, real >, dim > &  unit_axis)
overridevirtual

Set the unit axes of the element.

Requires an ordered array of dealii::Tensor with vector axes. If none unit values are provided, will be rescaled and factored in to anisotropic ratios. Note: does nothing in the isotropic case.

Implements PHiLiP::GridRefinement::Element< dim, real >.

Definition at line 117 of file field.cpp.

Friends And Related Function Documentation

◆ operator<<

template<int dim, typename real>
std::ostream& operator<< ( std::ostream &  os,
const ElementIsotropic< dim, real > &  element 
)
friend

Write properties of element to ostream.

Used in Field.serialize(os) to provide summary of field.

Definition at line 266 of file field.h.


The documentation for this class was generated from the following files: