Rheolef  7.2
an efficient C++ finite element environment
 
Loading...
Searching...
No Matches
geo_element_contains.cc
Go to the documentation of this file.
1
20// x in K ?
21// where x is a point and K a geo_element
22//
23// TODO: omega.order == 1 only: not valid for curved elements
24// ---------------------------------------------------------------------
25#include "rheolef/geo_element_contains.h"
26
27namespace rheolef { namespace details {
28
29template <class T, class M>
30bool
32 const geo_element& K,
33 const disarray<point_basic<T>,M>& node,
34 const point_basic<T>& x)
35{
36 const point_basic<T>& p0 = node.dis_at(K[0]);
37 const point_basic<T>& p1 = node.dis_at(K[1]);
38 return x[0] >= p0[0] && x[0] <= p1[0];
39}
40template <class T, class M>
41bool
43 const geo_element& K,
44 const disarray<point_basic<T>,M>& node,
45 const point_basic<T>& x)
46{
47 static const T eps = 1e3*std::numeric_limits<T>::epsilon();
48 const point_basic<T>& p0 = node.dis_at(K[0]);
49 const point_basic<T>& p1 = node.dis_at(K[1]);
50 const point_basic<T>& p2 = node.dis_at(K[2]);
51 if (orient2d( x, p1, p2) < -eps) return false;
52 if (orient2d(p0, x, p2) < -eps) return false;
53 if (orient2d(p0, p1, x) < -eps) return false;
54 return true;
55}
56template <class T, class M>
57bool
59 const geo_element& K,
60 const disarray<point_basic<T>,M>& node,
61 const point_basic<T>& x)
62{
63 static const T eps = 1e3*std::numeric_limits<T>::epsilon();
64 const point_basic<T>& p0 = node.dis_at(K[0]);
65 const point_basic<T>& p1 = node.dis_at(K[1]);
66 const point_basic<T>& p2 = node.dis_at(K[2]);
67 const point_basic<T>& p3 = node.dis_at(K[3]);
68 if (orient2d(x, p0, p1) < -eps) return false;
69 if (orient2d(x, p1, p2) < -eps) return false;
70 if (orient2d(x, p2, p3) < -eps) return false;
71 if (orient2d(x, p3, p0) < -eps) return false;
72 return true;
73}
74template <class T, class M>
75bool
77 const geo_element& K,
78 const disarray<point_basic<T>,M>& node,
79 const point_basic<T>& x)
80{
81 static const T eps = 1e3*std::numeric_limits<T>::epsilon();
82 const point_basic<T>& p0 = node.dis_at(K[0]);
83 const point_basic<T>& p1 = node.dis_at(K[1]);
84 const point_basic<T>& p2 = node.dis_at(K[2]);
85 const point_basic<T>& p3 = node.dis_at(K[3]);
86 if (orient3d( x, p1, p2, p3) < -eps) return false;
87 if (orient3d(p0, x, p2, p3) < -eps) return false;
88 if (orient3d(p0, p1, x, p3) < -eps) return false;
89 if (orient3d(p0, p1, p2, x) < -eps) return false;
90 return true;
91}
92template <class T, class M>
93bool
95 const geo_element& K,
96 const disarray<point_basic<T>,M>& node,
97 const point_basic<T>& x)
98{
100 static const T eps = 1e3*std::numeric_limits<T>::epsilon();
101 const point_basic<T>& p0 = node.dis_at(K[0]);
102 const point_basic<T>& p1 = node.dis_at(K[1]);
103 const point_basic<T>& p2 = node.dis_at(K[2]);
104 const point_basic<T>& p3 = node.dis_at(K[3]);
105 const point_basic<T>& p4 = node.dis_at(K[4]);
106 const point_basic<T>& p5 = node.dis_at(K[5]);
107 const point_basic<T>& p6 = node.dis_at(K[6]);
108 const point_basic<T>& p7 = node.dis_at(K[7]);
109 const point_basic<T>* q[8] = {&p0, &p1, &p2, &p3, &p4, &p5, &p6, &p7};
110 for (size_type loc_isid = 0, loc_nsid = 6; loc_isid < loc_nsid; loc_isid++) {
111 size_type j0 = reference_element_H::subgeo_local_node (1, 2, loc_isid, 0);
112 size_type j1 = reference_element_H::subgeo_local_node (1, 2, loc_isid, 1);
113 size_type j2 = reference_element_H::subgeo_local_node (1, 2, loc_isid, 2);
114 if (orient3d(x, *(q[j0]), *(q[j1]), *(q[j2])) < -eps) return false;
115 }
116 return true;
117}
118template <class T, class M>
119bool
121 const geo_element& K,
122 const disarray<point_basic<T>,M>& node,
123 const point_basic<T>& x)
124{
126 static const T eps = 1e3*std::numeric_limits<T>::epsilon();
127 const point_basic<T>& p0 = node.dis_at(K[0]);
128 const point_basic<T>& p1 = node.dis_at(K[1]);
129 const point_basic<T>& p2 = node.dis_at(K[2]);
130 const point_basic<T>& p3 = node.dis_at(K[3]);
131 const point_basic<T>& p4 = node.dis_at(K[4]);
132 const point_basic<T>& p5 = node.dis_at(K[5]);
133 const point_basic<T>* q[6] = {&p0, &p1, &p2, &p3, &p4, &p5};
134 for (size_type loc_isid = 0, loc_nsid = 5; loc_isid < loc_nsid; loc_isid++) {
135 size_type j0 = reference_element_P::subgeo_local_node (1, 2, loc_isid, 0);
136 size_type j1 = reference_element_P::subgeo_local_node (1, 2, loc_isid, 1);
137 size_type j2 = reference_element_P::subgeo_local_node (1, 2, loc_isid, 2);
138 if (orient3d(x, *(q[j0]), *(q[j1]), *(q[j2])) < -eps) return false;
139 }
140 return true;
141}
142template <class T, class M>
143bool
145 const geo_element& K,
146 const disarray<point_basic<T>,M>& node,
147 const point_basic<T>& x)
148{
149 switch (K.variant()) {
150 case reference_element::e: return point_belongs_to_e (K, node, x);
151 case reference_element::t: return point_belongs_to_t (K, node, x);
152 case reference_element::q: return point_belongs_to_q (K, node, x);
153 case reference_element::T: return point_belongs_to_T (K, node, x);
154 case reference_element::P: return point_belongs_to_P (K, node, x);
155 case reference_element::H: return point_belongs_to_H (K, node, x);
156 default: error_macro ("unsupported element type '" << K.name() << "'"); return false;
157 }
158}
159// ----------------------------------------------------------------------------
160// instanciation in library
161// ----------------------------------------------------------------------------
162#define _RHEOLEF_instanciation(T,M) \
163template bool contains ( \
164 const geo_element& , \
165 const disarray<point_basic<T>,M>&, \
166 const point_basic<T>&);
167
169#ifdef _RHEOLEF_HAVE_MPI
170_RHEOLEF_instanciation(Float,distributed)
171#endif // _RHEOLEF_HAVE_MPI
172#undef _RHEOLEF_instanciation
173
174}} // namespace rheolef::details
#define _RHEOLEF_instanciation(T, M, A)
Definition asr.cc:223
field::size_type size_type
Definition branch.cc:430
see the Float page for the full documentation
see the disarray page for the full documentation
Definition disarray.h:497
see the geo_element page for the full documentation
variant_type variant() const
static size_type subgeo_local_node(size_type order, size_type side_dim, size_type loc_isid, size_type loc_jsidnod)
static size_type subgeo_local_node(size_type order, size_type side_dim, size_type loc_isid, size_type loc_jsidnod)
static const variant_type H
static const variant_type q
static const variant_type e
std::vector< int >::size_type size_type
static const variant_type T
static const variant_type P
static const variant_type t
#define error_macro(message)
Definition dis_macros.h:49
Expr1::float_type T
Definition field_expr.h:230
bool contains(const geo_element &K, const disarray< point_basic< T >, M > &node, const point_basic< T > &x)
bool point_belongs_to_t(const geo_element &K, const disarray< point_basic< T >, M > &node, const point_basic< T > &x)
bool point_belongs_to_H(const geo_element &K, const disarray< point_basic< T >, M > &node, const point_basic< T > &x)
bool point_belongs_to_q(const geo_element &K, const disarray< point_basic< T >, M > &node, const point_basic< T > &x)
bool point_belongs_to_T(const geo_element &K, const disarray< point_basic< T >, M > &node, const point_basic< T > &x)
bool point_belongs_to_P(const geo_element &K, const disarray< point_basic< T >, M > &node, const point_basic< T > &x)
bool point_belongs_to_e(const geo_element &K, const disarray< point_basic< T >, M > &node, const point_basic< T > &x)
This file is part of Rheolef.
T orient2d(const point_basic< T > &a, const point_basic< T > &b, const point_basic< T > &c)
T orient3d(const point_basic< T > &a, const point_basic< T > &b, const point_basic< T > &c, const point_basic< T > &d)
Expr1::memory_type M