MAIA bb96820c
Multiphysics at AIA
Loading...
Searching...
No Matches
Public Member Functions | Private Member Functions | List of all members
CHECKNORMAL< nDim_ > Class Template Reference

#include <lscartesiansolver.h>

Inheritance diagram for CHECKNORMAL< nDim_ >:
[legend]

Public Member Functions

void rotation (const MFloat q[3], MFloat r[3], MFloat tr[3][3])
 
MInt PointInsideTriangle (GeometryElement< nDim_ > *el, MFloat q[3], MFloat tr[3][3])
 

Private Member Functions

MFloat dot (const MFloat p[3], const MFloat q[3])
 
void cross (const MFloat p[3], const MFloat q[3], MFloat r[3])
 
void transformationmatrix (GeometryElement< nDim_ > *el, MFloat tr[3][3])
 

Detailed Description

template<MInt nDim_>
class CHECKNORMAL< nDim_ >

Definition at line 38 of file lscartesiansolver.h.

Member Function Documentation

◆ cross()

template<MInt nDim_>
void CHECKNORMAL< nDim_ >::cross ( const MFloat  p[3],
const MFloat  q[3],
MFloat  r[3] 
)
inlineprivate

Definition at line 43 of file lscartesiansolver.h.

43 {
44 r[0] = (p[1] * q[2]) - (q[1] * p[2]);
45 r[1] = (p[2] * q[0]) - (q[2] * p[0]);
46 r[2] = (p[0] * q[1]) - (q[0] * p[1]);
47 }
maia::fc::cell::p
constexpr std::underlying_type< FcCell >::type p(const FcCell property)
Converts property name to underlying integer value.
Definition: fccellproperties.h:32

◆ dot()

template<MInt nDim_>
MFloat CHECKNORMAL< nDim_ >::dot ( const MFloat  p[3],
const MFloat  q[3] 
)
inlineprivate

Definition at line 41 of file lscartesiansolver.h.

41{ return p[0] * q[0] + p[1] * q[1] + p[2] * q[2]; }

◆ PointInsideTriangle()

template<MInt nDim_>
MInt CHECKNORMAL< nDim_ >::PointInsideTriangle ( GeometryElement< nDim_ > *  el,
MFloat  q[3],
MFloat  tr[3][3] 
)
inline

Definition at line 132 of file lscartesiansolver.h.

132 {
133 // Set transformation matrix
134 transformationmatrix(el, tr);
135
136 // Copy vertices to vv[][]
137 MFloat vv[3][3] = {{F0, F0, F0}, {F0, F0, F0}, {F0, F0, F0}};
138
139 MFloat noVertices;
140 noVertices = nDim_;
141
142 for(MInt i = 0; i < noVertices; i++)
143 for(MInt j = 0; j < nDim_; j++) {
144 vv[i][j] = el->m_vertices[i][j];
145 }
146
147 // Rotation of vertices
148 MFloat ro[3] = {0.0, 0.0, 0.0}, uv[6] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0}, quv[2] = {0.0, 0.0};
149 for(MInt i = 0; i < noVertices; i++) {
150 rotation(vv[i], ro, tr);
151 for(MInt j = 0; j < nDim_ - 1; j++) {
152 uv[(nDim_ - 1) * i + j] = ro[j];
153 }
154 }
155
156 // Rotation of q
157 rotation(q, ro, tr);
158 for(MInt j = 0; j < nDim_ - 1; j++)
159 quv[j] = ro[j];
160
161 // Move origin to q
162 for(MInt j = 0; j < nDim_ - 1; j++) {
163 for(MInt i = 0; i < noVertices; i++)
164 uv[(nDim_ - 1) * i + j] -= quv[j];
165 }
166
167 // Define intersection
168 IF_CONSTEXPR(nDim_ == 2) {
169 MInt SH = 0, NSH = 0;
170
171 if(uv[0] > 0)
172 SH = 1;
173 else if(approx(uv[0], 0.0, MFloatEps))
174 SH = 0;
175 else
176 SH = -1;
177
178 if(uv[1] > 0)
179 NSH = 1;
180 else if(approx(uv[1], 0.0, MFloatEps))
181 NSH = 0;
182 else
183 NSH = -1;
184
185 if((SH != NSH) || (SH == 0) || (NSH == 0)) return 1;
186 }
187 else {
188 MInt NC = 0;
189
190 for(MInt i = 0; i < noVertices; i++) {
191 MInt ni = (i + 1) % 3; // next i
192 // cast ray in +u direction
193 MInt SH = 0, NSH = 0;
194
195 if(uv[i * 2 + 1] > 0)
196 SH = 1;
197 else if(approx(uv[i * 2 + 1], 0.0, MFloatEps))
198 SH = 0;
199 else
200 SH = -1;
201
202 if(uv[ni * 2 + 1] > 0)
203 NSH = 1;
204 else if(approx(uv[ni * 2 + 1], 0.0, MFloatEps))
205 NSH = 0;
206 else
207 NSH = -1;
208
209 // predict intersection
210 if(SH == 0 && NSH == 0) {
211 NC = 1;
212 break;
213 } // horizontal line through origin
214 else if((SH == 0 && approx(uv[i * 2 + 0], 0.0, MFloatEps))
215 || (NSH == 0 && approx(uv[ni * 2 + 0], 0.0, MFloatEps))) {
216 NC = 1;
217 break;
218 } else if(approx((uv[i * 2 + 0]
219 - uv[i * 2 + 1] * (uv[ni * 2 + 0] - uv[i * 2 + 0]) / (uv[ni * 2 + 1] - uv[i * 2 + 1])),
220 0.0, MFloatEps)) {
221 NC = 1;
222 break;
223 } // line through origin
224 else if(SH != NSH) {
225 if((uv[i * 2 + 0] > 0) && (uv[ni * 2 + 0] > 0))
226 NC++; // there is intersection in +u
227 else if((uv[i * 2 + 0] > 0) || (uv[ni * 2 + 0] > 0)) { // there can be an intersection in +u
228 // need to calculate
229 if((uv[i * 2 + 0] - uv[i * 2 + 1] * (uv[ni * 2 + 0] - uv[i * 2 + 0]) / (uv[ni * 2 + 1] - uv[i * 2 + 1]))
230 > 0)
231 NC++;
232 }
233 }
234 }
235 if((NC % 2) != 0) {
236 return 1;
237 }
238 }
239 return 0;
240 }
CHECKNORMAL::transformationmatrix
void transformationmatrix(GeometryElement< nDim_ > *el, MFloat tr[3][3])
Definition: lscartesiansolver.h:50
CHECKNORMAL::rotation
void rotation(const MFloat q[3], MFloat r[3], MFloat tr[3][3])
Definition: lscartesiansolver.h:122
GeometryElement::m_vertices
MFloat ** m_vertices
Definition: geometryelement.h:25
approx
MBool approx(const T &, const U &, const T)
Definition: functions.h:272
MInt
int32_t MInt
Definition: maiatypes.h:62
MFloat
double MFloat
Definition: maiatypes.h:52

◆ rotation()

template<MInt nDim_>
void CHECKNORMAL< nDim_ >::rotation ( const MFloat  q[3],
MFloat  r[3],
MFloat  tr[3][3] 
)
inline

Definition at line 122 of file lscartesiansolver.h.

122 {
123 for(MInt i = 0; i < 3; i++) {
124 r[i] = 0;
125 for(MInt j = 0; j < 3; j++) {
126 r[i] += tr[i][j] * q[j];
127 }
128 }
129 }

◆ transformationmatrix()

template<MInt nDim_>
void CHECKNORMAL< nDim_ >::transformationmatrix ( GeometryElement< nDim_ > *  el,
MFloat  tr[3][3] 
)
inlineprivate

Definition at line 50 of file lscartesiansolver.h.

50 {
51 // Step 1: get vertices in global coordinates
52 MFloat v[3][3] = {{F0, F0, F0}, {F0, F0, F0}, {F0, F0, F0}};
53 MInt noVertices = nDim_;
54
55 for(MInt i = 0; i < noVertices; i++) {
56 for(MInt j = 0; j < nDim_; j++) {
57 v[i][j] = el->m_vertices[i][j];
58 }
59 }
60
61 // Step 2: define the transformation matrix for the current geometrical element
62 IF_CONSTEXPR(nDim_ == 2) {
63 MFloat theta_rad = 0;
64 MFloat dx = v[1][0] - v[0][0];
65 if(approx(dx, 0.0, MFloatEps))
66 theta_rad = 3.141592654 / 2.0;
67 else
68 theta_rad = atan((v[1][1] - v[0][1]) / dx);
69
70 tr[0][0] = cos(theta_rad);
71 tr[0][1] = sin(theta_rad);
72 tr[2][0] = 0;
73 tr[1][0] = -sin(theta_rad);
74 tr[1][1] = cos(theta_rad);
75 tr[2][1] = 0;
76 tr[0][2] = 0;
77 tr[1][2] = 0;
78 tr[2][2] = 0;
79 }
80 else {
81 // Euclidean frame
82 MFloat xm[3] = {1, 0, 0};
83 MFloat ym[3] = {0, 1, 0};
84 MFloat zm[3] = {0, 0, 1};
85 // n subscript means rotated coordinates
86 // set zn = normal vector
87 MFloat zn[3];
88 for(MInt i = 0; i < 3; i++) {
89 zn[i] = el->m_normal[i];
90 }
91
92 // set xn to the direction of first edge
93 MFloat xn[3];
94 for(MInt i = 0; i < 3; i++) {
95 xn[i] = v[1][i] - v[0][i];
96 }
97
98 // Normalize xn
99 MFloat mag = sqrt(pow(xn[0], 2) + pow(xn[1], 2) + pow(xn[2], 2));
100 for(MFloat& i : xn) {
101 i = i / mag;
102 }
103
104 // yn is the outer product of zn x xn
105 MFloat yn[3];
106 cross(zn, xn, yn);
107 // defining transformation matrix
108 tr[0][0] = dot(xn, xm);
109 tr[0][1] = dot(xn, ym);
110 tr[0][2] = dot(xn, zm);
111 tr[1][0] = dot(yn, xm);
112 tr[1][1] = dot(yn, ym);
113 tr[1][2] = dot(yn, zm);
114 tr[2][0] = dot(zn, xm);
115 tr[2][1] = dot(zn, ym);
116 tr[2][2] = dot(zn, zm);
117 }
118 }
CHECKNORMAL::dot
MFloat dot(const MFloat p[3], const MFloat q[3])
Definition: lscartesiansolver.h:41
CHECKNORMAL::cross
void cross(const MFloat p[3], const MFloat q[3], MFloat r[3])
Definition: lscartesiansolver.h:43
maia::filter::slope::cos
T cos(const T a, const T b, const T x)
Cosine slope filter.
Definition: filter.h:125
The documentation for this class was generated from the following file: