geo_element.cc
Go to the documentation of this file.
#include "rheolef/geo_element.h"
#include "geo_element_aux.icc"
namespace rheolef {
// --------------------------------------------------------------------------
// i/o
// --------------------------------------------------------------------------
char
skip_blancs_and_tabs (std::istream& is)
{
if (!is.good()) return 0;
char c = is.peek();
while (is.good() && (c == ' ' || c == '\t')) {
is >> c;
} while (is.good() && (c == ' ' || c == '\t'));
return c;
}
void
geo_element::get (std::istream& is)
{
geo_element& K = *this; // for more readability
typedef geo_element::size_type size_type;
extern char skip_blancs_and_tabs (std::istream&);
char c;
is >> std::ws >> c;
if (!isdigit(c)) {
// have an element type: 'e', 't' etc
geo_element::variant_type variant = reference_element::variant (c);
check_macro (variant != reference_element::max_variant, "undefined element variant `"<<c<<"'");
size_type order = 1;
char c2;
is >> std::ws >> c2;
if (c2 == 'p') {
// have a high order spec as "p2"
is >> order;
} else {
is.unget(); // or is.putback(c); ?
}
K.reset (variant, order);
for (size_type i = 0, n = K.n_node(); i < n; i++) {
is >> K[i];
}
return;
}
// here, starts with a digit: triangle (e.g. "21 22 23"), edge (e.g. "12 13") or point (e.g. "11")
// end of element input: with end-of-line or non-digit input char
// Note: with 4 vertices, there is an ambiguity between quadrangle and thetraedra
size_type tmp [3];
size_type nvertex = 0;
while (is.good() && isdigit(c) && nvertex < 3) {
is.unget(); // or is.putback(c); ?
is >> tmp [nvertex];
nvertex++;
c = skip_blancs_and_tabs (is);
}
size_type variant = reference_element::max_variant;
switch (nvertex) {
case 1: variant = reference_element::p; break;
case 2: variant = reference_element::e; break;
case 3: variant = reference_element::t; break;
default: error_macro ("unexpected element with "<<nvertex<<" vertices");
}
K.reset (variant, 1);
for (size_type i = 0, n = K.n_node(); i < n; i++) {
K[i] = tmp[i];
}
return;
}
void
geo_element::put (std::ostream& os) const
{
typedef geo_element::size_type size_type;
os << name() << "\t";
if (order() > 1) {
os << "p" << order() << " ";
}
for (size_type loc_inod = 0, loc_nnod = n_node(); loc_inod < loc_nnod; loc_inod++) {
if (loc_inod != 0) os << " ";
os << operator[](loc_inod);
}
}
// --------------------------------------------------------------------------
bool
geo_element::get_orientation_and_shift (
const geo_element& S,
orientation_type& orient,
shift_type& shift) const
{
check_macro (S.dimension() == dimension(),
"get_orientation_and_shift: elements have different dimensions "<<dimension()<<" and "<<S.dimension());
check_macro (S.size() == size(),
"get_orientation_and_shift: elements have different sizes "<<size()<<" and "<<S.size());
if (S.dimension() == 0) {
orient = 1;
shift = 0;
return true;
}
if (S.dimension() == 1) {
orient = (operator[](0) == S[0]) ? 1 : -1;
shift = 0;
return true;
}
if (S.dimension() == 2) {
size_type n = size();
// 3d face: triangle or quadrangle
for (shift = 0; shift < shift_type(n); shift++) {
if (operator[](0) == S[shift]) break;
}
if (shift == shift_type(n)) {
orient = 0;
shift = std::numeric_limits<shift_type>::max();
return false;
}
orient = (operator[](1) == S[(shift+1)%n]) ? 1 : -1;
return true;
}
// S.dimension() == 3: TODO volumic domain: tetra can be rotated ?
orient = 1;
shift = 0;
return true;
}
// if edge (dis_iv0,dis_iv1) has the same orientation as current edge
// assume that current geo_element is an edge that contains the same vertices
geo_element::orientation_type
geo_element::get_edge_orientation (
size_type dis_iv0, size_type dis_iv1) const
{
return (operator[](0) == dis_iv0) ? 1 : -1;
}
// for a triangular side (dis_iv0,dis_iv1,dis_iv2)
// assume that current geo_element is a triangle that contains the same vertices
void
geo_element::get_orientation_and_shift (
size_type dis_iv0, size_type dis_iv1, size_type dis_iv2,
orientation_type& orient,
shift_type& shift) const
{
geo_element_get_orientation_and_shift (*this, dis_iv0, dis_iv1, dis_iv2, orient, shift);
}
// for a quadrangular side (dis_iv0,dis_iv1,dis_iv2,dis_iv3)
// assume that current geo_element is a quadrangular that contains the same vertices
void
geo_element::get_orientation_and_shift (
size_type dis_iv0, size_type dis_iv1, size_type dis_iv2, size_type dis_iv3,
orientation_type& orient,
shift_type& shift) const
{
geo_element_get_orientation_and_shift (*this, dis_iv0, dis_iv1, dis_iv2, dis_iv3, orient, shift);
}
// let K=(*this) one element and S one side of K
// S could have the good or the opposite orientation on K:
// gives its sign, local side index and node shift (in 3D)
// Note: assume that S is a valid side of K
geo_element::orientation_type
geo_element::get_side_informations (
const geo_element& S,
size_type& loc_isid,
size_type& shift) const
{
loc_isid = shift = 0;
check_macro (S.dimension() + 1 == dimension(),
"get_side_orientation: side have unexpected dimension "<<S.dimension()<<": "
<<dimension()<<"-1 was expected");
const geo_element& K = *this;
size_type side_dim = S.dimension();
for (size_type loc_nsid = K.n_subgeo(side_dim); loc_isid < loc_nsid; loc_isid++) {
size_type sid_nloc = K.subgeo_size (side_dim, loc_isid);
if (sid_nloc != S.size()) continue;
for (shift = 0; shift < sid_nloc; shift++) {
size_type loc_jv = K.subgeo_local_vertex (side_dim, loc_isid, shift);
if (K[loc_jv] != S[0]) continue;
// one node matches with S[0] on K: loc_isid and shift
// check others nodes in a first rotation direction:
bool matches = true;
for (size_type sid_kloc = 1; sid_kloc < sid_nloc; sid_kloc++) {
size_type loc_kv = K.subgeo_local_vertex (side_dim, loc_isid, (shift+sid_kloc) % sid_nloc);
if (K[loc_kv] != S[sid_kloc]) { matches = false; break; }
}
if (matches) {
if (side_dim == 1 && shift == 1) {
// shift=1 for an edge means a change of orientation of this edge (for DG)
shift = 0;
return -1;
}
return 1;
}
// check others nodes in the opposite rotation direction:
matches = true;
for (size_type sid_kloc = 1; sid_kloc < sid_nloc; sid_kloc++) {
size_type loc_kv = K.subgeo_local_vertex (side_dim, loc_isid, (shift+sid_nloc-sid_kloc) % sid_nloc);
if (K[loc_kv] != S[sid_kloc]) { matches = false; break; }
}
if (matches) { return -1; }
}
}
fatal_macro ("get_side_orientation: side is not part of the element");
return 0; // not reached
}
void
geo_element::get_side_informations (
const geo_element& S,
side_information_type& sid) const
{
sid.orient = get_side_informations (S, sid.loc_isid, sid.shift);
sid.dim = S.dimension();
sid.n_vertex = subgeo_size (sid.dim, sid.loc_isid);
sid.hat.set_variant (sid.n_vertex, sid.dim);
}
geo_element::orientation_type
geo_element::get_side_orientation (const geo_element& S) const
{
size_type loc_isid, shift;
return get_side_informations (S, loc_isid, shift);
}
// =========================================================================
// fix rotation and orientation on a 2d edge or 3d face
// =========================================================================
// --------------
// 1) edges
// --------------
geo_element::size_type
geo_element::fix_edge_indirect (
const geo_element& K,
size_type loc_iedg,
size_type loc_iedg_j,
size_type order)
{
return geo_element_fix_edge_indirect (K, loc_iedg, loc_iedg_j, order);
}
// --------------
// 2) triangles
// --------------
void
geo_element::loc_tri_inod2lattice (
size_type loc_tri_inod,
size_type order,
point_basic<size_type>& ij_lattice)
{
geo_element_loc_tri_inod2lattice (loc_tri_inod, order, ij_lattice);
}
geo_element::size_type
geo_element::fix_triangle_indirect (
orientation_type orient,
shift_type shift,
size_type order,
size_type loc_itri_j)
{
return geo_element_fix_triangle_indirect (orient, shift, order, loc_itri_j);
}
geo_element::size_type
geo_element::fix_triangle_indirect (
const geo_element& K,
size_type loc_ifac,
size_type loc_itri_j,
size_type order)
{
return geo_element_fix_triangle_indirect (K, loc_ifac, loc_itri_j, order);
}
// --------------
// 3) quadrangles
// --------------
void
geo_element::loc_qua_inod2lattice (
size_type loc_qua_inod,
size_type order,
point_basic<size_type>& ij_lattice)
{
geo_element_loc_qua_inod2lattice (loc_qua_inod, order, ij_lattice);
}
geo_element::size_type
geo_element::fix_quadrangle_indirect (
orientation_type orient,
shift_type shift,
size_type order,
size_type loc_iqua_j)
{
return geo_element_fix_quadrangle_indirect (orient, shift, order, loc_iqua_j);
}
geo_element::size_type
geo_element::fix_quadrangle_indirect (
const geo_element& K,
size_type loc_ifac,
size_type loc_iqua_j,
size_type order)
{
return geo_element_fix_quadrangle_indirect (K, loc_ifac, loc_iqua_j, order);
}
geo_element::size_type
geo_element::fix_indirect (
const geo_element& K,
size_type subgeo_variant,
size_type loc_ige,
size_type loc_comp_idof_on_subgeo,
size_type order)
{
switch (subgeo_variant) {
case reference_element::e: return geo_element::fix_edge_indirect (K, loc_ige, loc_comp_idof_on_subgeo, order);
case reference_element::t: return geo_element::fix_triangle_indirect (K, loc_ige, loc_comp_idof_on_subgeo, order);
case reference_element::q: return geo_element::fix_quadrangle_indirect (K, loc_ige, loc_comp_idof_on_subgeo, order);
default: return loc_comp_idof_on_subgeo;
}
}
} // namespace rheolef
rheolef::geo_element
see the geo_element page for the full documentation
Definition: geo_element.h:102
rheolef::geo_element::operator[]
size_type & operator[](size_type loc_inod)
Definition: geo_element.h:183
rheolef::geo_element::orientation_type
geo_element_indirect::orientation_type orientation_type
Definition: geo_element.h:134
rheolef::geo_element::get_orientation_and_shift
bool get_orientation_and_shift(const geo_element &S, orientation_type &orient, shift_type &shift) const
return orientation and shift between *this element and S
Definition: geo_element.cc:114
rheolef::geo_element::n_node
size_type n_node() const
Definition: geo_element.h:170
rheolef::geo_element::subgeo_local_vertex
size_type subgeo_local_vertex(size_type subgeo_dim, size_type i_subgeo, size_type i_subgeo_vertex) const
Definition: geo_element.h:223
rheolef::geo_element::fix_edge_indirect
static size_type fix_edge_indirect(const geo_element &K, size_type loc_iedg, size_type loc_iedg_j, size_type order)
Definition: geo_element.cc:253
rheolef::geo_element::size
size_type size() const
Definition: geo_element.h:168
rheolef::geo_element::fix_triangle_indirect
static size_type fix_triangle_indirect(const geo_element &K, size_type loc_itri, size_type loc_itri_j, size_type order)
Definition: geo_element.cc:282
rheolef::geo_element::size_type
reference_element::size_type size_type
Definition: geo_element.h:125
rheolef::geo_element::put
void put(std::ostream &is) const
Definition: geo_element.cc:93
rheolef::geo_element::dimension
size_type dimension() const
Definition: geo_element.h:167
rheolef::geo_element::get
void get(std::istream &os)
Definition: geo_element.cc:42
rheolef::geo_element::subgeo_size
size_type subgeo_size(size_type subgeo_dim, size_type loc_isid) const
Definition: geo_element.h:221
rheolef::geo_element::get_side_orientation
orientation_type get_side_orientation(const geo_element &S) const
Definition: geo_element.cc:240
rheolef::geo_element::fix_indirect
static size_type fix_indirect(const geo_element &K, size_type subgeo_variant, size_type loc_ige, size_type loc_comp_idof_on_subgeo, size_type order)
Definition: geo_element.cc:320
rheolef::geo_element::variant
variant_type variant() const
Definition: geo_element.h:161
rheolef::geo_element::loc_qua_inod2lattice
static void loc_qua_inod2lattice(size_type loc_qua_inod, size_type order, point_basic< size_type > &ij_lattice)
Definition: geo_element.cc:294
rheolef::geo_element::get_edge_orientation
orientation_type get_edge_orientation(size_type dis_iv0, size_type dis_iv1) const
Definition: geo_element.cc:155
rheolef::geo_element::loc_tri_inod2lattice
static void loc_tri_inod2lattice(size_type loc_tri_inod, size_type order, point_basic< size_type > &ij_lattice)
Definition: geo_element.cc:265
rheolef::geo_element::shift_type
geo_element_indirect::shift_type shift_type
Definition: geo_element.h:135
rheolef::geo_element::get_side_informations
orientation_type get_side_informations(const geo_element &S, size_type &loc_isid, size_type &shift) const
Definition: geo_element.cc:185
rheolef::geo_element::n_subgeo
size_type n_subgeo(size_type subgeo_dim) const
Definition: geo_element.h:212
rheolef::geo_element::fix_quadrangle_indirect
static size_type fix_quadrangle_indirect(const geo_element &K, size_type loc_iqua, size_type loc_iqua_j, size_type order)
Definition: geo_element.cc:311
rheolef::geo_element::reset
virtual void reset(variant_type variant, size_type order)=0
rheolef::geo_element::name
char name() const
Definition: geo_element.h:169
rheolef::geo_element::variant_type
reference_element::variant_type variant_type
Definition: geo_element.h:126
rheolef::geo_element::order
size_type order() const
Definition: geo_element.h:162
rheolef::point_basic
Definition: point.h:87
rheolef::reference_element::q
static const variant_type q
Definition: reference_element.h:78
rheolef::reference_element::e
static const variant_type e
Definition: reference_element.h:76
rheolef::reference_element::set_variant
void set_variant(variant_type x)
Definition: reference_element.h:93
rheolef::reference_element::max_variant
static const variant_type max_variant
Definition: reference_element.h:82
rheolef::reference_element::p
static const variant_type p
Definition: reference_element.h:75
rheolef::reference_element::variant
variant_type variant() const
Definition: reference_element.h:99
rheolef::reference_element::t
static const variant_type t
Definition: reference_element.h:77
error_macro
#define error_macro(message)
Definition: dis_macros.h:49
fatal_macro
#define fatal_macro(message)
Definition: dis_macros.h:33
check_macro
check_macro(expr1.have_homogeneous_space(Xh1), "dual(expr1,expr2); expr1 should have homogeneous space. HINT: use dual(interpolate(Xh, expr1),expr2)")
geo_element_aux.icc
mkgeo_ball.c
int c
Definition: mkgeo_ball.sh:153
mkgeo_ball.n
int n
Definition: mkgeo_ball.sh:150
mkgeo_ball.order
order
Definition: mkgeo_ball.sh:343
mkgeo_ball.tmp
tmp
Definition: mkgeo_ball.sh:380
mkgeo_contraction.c2
c2
Definition: mkgeo_contraction.sh:199
rheolef
This file is part of Rheolef.
Definition: compiler_eigen.h:37
rheolef::skip_blancs_and_tabs
char skip_blancs_and_tabs(std::istream &is)
Definition: geo_element.cc:32
rheolef::side_information_type
Definition: reference_element_face_transformation.h:37
rheolef::side_information_type::hat
reference_element hat
Definition: reference_element_face_transformation.h:43
rheolef::side_information_type::dim
size_t dim
Definition: reference_element_face_transformation.h:41
rheolef::side_information_type::n_vertex
size_t n_vertex
Definition: reference_element_face_transformation.h:42
rheolef::side_information_type::orient
geo_element_indirect::orientation_type orient
Definition: reference_element_face_transformation.h:40
rheolef::side_information_type::loc_isid
size_t loc_isid
Definition: reference_element_face_transformation.h:38
rheolef::side_information_type::shift
size_t shift
Definition: reference_element_face_transformation.h:39