TrioCFD 1.9.9_beta
TrioCFD documentation
Loading...
Searching...
No Matches
Tetra_VEF.cpp
1/****************************************************************************
2* Copyright (c) 2026, CEA
3* All rights reserved.
4*
5* Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
6* 1. Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
7* 2. Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.
8* 3. Neither the name of the copyright holder nor the names of its contributors may be used to endorse or promote products derived from this software without specific prior written permission.
9*
10* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
11* IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS;
12* OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
13*
14*****************************************************************************/
15
16#include <Tetra_VEF.h>
17#include <Domaine.h>
18#include <Domaine_VEF.h>
19#include <Champ_P1NC.h>
20
21Implemente_instanciable_sans_constructeur(Tetra_VEF,"Tetra_VEF",Elem_VEF_base);
22
23// printOn and readOn
24
25
27{
28 return s << que_suis_je() << finl;
29}
30
32{
33 return s ;
34}
35/*! @brief Returns for sub-facet fa7: for j=0,j=1: the local indices of the 2 faces surrounding fa7.
36 *
37 * For j=2,j=3: the local indices of the tetrahedron vertices belonging to fa7.
38 *
39 */
41{
42 int tmp[4][6]=
43 {
44 {0, 0, 0, 1, 1, 2},
45 {1, 2, 3, 2, 3, 3},
46 {2, 1, 1, 3, 2, 1},
47 {3, 3, 2, 0, 0, 0}
48 };
49 KEL_.resize(4,6);
50 for (int i=0; i<4; i++)
51 for (int j=0; j<6; j++)
52 KEL_(i,j)=tmp[i][j];
53}
54
55void Tetra_VEF::creer_face_normales(DoubleTab& tab_Face_normales,
56 const IntTab& tab_Face_sommets,
57 const IntTab& tab_Face_voisins,
58 const IntTab& tab_elem_faces,
59 const Domaine& domaine_geom) const
60{
61 int nb_face_tot = tab_Face_normales.dimension_tot(0);
62 CDoubleTabView les_coords = domaine_geom.coord_sommets().view_ro();
63 CIntTabView Face_sommets = tab_Face_sommets.view_ro();
64 CIntTabView Face_voisins = tab_Face_voisins.view_ro();
65 CIntTabView elem_faces = tab_elem_faces.view_ro();
66 DoubleTabView Face_normales = tab_Face_normales.view_rw();
67 Kokkos::parallel_for(start_gpu_timer(__KERNEL_NAME__), range_1D(0, nb_face_tot), KOKKOS_LAMBDA(const int num_Face)
68 {
69 int n0 = Face_sommets(num_Face, 0);
70 int n1 = Face_sommets(num_Face, 1);
71 int n2 = Face_sommets(num_Face, 2);
72
73 double x1 = les_coords(n0, 0) - les_coords(n1, 0);
74 double y1 = les_coords(n0, 1) - les_coords(n1, 1);
75 double z1 = les_coords(n0, 2) - les_coords(n1, 2);
76
77 double x2 = les_coords(n2, 0) - les_coords(n1, 0);
78 double y2 = les_coords(n2, 1) - les_coords(n1, 1);
79 double z2 = les_coords(n2, 2) - les_coords(n1, 2);
80
81 double nx = (y1 * z2 - y2 * z1) / 2;
82 double ny = (-x1 * z2 + x2 * z1) / 2;
83 double nz = (x1 * y2 - x2 * y1) / 2;
84
85 // Orientation of the normal from elem1 to elem2:
86 // find the vertex of elem1 that does not lie on the face
87 int elem1 = Face_voisins(num_Face, 0);
88 int f0 = elem_faces(elem1, 0);
89 if (f0 == num_Face)
90 f0 = elem_faces(elem1, 1);
91
92 int no4 = Face_sommets(f0, 0);
93 if (no4 == n0 || no4 == n1 || no4 == n2)
94 {
95 no4 = Face_sommets(f0, 1);
96 if (no4 == n0 || no4 == n1 || no4 == n2)
97 no4 = Face_sommets(f0, 2);
98 }
99
100 x1 = les_coords(no4, 0) - les_coords(n0, 0);
101 y1 = les_coords(no4, 1) - les_coords(n0, 1);
102 z1 = les_coords(no4, 2) - les_coords(n0, 2);
103
104 double sign = ((nx * x1 + ny * y1 + nz * z1) > 0) ? -1. : 1.;
105 Face_normales(num_Face, 0) = sign * nx;
106 Face_normales(num_Face, 1) = sign * ny;
107 Face_normales(num_Face, 2) = sign * nz;
108 });
109 end_gpu_timer(__KERNEL_NAME__);
110}
111
112/*! @brief Fills the face_normales array in the Domaine_VEF.
113 *
114 */
116 const IntVect& tab_rang_elem_non_std) const
117{
118 const Domaine& domaine_geom = dom_VEF.domaine();
119 auto& facette_normales = const_cast<Domaine_VEF&>(dom_VEF).facette_normales();
120 int nb_elem_tot = domaine_geom.nb_elem_tot();
121
122 if (facette_normales.dimension(0) != nb_elem_tot)
123 facette_normales.resize(nb_elem_tot,6,3);
124
125 CDoubleTabView les_coords = domaine_geom.coord_sommets().view_ro();
126 CIntTabView les_Polys = domaine_geom.les_elems().view_ro();
127 CIntArrView rang_elem_non_std = tab_rang_elem_non_std.view_ro();
128 CIntTabView KEL = KEL_.view_ro();
129 DoubleTabView3 facette_normale = facette_normales.view_rw<3>();
130 Kokkos::parallel_for(start_gpu_timer(__KERNEL_NAME__), range_1D(0, nb_elem_tot), KOKKOS_LAMBDA(const int i)
131 {
132 if (rang_elem_non_std(i) == -1)
133 {
134 int num_som[4];
135 double x[4][3];
136 num_som[0] = les_Polys(i, 0);
137 num_som[1] = les_Polys(i, 1);
138 num_som[2] = les_Polys(i, 2);
139 num_som[3] = les_Polys(i, 3);
140 for (int s = 0; s < 4; s++)
141 for (int d = 0; d < 3; d++)
142 x[s][d] = les_coords(num_som[s], d);
143 double xg[3];
144 xg[0] = 0.25*(x[0][0]+x[1][0]+x[2][0]+x[3][0]);
145 xg[1] = 0.25*(x[0][1]+x[1][1]+x[2][1]+x[3][1]);
146 xg[2] = 0.25*(x[0][2]+x[1][2]+x[2][2]+x[3][2]);
147 for (int fa7 = 0; fa7 < 6; fa7++)
148 {
149 // fa7 has vertices kel(2,fa7), kel(3,fa7), "G"
150 double u[3], v[3], pv[3], xj0[3];
151 u[0] = x[KEL(2,fa7)][0]-xg[0];
152 u[1] = x[KEL(2,fa7)][1]-xg[1];
153 u[2] = x[KEL(2,fa7)][2]-xg[2];
154 v[0] = x[KEL(3,fa7)][0]-xg[0];
155 v[1] = x[KEL(3,fa7)][1]-xg[1];
156 v[2] = x[KEL(3,fa7)][2]-xg[2];
157 // inline prodvect (not KOKKOS_INLINE_FUNCTION)
158 pv[0] = u[1]*v[2]-u[2]*v[1];
159 pv[1] = u[2]*v[0]-u[0]*v[2];
160 pv[2] = u[0]*v[1]-u[1]*v[0];
161 // Normal orientation:
162 xj0[0] = x[KEL(0,fa7)][0]-x[KEL(1,fa7)][0];
163 xj0[1] = x[KEL(0,fa7)][1]-x[KEL(1,fa7)][1];
164 xj0[2] = x[KEL(0,fa7)][2]-x[KEL(1,fa7)][2];
165 double psc = xj0[0]*pv[0]+xj0[1]*pv[1]+xj0[2]*pv[2];
166 double sign = (psc < 0) ? -0.5 : 0.5;
167 facette_normale(i,fa7,0) = sign*pv[0];
168 facette_normale(i,fa7,1) = sign*pv[1];
169 facette_normale(i,fa7,2) = sign*pv[2];
170 }
171 }
172 });
173 end_gpu_timer(__KERNEL_NAME__);
174}
175
176
177/*! @brief Fills the normales_facettes_Cl array in the Domaine_Cl_VEF for sub-facet fa7 of element num_elem.
178 *
179 */
180void Tetra_VEF::creer_normales_facettes_Cl(DoubleTab& normales_facettes_Cl,
181 int fa7,
182 int num_elem,const DoubleTab& x,
183 const DoubleVect& xg, const Domaine& domaine_geom) const
184{
185 // x contains the coordinates of the tetrahedron vertices
186 // xg contains the coordinates of the tetrahedron "center"
187 double u[3];
188 double v[3];
189 double xj0[3];
190 double psc;
191 // Cerr << "fa7 " << fa7 << " et le num_elem " << num_elem << finl;
192 int i0 = KEL_(0,fa7);
193 int i1 = KEL_(1,fa7);
194 int i2 = KEL_(2,fa7);
195 int i3 = KEL_(3,fa7);
196 // Cerr << "i0 " << i0 << " i1 " << i1 << " i2 " << i2 << " i3 " << i3 << finl;
197 u[0]=x(i2,0)-xg[0];
198 u[1]=x(i2,1)-xg[1];
199 u[2]=x(i2,2)-xg[2];
200 v[0]=x(i3,0)-xg[0];
201 v[1]=x(i3,1)-xg[1];
202 v[2]=x(i3,2)-xg[2];
203 double pv[3];
204 prodvect(u,v,pv);
205
206 // Normal orientation:
207 xj0[0]= x(i0,0)-x(i1,0);
208 xj0[1]= x(i0,1)-x(i1,1);
209 xj0[2]= x(i0,2)-x(i1,2);
210
211 psc=xj0[0]*pv[0] + xj0[1]*pv[1] + xj0[2]*pv[2] ;
212
213 if (psc < 0)
214 {
215 normales_facettes_Cl(num_elem,fa7,0) = -pv[0]/2;
216 normales_facettes_Cl(num_elem,fa7,1) = -pv[1]/2;
217 normales_facettes_Cl(num_elem,fa7,2) = -pv[2]/2;
218 }
219 else
220 {
221 normales_facettes_Cl(num_elem,fa7,0) = pv[0]/2;
222 normales_facettes_Cl(num_elem,fa7,1) = pv[1]/2;
223 normales_facettes_Cl(num_elem,fa7,2) = pv[2]/2;
224 }
225}
226
227
229 const Domaine_VEF& le_dom_VEF,
230 DoubleVect& volumes_entrelaces_Cl,
231 int type_cl) const
232{
233 double vol_mod;
234 const DoubleVect& volumes_entrelaces = le_dom_VEF.volumes_entrelaces();
235 const IntTab& elem_faces = le_dom_VEF.elem_faces();
236
237 switch(type_cl)
238 {
239
240 // no Dirichlet face: impossible
241 case 0:
242 {
243 Cerr << "Tetra_VEF::modif_volumes_entrelaces() type 0 not possible!\n";
244 break;
245 }
246
247 case 1: // one Dirichlet face: Face 3
248 {
249 vol_mod = volumes_entrelaces[j]/3 ;
250 volumes_entrelaces_Cl[elem_faces(elem,0)] += vol_mod;
251 volumes_entrelaces_Cl[elem_faces(elem,1)] += vol_mod;
252 volumes_entrelaces_Cl[elem_faces(elem,2)] += vol_mod;
253 break;
254 }
255
256 case 2: // one Dirichlet face: Face 2
257 {
258 vol_mod = volumes_entrelaces[j]/3 ;
259 volumes_entrelaces_Cl[elem_faces(elem,0)] += vol_mod;
260 volumes_entrelaces_Cl[elem_faces(elem,1)] += vol_mod;
261 volumes_entrelaces_Cl[elem_faces(elem,3)] += vol_mod;
262 break;
263 }
264
265 case 4: // one Dirichlet face: Face 1
266 {
267 vol_mod = volumes_entrelaces[j]/3 ;
268 volumes_entrelaces_Cl[elem_faces(elem,0)] += vol_mod;
269 volumes_entrelaces_Cl[elem_faces(elem,2)] += vol_mod;
270 volumes_entrelaces_Cl[elem_faces(elem,3)] += vol_mod;
271 break;
272 }
273
274 case 8: // one Dirichlet face: Face 0
275 {
276 vol_mod = volumes_entrelaces[j]/3 ;
277 volumes_entrelaces_Cl[elem_faces(elem,1)] += vol_mod;
278 volumes_entrelaces_Cl[elem_faces(elem,2)] += vol_mod;
279 volumes_entrelaces_Cl[elem_faces(elem,3)] += vol_mod;
280 break;
281 }
282
283 case 3: // two Dirichlet faces: faces 2,3
284 {
285 vol_mod = (volumes_entrelaces[elem_faces(elem,2)]
286 + volumes_entrelaces[elem_faces(elem,3)])/2;
287 volumes_entrelaces_Cl[elem_faces(elem,0)] += vol_mod;
288 volumes_entrelaces_Cl[elem_faces(elem,1)] += vol_mod;
289 break;
290 }
291
292 case 5: // two Dirichlet faces: faces 1,3
293 {
294 vol_mod = (volumes_entrelaces[elem_faces(elem,1)]
295 + volumes_entrelaces[elem_faces(elem,3)])/2;
296 volumes_entrelaces_Cl[elem_faces(elem,0)] += vol_mod;
297 volumes_entrelaces_Cl[elem_faces(elem,2)] += vol_mod;
298 break;
299 }
300
301 case 6: // the tetrahedron has two Dirichlet faces: faces 1,2
302 {
303 vol_mod = (volumes_entrelaces[elem_faces(elem,1)]
304 + volumes_entrelaces[elem_faces(elem,2)])/2;
305 volumes_entrelaces_Cl[elem_faces(elem,0)] += vol_mod;
306 volumes_entrelaces_Cl[elem_faces(elem,3)] += vol_mod;
307 break;
308 }
309
310 case 9: // two Dirichlet faces: faces 0,3
311 {
312 vol_mod = (volumes_entrelaces[elem_faces(elem,0)]
313 + volumes_entrelaces[elem_faces(elem,3)])/2;
314 volumes_entrelaces_Cl[elem_faces(elem,1)] += vol_mod;
315 volumes_entrelaces_Cl[elem_faces(elem,2)] += vol_mod;
316 break;
317 }
318
319 case 10: // two Dirichlet faces: faces 0,2
320 {
321 vol_mod = (volumes_entrelaces[elem_faces(elem,0)]
322 + volumes_entrelaces[elem_faces(elem,2)])/2;
323 volumes_entrelaces_Cl[elem_faces(elem,1)] += vol_mod;
324 volumes_entrelaces_Cl[elem_faces(elem,3)] += vol_mod;
325 break;
326 }
327
328 case 12: // two Dirichlet faces: faces 0,1
329 {
330 vol_mod = (volumes_entrelaces[elem_faces(elem,0)]
331 + volumes_entrelaces[elem_faces(elem,1)])/2;
332 volumes_entrelaces_Cl[elem_faces(elem,2)] += vol_mod;
333 volumes_entrelaces_Cl[elem_faces(elem,3)] += vol_mod;
334 break;
335 }
336
337
338
339 case 7: // three Dirichlet faces: faces 1,2,3
340 {
341 vol_mod = volumes_entrelaces[elem_faces(elem,1)]
342 + volumes_entrelaces[elem_faces(elem,2)]
343 + volumes_entrelaces[elem_faces(elem,3)];
344 volumes_entrelaces_Cl[elem_faces(elem,0)] += vol_mod;
345 break;
346 }
347
348 case 11: // three Dirichlet faces: faces 0,2,3
349 {
350 vol_mod = volumes_entrelaces[elem_faces(elem,0)]
351 + volumes_entrelaces[elem_faces(elem,2)]
352 + volumes_entrelaces[elem_faces(elem,3)];
353 volumes_entrelaces_Cl[elem_faces(elem,1)] += vol_mod;
354 break;
355 }
356
357 case 13: // three Dirichlet faces: faces 0,1,3
358 {
359
360 vol_mod = volumes_entrelaces[elem_faces(elem,0)]
361 + volumes_entrelaces[elem_faces(elem,1)]
362 + volumes_entrelaces[elem_faces(elem,3)];
363 volumes_entrelaces_Cl[elem_faces(elem,2)] += vol_mod;
364 break;
365 }
366
367 case 14: // three Dirichlet faces: faces 0,1,2
368 {
369 vol_mod = volumes_entrelaces[elem_faces(elem,0)]
370 + volumes_entrelaces[elem_faces(elem,1)]
371 + volumes_entrelaces[elem_faces(elem,2)];
372 volumes_entrelaces_Cl[elem_faces(elem,3)] += vol_mod;
373 break;
374 }
375 default :
376 {
377 Cerr << "\n unknown type in Tetra_VEF::modif_volumes_entrelaces: " << type_cl;
378 exit();
379 }
380
381 } // end of switch
382}
383
385 const Domaine_VEF& le_dom_VEF,
386 DoubleVect& volumes_entrelaces_Cl,
387 int type_cl) const
388{
389 double vol_mod;
390 const DoubleVect& volumes_entrelaces = le_dom_VEF.volumes_entrelaces();
391 const IntTab& elem_faces = le_dom_VEF.elem_faces();
392
393 int face;
394 int nb_faces_cl = volumes_entrelaces_Cl.size();
395 switch(type_cl)
396 {
397
398 // no Dirichlet face: impossible
399 case 0:
400 {
401 Cerr << "Tetra_VEF::modif_volumes_entrelaces() type 0 not possible!\n";
402 break;
403 }
404
405 case 1: // one Dirichlet face: Face 3
406 {
407 vol_mod = volumes_entrelaces[j]/3 ;
408 face=elem_faces(elem,0);
409 if(face<nb_faces_cl)
410 volumes_entrelaces_Cl[face] += vol_mod;
411 face=elem_faces(elem,1);
412 if(face<nb_faces_cl)
413 volumes_entrelaces_Cl[face] += vol_mod;
414 face=elem_faces(elem,2);
415 if(face<nb_faces_cl)
416 volumes_entrelaces_Cl[face] += vol_mod;
417 break;
418 }
419
420 case 2: // one Dirichlet face: Face 2
421 {
422 vol_mod = volumes_entrelaces[j]/3 ;
423 face=elem_faces(elem,0);
424 if(face<nb_faces_cl)
425 volumes_entrelaces_Cl[face] += vol_mod;
426 face=elem_faces(elem,1);
427 if(face<nb_faces_cl)
428 volumes_entrelaces_Cl[face] += vol_mod;
429 face=elem_faces(elem,3);
430 if(face<nb_faces_cl)
431 volumes_entrelaces_Cl[face] += vol_mod;
432 break;
433 }
434
435 case 4: // one Dirichlet face: Face 1
436 {
437 vol_mod = volumes_entrelaces[j]/3 ;
438 face=elem_faces(elem,0);
439 if(face<nb_faces_cl)
440 volumes_entrelaces_Cl[face] += vol_mod;
441 face=elem_faces(elem,2);
442 if(face<nb_faces_cl)
443 volumes_entrelaces_Cl[face] += vol_mod;
444 face=elem_faces(elem,3);
445 if(face<nb_faces_cl)
446 volumes_entrelaces_Cl[face] += vol_mod;
447 break;
448 }
449
450 case 8: // one Dirichlet face: Face 0
451 {
452 vol_mod = volumes_entrelaces[j]/3 ;
453 face=elem_faces(elem,1);
454 if(face<nb_faces_cl)
455 volumes_entrelaces_Cl[face] += vol_mod;
456 face=elem_faces(elem,2);
457 if(face<nb_faces_cl)
458 volumes_entrelaces_Cl[face] += vol_mod;
459 face=elem_faces(elem,3);
460 if(face<nb_faces_cl)
461 volumes_entrelaces_Cl[face] += vol_mod;
462 break;
463 }
464
465 case 3: // two Dirichlet faces: faces 2,3
466 {
467 vol_mod = (volumes_entrelaces[elem_faces(elem,2)]
468 + volumes_entrelaces[elem_faces(elem,3)])/2;
469 face=elem_faces(elem,0);
470 if(face<nb_faces_cl)
471 volumes_entrelaces_Cl[face] += vol_mod;
472 face=elem_faces(elem,1);
473 if(face<nb_faces_cl)
474 volumes_entrelaces_Cl[face] += vol_mod;
475 break;
476 }
477
478 case 5: // two Dirichlet faces: faces 1,3
479 {
480 vol_mod = (volumes_entrelaces[elem_faces(elem,1)]
481 + volumes_entrelaces[elem_faces(elem,3)])/2;
482 face=elem_faces(elem,0);
483 if(face<nb_faces_cl)
484 volumes_entrelaces_Cl[face] += vol_mod;
485 face=elem_faces(elem,2);
486 if(face<nb_faces_cl)
487 volumes_entrelaces_Cl[face] += vol_mod;
488 break;
489 }
490
491 case 6: // the tetrahedron has two Dirichlet faces: faces 1,2
492 {
493 vol_mod = (volumes_entrelaces[elem_faces(elem,1)]
494 + volumes_entrelaces[elem_faces(elem,2)])/2;
495 face=elem_faces(elem,0);
496 if(face<nb_faces_cl)
497 volumes_entrelaces_Cl[face] += vol_mod;
498 face=elem_faces(elem,3);
499 if(face<nb_faces_cl)
500 volumes_entrelaces_Cl[face] += vol_mod;
501 break;
502 }
503
504 case 9: // two Dirichlet faces: faces 0,3
505 {
506 vol_mod = (volumes_entrelaces[elem_faces(elem,0)]
507 + volumes_entrelaces[elem_faces(elem,3)])/2;
508 face=elem_faces(elem,1);
509 if(face<nb_faces_cl)
510 volumes_entrelaces_Cl[face] += vol_mod;
511 face=elem_faces(elem,2);
512 if(face<nb_faces_cl)
513 volumes_entrelaces_Cl[face] += vol_mod;
514 break;
515 }
516
517 case 10: // two Dirichlet faces: faces 0,2
518 {
519 vol_mod = (volumes_entrelaces[elem_faces(elem,0)]
520 + volumes_entrelaces[elem_faces(elem,2)])/2;
521 face=elem_faces(elem,1);
522 if(face<nb_faces_cl)
523 volumes_entrelaces_Cl[face] += vol_mod;
524 face=elem_faces(elem,3);
525 if(face<nb_faces_cl)
526 volumes_entrelaces_Cl[face] += vol_mod;
527 break;
528 }
529
530 case 12: // two Dirichlet faces: faces 0,1
531 {
532 vol_mod = (volumes_entrelaces[elem_faces(elem,0)]
533 + volumes_entrelaces[elem_faces(elem,1)])/2;
534 face=elem_faces(elem,2);
535 if(face<nb_faces_cl)
536 volumes_entrelaces_Cl[face] += vol_mod;
537 face=elem_faces(elem,3);
538 if(face<nb_faces_cl)
539 volumes_entrelaces_Cl[face] += vol_mod;
540 break;
541 }
542
543
544
545 case 7: // three Dirichlet faces: faces 1,2,3
546 {
547 vol_mod = volumes_entrelaces[elem_faces(elem,1)]
548 + volumes_entrelaces[elem_faces(elem,2)]
549 + volumes_entrelaces[elem_faces(elem,3)];
550 face=elem_faces(elem,0);
551 if(face<nb_faces_cl)
552 volumes_entrelaces_Cl[face] += vol_mod;
553 break;
554 }
555
556 case 11: // three Dirichlet faces: faces 0,2,3
557 {
558 vol_mod = volumes_entrelaces[elem_faces(elem,0)]
559 + volumes_entrelaces[elem_faces(elem,2)]
560 + volumes_entrelaces[elem_faces(elem,3)];
561 face=elem_faces(elem,1);
562 if(face<nb_faces_cl)
563 volumes_entrelaces_Cl[face] += vol_mod;
564 break;
565 }
566
567 case 13: // three Dirichlet faces: faces 0,1,3
568 {
569
570 vol_mod = volumes_entrelaces[elem_faces(elem,0)]
571 + volumes_entrelaces[elem_faces(elem,1)]
572 + volumes_entrelaces[elem_faces(elem,3)];
573 face=elem_faces(elem,2);
574 if(face<nb_faces_cl)
575 volumes_entrelaces_Cl[face] += vol_mod;
576 break;
577 }
578
579 case 14: // three Dirichlet faces: faces 0,1,2
580 {
581 vol_mod = volumes_entrelaces[elem_faces(elem,0)]
582 + volumes_entrelaces[elem_faces(elem,1)]
583 + volumes_entrelaces[elem_faces(elem,2)];
584 face=elem_faces(elem,3);
585 if(face<nb_faces_cl)
586 volumes_entrelaces_Cl[face] += vol_mod;
587 break;
588 }
589 default :
590 {
591 Cerr << "\n unknown type in Tetra_VEF::modif_volumes_entrelaces: " << type_cl;
592 exit();
593 }
594
595 } // end of switch
596}
597
598void Tetra_VEF::calcul_vc(const ArrOfInt& Face,ArrOfDouble& vc,
599 const ArrOfDouble& vs,const DoubleTab& vsom,
600 const Champ_Inc_base& vitesse,int type_cl, const DoubleVect& porosite_face) const
601{
602 DoubleTab poro(4);
603 DoubleTab vfa(4,3);
604 for (int i=0; i<4; i++)
605 for (int j=0; j<3; j++)
606 {
607 vfa(i,j) = vitesse.valeurs()(Face[i],j);
608 }
609
610 for (int i=0; i<4; i++)
611 {
612 poro(i) = porosite_face(Face[i]);
613 }
614
615 calcul_vc_tetra(Face.addr(), vc.addr(), vs.addr(), vsom.addr(), vfa.addr(), (int)type_cl, poro.addr());
616}
617
618/*! @brief Computes the coordinates xg of the center of a non-standard element.
619 *
620 * Also computes idirichlet = number of Dirichlet faces of the element.
621 *
622 */
623void Tetra_VEF::calcul_xg(DoubleVect& xg,const DoubleTab& x, const int type_elem_Cl,
624 int& idirichlet,int& n1,int& n2,int& n3) const
625{
626 calcul_xg_tetra(xg.addr(), x.addr(), (int)type_elem_Cl, idirichlet, n1, n2, n3);
627}
628
629void Tetra_VEF::modif_normales_facettes_Cl(DoubleTab& normales_facettes_Cl,
630 int fa7,int num_elem,int idirichlet,
631 int n1,int n2,int n3) const
632{
633 switch (idirichlet)
634 {
635
636 case 0:
637 break;
638
639 case 1:
640 break;
641
642 case 2: // one null sub-facet n1
643 {
644 fa7=n1;
645 normales_facettes_Cl(num_elem,fa7,0) = 0;
646 normales_facettes_Cl(num_elem,fa7,1) = 0;
647 normales_facettes_Cl(num_elem,fa7,2) = 0;
648 break;
649 }
650
651 case 3:
652 {
653 fa7=n1;
654 normales_facettes_Cl(num_elem,fa7,0) = 0;
655 normales_facettes_Cl(num_elem,fa7,1) = 0;
656 normales_facettes_Cl(num_elem,fa7,2) = 0;
657 fa7=n2;
658 normales_facettes_Cl(num_elem,fa7,0) = 0;
659 normales_facettes_Cl(num_elem,fa7,1) = 0;
660 normales_facettes_Cl(num_elem,fa7,2) = 0;
661 fa7=n3;
662 normales_facettes_Cl(num_elem,fa7,0) = 0;
663 normales_facettes_Cl(num_elem,fa7,1) = 0;
664 normales_facettes_Cl(num_elem,fa7,2) = 0;
665 break;
666 }
667
668 } // end of switch
669}
670
671
672
Class Champ_Inc_base.
DoubleTab & valeurs() override
Returns the array of field values at the current time.
int_t nb_elem_tot() const
Definition Domaine.h:132
IntTab_t & les_elems()
Definition Domaine.h:129
const DoubleTab_t & coord_sommets() const
Definition Domaine.h:112
class Domaine_VEF
Definition Domaine_VEF.h:53
DoubleVect & volumes_entrelaces()
Definition Domaine_VF.h:99
int elem_faces(int i, int j) const
Returns the index of the i-th face of element num_elem; the face numbering convention is.
Definition Domaine_VF.h:542
const Domaine & domaine() const
const IntTab & KEL() const
Class defining operators and methods for all reading operation in an input flow (file,...
Definition Entree.h:42
const Nom & que_suis_je() const
Returns the string identifying the class.
Definition Objet_U.cpp:104
virtual Entree & readOn(Entree &)
Reads an Objet_U from an input stream. Virtual method to override.
Definition Objet_U.cpp:289
virtual Sortie & printOn(Sortie &) const
Writes the object to an output stream. Virtual method to override.
Definition Objet_U.cpp:278
static void exit(int exit_code=-1)
Exit routine for TRUST within a Kokkos region.
Definition Process.cpp:466
Base class for output streams.
Definition Sortie.h:52
_TYPE_ * addr()
_SIZE_ dimension_tot(int) const override
Definition TRUSTTab.tpp:160
std::enable_if_t< is_default_exec_space< EXEC_SPACE >, ConstView< _TYPE_, _SHAPE_ > > view_ro() const
Definition TRUSTTab.h:261
std::enable_if_t< is_default_exec_space< EXEC_SPACE >, View< _TYPE_, _SHAPE_ > > view_rw()
Definition TRUSTTab.h:291
_SIZE_ size() const
Definition TRUSTVect.tpp:45
void modif_volumes_entrelaces_faces_joints(int, int, const Domaine_VEF &, DoubleVect &, int) const override
void modif_volumes_entrelaces(int, int, const Domaine_VEF &, DoubleVect &, int) const override
void creer_normales_facettes_Cl(DoubleTab &, int, int, const DoubleTab &, const DoubleVect &, const Domaine &) const override
Fills the normales_facettes_Cl array in the Domaine_Cl_VEF for sub-facet fa7 of element num_elem.
void calcul_vc(const ArrOfInt &, ArrOfDouble &, const ArrOfDouble &, const DoubleTab &, const Champ_Inc_base &, int, const DoubleVect &) const override
Tetra_VEF()
Returns for sub-facet fa7: for j=0,j=1: the local indices of the 2 faces surrounding fa7.
Definition Tetra_VEF.cpp:40
void modif_normales_facettes_Cl(DoubleTab &, int, int, int, int, int, int) const override
void calcul_xg(DoubleVect &, const DoubleTab &, const int, int &, int &, int &, int &) const override
Computes the coordinates xg of the center of a non-standard element.
void creer_facette_normales(const Domaine_VEF &, const IntVect &) const override
Fills the face_normales array in the Domaine_VEF.
void creer_face_normales(DoubleTab &, const IntTab &, const IntTab &, const IntTab &, const Domaine &) const override
Definition Tetra_VEF.cpp:55