TrioCFD 1.9.8
TrioCFD documentation
Loading...
Searching...
No Matches
Tetraedriser_homogene_fin.cpp
1/****************************************************************************
2* Copyright (c) 2023, 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 <Tetraedriser_homogene_fin.h>
17#include <Domaine.h>
18
19Implemente_instanciable(Tetraedriser_homogene_fin, "Tetraedriser_homogene_fin", Triangulation_base);
20// XD tetraedriser_homogene_fin tetraedriser tetraedriser_homogene_fin INHERITS_BRACE Tetraedriser_homogene_fin is the
21// XD_CONT recommended option to tetrahedralise blocks. As an extension (subdivision) of Tetraedriser_homogene_compact,
22// XD_CONT this last one cut each initial block in 48 tetrahedra (against 24, previously). This cutting ensures : NL2 -
23// XD_CONT a correct cutting in the corners (in respect to pressure discretization PreP1B), NL2 - a better isotropy of
24// XD_CONT elements than with Tetraedriser_homogene_compact, NL2 - a better alignment of summits (this could have a
25// XD_CONT benefit effect on calculation near walls since first elements in contact with it are all contained in the
26// XD_CONT same constant thickness and ii/ by the way, a 3D cartesian grid based on summits can be engendered and used
27// XD_CONT to realise spectral analysis in HIT for instance). Initial block is divided in 48 tetrahedra:
28// XD_CONT \includeimage{{tetraedriserhomogenefin.jpeg}}
29
30
31Sortie& Tetraedriser_homogene_fin::printOn(Sortie& os) const { return Interprete::printOn(os); }
32
33Entree& Tetraedriser_homogene_fin::readOn(Entree& is) { return Interprete::readOn(is); }
34
35// Traitement des faces
36//
37void Tetraedriser_homogene_fin::decoupe(Domaine& domaine, Faces& faces, IntTab& new_soms_old_elems, IntTab& fait_sommet, int nface, IntTab& fait_sommet_arete, int narete) const
38{
39 IntTab& sommets = faces.les_sommets();
40 int nb_faces = sommets.dimension(0);
41 assert(sommets.dimension(1) == 4);
42 IntTab nouveaux(8 * nb_faces, 3);
43 faces.voisins().resize(8 * nb_faces, 2);
44 faces.voisins() = -1;
45 IntVect indice(5);
46
47 for (int i = 0; i < nb_faces; i++)
48 {
49 indice(0) = sommets(i, 0);
50 indice(1) = sommets(i, 1);
51 indice(2) = sommets(i, 2);
52 indice(3) = sommets(i, 3);
53
54 // on remet un peu d'ordre dans les indices
55
56 int permut, ind;
57
58 do
59 {
60 permut = 0;
61 for (int j = 1; j < 4; j++)
62 {
63 if (indice[j] < indice[j - 1])
64 {
65 ind = indice[j - 1];
66 indice[j - 1] = indice[j];
67 indice[j] = ind;
68 permut++;
69 }
70 }
71 }
72 while (permut);
73
74 int i0 = indice(0);
75 int i1 = indice(1);
76 int i2 = indice(2);
77 int i3 = indice(3);
78 int i4 = -1;
79
80 for (int ii = 0; ii < nface; ii++)
81 {
82 if (fait_sommet(ii, 0) == i0 && fait_sommet(ii, 1) == i1 && fait_sommet(ii, 2) == i2)
83 {
84 i4 = fait_sommet(ii, 3);
85 break;
86 }
87 }
88
89 if (i4 < 0)
90 {
91 Cerr << " FALSE " << finl;
92 exit();
93 }
94 assert(i4 >= 0);
95
96 int jj = 0;
97
98 for (int ii = 0; ii < narete; ii++)
99 {
100 if (fait_sommet_arete(ii, 0) == i0 && fait_sommet_arete(ii, 1) == i1)
101 {
102 nouveaux(jj * nb_faces + i, 0) = i0;
103 nouveaux(jj * nb_faces + i, 1) = fait_sommet_arete(ii, 2);
104 nouveaux(jj * nb_faces + i, 2) = i4;
105 jj++;
106
107 nouveaux(jj * nb_faces + i, 0) = i1;
108 nouveaux(jj * nb_faces + i, 1) = fait_sommet_arete(ii, 2);
109 nouveaux(jj * nb_faces + i, 2) = i4;
110 jj++;
111 }
112
113 if (fait_sommet_arete(ii, 0) == i0 && fait_sommet_arete(ii, 1) == i2)
114 {
115 nouveaux(jj * nb_faces + i, 0) = i0;
116 nouveaux(jj * nb_faces + i, 1) = fait_sommet_arete(ii, 2);
117 nouveaux(jj * nb_faces + i, 2) = i4;
118 jj++;
119
120 nouveaux(jj * nb_faces + i, 0) = i2;
121 nouveaux(jj * nb_faces + i, 1) = fait_sommet_arete(ii, 2);
122 nouveaux(jj * nb_faces + i, 2) = i4;
123 jj++;
124 }
125
126 if (fait_sommet_arete(ii, 0) == i0 && fait_sommet_arete(ii, 1) == i3)
127 {
128 nouveaux(jj * nb_faces + i, 0) = i0;
129 nouveaux(jj * nb_faces + i, 1) = fait_sommet_arete(ii, 2);
130 nouveaux(jj * nb_faces + i, 2) = i4;
131 jj++;
132
133 nouveaux(jj * nb_faces + i, 0) = i3;
134 nouveaux(jj * nb_faces + i, 1) = fait_sommet_arete(ii, 2);
135 nouveaux(jj * nb_faces + i, 2) = i4;
136 jj++;
137 }
138
139 if (fait_sommet_arete(ii, 0) == i1 && fait_sommet_arete(ii, 1) == i2)
140 {
141 nouveaux(jj * nb_faces + i, 0) = i1;
142 nouveaux(jj * nb_faces + i, 1) = fait_sommet_arete(ii, 2);
143 nouveaux(jj * nb_faces + i, 2) = i4;
144 jj++;
145
146 nouveaux(jj * nb_faces + i, 0) = i2;
147 nouveaux(jj * nb_faces + i, 1) = fait_sommet_arete(ii, 2);
148 nouveaux(jj * nb_faces + i, 2) = i4;
149 jj++;
150 }
151
152 if (fait_sommet_arete(ii, 0) == i1 && fait_sommet_arete(ii, 1) == i3)
153 {
154 nouveaux(jj * nb_faces + i, 0) = i1;
155 nouveaux(jj * nb_faces + i, 1) = fait_sommet_arete(ii, 2);
156 nouveaux(jj * nb_faces + i, 2) = i4;
157 jj++;
158
159 nouveaux(jj * nb_faces + i, 0) = i3;
160 nouveaux(jj * nb_faces + i, 1) = fait_sommet_arete(ii, 2);
161 nouveaux(jj * nb_faces + i, 2) = i4;
162 jj++;
163 }
164
165 if (fait_sommet_arete(ii, 0) == i2 && fait_sommet_arete(ii, 1) == i3)
166 {
167 nouveaux(jj * nb_faces + i, 0) = i2;
168 nouveaux(jj * nb_faces + i, 1) = fait_sommet_arete(ii, 2);
169 nouveaux(jj * nb_faces + i, 2) = i4;
170 jj++;
171
172 nouveaux(jj * nb_faces + i, 0) = i3;
173 nouveaux(jj * nb_faces + i, 1) = fait_sommet_arete(ii, 2);
174 nouveaux(jj * nb_faces + i, 2) = i4;
175 jj++;
176 }
177 }
178
179 // if(jj!=8)
180 // {
181 // Cerr << " il manque des nouvelles faces " << finl;
182 // exit();
183 // }
184 assert(jj == 8);
185
186 }
187 sommets.ref(nouveaux);
188}
189
190int Tetraedriser_homogene_fin::creer_sommet(Domaine& domaine, DoubleTab& new_soms, IntTab& elem_traite, IntTab& new_soms_old_elems, int NbSom, IntTab& sommets, int& compteur, int oldnbsom,
191 int& nbnewsoms, IntTab& fait_sommet, int& nface, IntTab& fait_sommet_arete, int& narete) const
192{
193
194 const DoubleTab& coord_sommets = domaine.coord_sommets();
195 //int nb_domaines = dom.nb_domaines();
196 int _out;
197
198 if (NbSom == 4) // On verifie que le sommet rattache a une face n'a pas ete cree ulterieurement
199 {
200 int permut, som;
201
202 do
203 {
204 permut = 0;
205 for (int j = 1; j < 4; j++)
206 {
207 if (sommets[j] < sommets[j - 1])
208 {
209 som = sommets[j - 1];
210 sommets[j - 1] = sommets[j];
211 sommets[j] = som;
212 permut++;
213 }
214 }
215 }
216 while (permut);
217
218 }
219
220 if (NbSom == 2) // On verifie que le sommet rattache a une arete n'a pas ete cree ulterieurement
221 {
222 int permut, som;
223
224 do
225 {
226 permut = 0;
227 for (int j = 1; j < 2; j++)
228 {
229 if (sommets[j] < sommets[j - 1])
230 {
231 som = sommets[j - 1];
232 sommets[j - 1] = sommets[j];
233 sommets[j] = som;
234 permut++;
235 }
236 }
237 }
238 while (permut);
239
240 }
241
242 double x = 0, y = 0, z = 0;
243 //x = 0.5*(coord_sommets(i1,0)+coord_sommets(i2,0));
244 //y = 0.5*(coord_sommets(i1,1)+coord_sommets(i2,1));
245 //z = 0.5*(coord_sommets(i1,2)+coord_sommets(i2,2));
246 for (int ii = 0; ii < NbSom; ii++)
247 {
248 x += coord_sommets(sommets(ii), 0);
249 y += coord_sommets(sommets(ii), 1);
250 z += coord_sommets(sommets(ii), 2);
251 }
252 x /= NbSom;
253 y /= NbSom;
254 z /= NbSom;
255
256 int trouve = -1;
257
258 trouve = domaine.chercher_elements(x, y, z);
259
260 if (trouve < 0)
261 {
262 for (int ii = 0; ii < NbSom; ii++)
263 {
264 Cerr << "i" << ii << "=" << sommets(ii) << " coord=" << coord_sommets(sommets(ii), 0) << " " << coord_sommets(sommets(ii), 1) << " " << coord_sommets(sommets(ii), 2) << finl;
265 }
266 Cerr << "x=" << x << " y=" << y << " z=" << z << finl;
267 Cerr << "-----" << finl;
268 assert(trouve >= 0);
269 }
270
271 if (elem_traite(trouve))
272 {
273 _out = -1;
274 int j = 0;
275 double epsilon = Objet_U::precision_geom;
276 do
277 {
278
279 if ((std::fabs(coord_sommets(new_soms_old_elems(trouve, j), 0) - x) < epsilon) && (std::fabs(coord_sommets(new_soms_old_elems(trouve, j), 1) - y) < epsilon)
280 && (std::fabs(coord_sommets(new_soms_old_elems(trouve, j), 2) - z) < epsilon))
281
282 {
283 _out = new_soms_old_elems(trouve, j);
284 }
285 j++;
286 }
287 while ((_out == -1) && (j < 19));
288 if (_out < 0)
289 {
290 Cerr << "Error, not found !" << finl;
291 Cerr << "found=" << trouve << " x=" << x << " y=" << y << " z=" << z << finl;
292 for (int ii = 0; ii < NbSom; ii++)
293 {
294 Cerr << "i" << ii << "=" << sommets(ii) << " coord=" << coord_sommets(sommets(ii), 0) << " " << coord_sommets(sommets(ii), 1) << " " << coord_sommets(sommets(ii), 2) << finl;
295 }
296 exit();
297 }
298 }
299 else
300 {
301
302 _out = oldnbsom + nbnewsoms;
303 new_soms(compteur, 0) = x;
304 new_soms(compteur, 1) = y;
305 new_soms(compteur, 2) = z;
306 nbnewsoms++;
307 compteur++;
308
309 if (NbSom == 4)
310 {
311 fait_sommet(nface, 0) = sommets[0];
312 fait_sommet(nface, 1) = sommets[1];
313 fait_sommet(nface, 2) = sommets[2];
314 fait_sommet(nface, 3) = _out;
315 // Cerr << sommets(0) << " " << sommets(1) << " " << sommets(2) << " " << _out << finl;
316 nface++;
317 }
318
319 if (NbSom == 2)
320 {
321 fait_sommet_arete(narete, 0) = sommets[0];
322 fait_sommet_arete(narete, 1) = sommets[1];
323 fait_sommet_arete(narete, 2) = _out;
324 // Cerr << sommets(0) << " " << sommets(1) << " " << _out << finl;
325 narete++;
326 }
327 }
328
329 if (_out < 0)
330 {
331 Cerr << "Problem at the creation of node" << finl;
332 exit();
333 }
334 return _out;
335
336}
337
338void Tetraedriser_homogene_fin::trianguler(Domaine& domaine) const
339{
340 if ((domaine.type_elem()->que_suis_je() == "Hexaedre") || (domaine.type_elem()->que_suis_je() == "Hexaedre_VEF"))
341 {
342 IntTab& les_elems = domaine.les_elems();
343 int oldsz = les_elems.dimension(0);
344 int nbs = domaine.nb_som();
345 IntTab elem_traite(oldsz);
346 int oldnbsom = domaine.nb_som();
347 IntTab new_elems(48 * oldsz, 4);
348 // pour chaque cube, liste des nouveaux sommets qu'il contient :
349 IntTab new_soms_old_elems(oldsz, 19);
350 IntTab sommets(8);
351 int compteur = 0;
352 int nbnewsoms = 0;
353 int nface = 0;
354 int narete = 0;
355 int i;
356
357 Cerr << " NB ELEM : " << oldsz << " NB NODE : " << nbs << finl;
358 IntTab fait_sommet(3 * nbs, 4); // PQ : 04-03 : dimensionnement de fait_sommet a 3*nbs (a revoir si insuffisant)
359 IntTab fait_sommet_arete(3 * nbs, 3); // PQ : 02-05 : dimensionnement de fait_sommet_arete a 3*nbs (a revoir si insuffisant)
360 fait_sommet = -1;
361 fait_sommet_arete = -1;
362
363 // Construction de l'Octree sur la grille "VDF" de base
364 domaine.construit_octree();
365
366 //On dimensionne une premiere fois le tableau des sommets avec la dimension maximum
367 //puis on redimensionnera seulement a la fin par la dimension exacte
368
369 DoubleTab& sommets_dom = domaine.les_sommets();
370 //int dim_som_max=19*oldsz*40;
371 //19 pour les nouveaux sommets et 8 pour les anciens sommets = 27
372 int dim_som_max = 27 * oldsz;
373 int dim_som_old = sommets_dom.dimension(0);
374 sommets_dom.resize(dim_som_max, 3);
375
376 IntTab indice(19);
377 DoubleTab new_soms(19, 3);
378
379 for (i = 0; i < oldsz; i++)
380 {
381 if (nface >= 3 * nbs)
382 {
383 Cerr << " The sizing of the array fait_sommet in Tetraedriser_homogene_fin::Trianguler is inappropriate " << finl;
384 exit();
385 }
386 if (narete >= 3 * nbs)
387 {
388 Cerr << " The sizing of the array fait_sommet_arete in Tetraedriser_homogene_fin::Trianguler is inappropriate " << finl;
389 exit();
390 }
391
392 int i0 = les_elems(i, 0);
393 int i1 = les_elems(i, 1);
394 int i2 = les_elems(i, 2);
395 int i3 = les_elems(i, 3);
396 int i4 = les_elems(i, 4);
397 int i5 = les_elems(i, 5);
398 int i6 = les_elems(i, 6);
399 int i7 = les_elems(i, 7);
400
401 compteur = 0;
402
403 // Definition des nouveaux sommets : creation des barycentres
404 //centre de l'hexaedre
405
406 sommets(0) = i0;
407 sommets(1) = i1;
408 sommets(2) = i2;
409 sommets(3) = i3;
410 sommets(4) = i4;
411 sommets(5) = i5;
412 sommets(6) = i6;
413 sommets(7) = i7;
414 indice(0) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 8, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
415
416 //centres des faces bas (0-1-3-2)
417 sommets(0) = i0;
418 sommets(1) = i1;
419 sommets(2) = i3;
420 sommets(3) = i2;
421 indice(1) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 4, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
422
423 //centres des faces haut (4-5-7-6)
424 sommets(0) = i4;
425 sommets(1) = i5;
426 sommets(2) = i7;
427 sommets(3) = i6;
428 indice(4) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 4, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
429
430 //centres des faces avant (0-1-5-4)
431 sommets(0) = i0;
432 sommets(1) = i1;
433 sommets(2) = i5;
434 sommets(3) = i4;
435 indice(3) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 4, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
436
437 //centres des faces arriere (2-3-7-6)
438 sommets(0) = i2;
439 sommets(1) = i3;
440 sommets(2) = i7;
441 sommets(3) = i6;
442 indice(6) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 4, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
443
444 //centres des faces gauche (0-2-6-4)
445 sommets(0) = i0;
446 sommets(1) = i2;
447 sommets(2) = i6;
448 sommets(3) = i4;
449 indice(2) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 4, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
450
451 //centres des faces droite (1-3-7-5)
452 sommets(0) = i1;
453 sommets(1) = i3;
454 sommets(2) = i7;
455 sommets(3) = i5;
456 indice(5) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 4, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
457
458 //centres des aretes (0-1)
459 sommets(0) = i0;
460 sommets(1) = i1;
461 indice(7) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 2, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
462
463 //centres des aretes (1-3)
464 sommets(0) = i1;
465 sommets(1) = i3;
466 indice(8) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 2, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
467
468 //centres des aretes (3-2)
469 sommets(0) = i3;
470 sommets(1) = i2;
471 indice(9) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 2, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
472
473 //centres des aretes (2-0)
474 sommets(0) = i2;
475 sommets(1) = i0;
476 indice(10) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 2, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
477
478 //centres des aretes (1-5)
479 sommets(0) = i1;
480 sommets(1) = i5;
481 indice(11) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 2, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
482
483 //centres des aretes (3-7)
484 sommets(0) = i3;
485 sommets(1) = i7;
486 indice(12) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 2, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
487
488 //centres des aretes (2-6)
489 sommets(0) = i2;
490 sommets(1) = i6;
491 indice(13) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 2, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
492
493 //centres des aretes (0-4)
494 sommets(0) = i0;
495 sommets(1) = i4;
496 indice(14) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 2, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
497
498 //centres des aretes (4-5)
499 sommets(0) = i4;
500 sommets(1) = i5;
501 indice(15) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 2, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
502
503 //centres des aretes (5-7)
504 sommets(0) = i5;
505 sommets(1) = i7;
506 indice(16) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 2, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
507
508 //centres des aretes (7-6)
509 sommets(0) = i7;
510 sommets(1) = i6;
511 indice(17) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 2, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
512
513 //centres des aretes (6-4)
514 sommets(0) = i6;
515 sommets(1) = i4;
516 indice(18) = creer_sommet(domaine, new_soms, elem_traite, new_soms_old_elems, 2, sommets, compteur, oldnbsom, nbnewsoms, fait_sommet, nface, fait_sommet_arete, narete);
517
518 // Liste des nouveaux sommets pour cet ancien cube
519 for (int t = 0; t < 19; t++)
520 {
521 new_soms_old_elems(i, t) = indice(t);
522 if (indice(t) < 0)
523 {
524 Cerr << "Error negative index" << finl;
525 exit();
526 }
527 for (int s = 0; s < t; s++)
528 if (indice(s) == indice(t))
529 {
530 Cerr << "Error repeated index" << finl;
531 exit();
532 }
533 }
534
535 for (int j = 0; j < compteur; j++)
536 for (int k = 0; k < 3; k++)
537 sommets_dom(dim_som_old + j, k) = new_soms(j, k);
538 dim_som_old += compteur;
539
540 // L'element en cours a ete "traite" en entier
541 elem_traite(i) = 1;
542
543 // pyramide bas (base : 0-1-3-2)
544
545 new_elems(i, 0) = les_elems(i, 0);
546 new_elems(i, 1) = indice(7);
547 new_elems(i, 2) = indice(1);
548 new_elems(i, 3) = indice(0);
549
550 new_elems(i + oldsz, 0) = indice(7);
551 new_elems(i + oldsz, 1) = les_elems(i, 1);
552 new_elems(i + oldsz, 2) = indice(1);
553 new_elems(i + oldsz, 3) = indice(0);
554 mettre_a_jour_sous_domaine(domaine, i, i + oldsz, 1);
555
556 new_elems(i + 2 * oldsz, 0) = les_elems(i, 1);
557 new_elems(i + 2 * oldsz, 1) = indice(8);
558 new_elems(i + 2 * oldsz, 2) = indice(1);
559 new_elems(i + 2 * oldsz, 3) = indice(0);
560 mettre_a_jour_sous_domaine(domaine, i, i + 2 * oldsz, 1);
561
562 new_elems(i + 3 * oldsz, 0) = indice(8);
563 new_elems(i + 3 * oldsz, 1) = les_elems(i, 3);
564 new_elems(i + 3 * oldsz, 2) = indice(1);
565 new_elems(i + 3 * oldsz, 3) = indice(0);
566 mettre_a_jour_sous_domaine(domaine, i, i + 3 * oldsz, 1);
567
568 new_elems(i + 4 * oldsz, 0) = les_elems(i, 3);
569 new_elems(i + 4 * oldsz, 1) = indice(9);
570 new_elems(i + 4 * oldsz, 2) = indice(1);
571 new_elems(i + 4 * oldsz, 3) = indice(0);
572 mettre_a_jour_sous_domaine(domaine, i, i + 4 * oldsz, 1);
573
574 new_elems(i + 5 * oldsz, 0) = indice(9);
575 new_elems(i + 5 * oldsz, 1) = les_elems(i, 2);
576 new_elems(i + 5 * oldsz, 2) = indice(1);
577 new_elems(i + 5 * oldsz, 3) = indice(0);
578 mettre_a_jour_sous_domaine(domaine, i, i + 5 * oldsz, 1);
579
580 new_elems(i + 6 * oldsz, 0) = les_elems(i, 2);
581 new_elems(i + 6 * oldsz, 1) = indice(10);
582 new_elems(i + 6 * oldsz, 2) = indice(1);
583 new_elems(i + 6 * oldsz, 3) = indice(0);
584 mettre_a_jour_sous_domaine(domaine, i, i + 6 * oldsz, 1);
585
586 new_elems(i + 7 * oldsz, 0) = indice(10);
587 new_elems(i + 7 * oldsz, 1) = les_elems(i, 0);
588 new_elems(i + 7 * oldsz, 2) = indice(1);
589 new_elems(i + 7 * oldsz, 3) = indice(0);
590 mettre_a_jour_sous_domaine(domaine, i, i + 7 * oldsz, 1);
591
592 // pyramide haut (base : 4-5-7-6)
593
594 new_elems(i + 8 * oldsz, 0) = les_elems(i, 4);
595 new_elems(i + 8 * oldsz, 1) = indice(15);
596 new_elems(i + 8 * oldsz, 2) = indice(4);
597 new_elems(i + 8 * oldsz, 3) = indice(0);
598 mettre_a_jour_sous_domaine(domaine, i, i + 8 * oldsz, 1);
599
600 new_elems(i + 9 * oldsz, 0) = indice(15);
601 new_elems(i + 9 * oldsz, 1) = les_elems(i, 5);
602 new_elems(i + 9 * oldsz, 2) = indice(4);
603 new_elems(i + 9 * oldsz, 3) = indice(0);
604 mettre_a_jour_sous_domaine(domaine, i, i + 9 * oldsz, 1);
605
606 new_elems(i + 10 * oldsz, 0) = les_elems(i, 5);
607 new_elems(i + 10 * oldsz, 1) = indice(16);
608 new_elems(i + 10 * oldsz, 2) = indice(4);
609 new_elems(i + 10 * oldsz, 3) = indice(0);
610 mettre_a_jour_sous_domaine(domaine, i, i + 10 * oldsz, 1);
611
612 new_elems(i + 11 * oldsz, 0) = indice(16);
613 new_elems(i + 11 * oldsz, 1) = les_elems(i, 7);
614 new_elems(i + 11 * oldsz, 2) = indice(4);
615 new_elems(i + 11 * oldsz, 3) = indice(0);
616 mettre_a_jour_sous_domaine(domaine, i, i + 11 * oldsz, 1);
617
618 new_elems(i + 12 * oldsz, 0) = les_elems(i, 7);
619 new_elems(i + 12 * oldsz, 1) = indice(17);
620 new_elems(i + 12 * oldsz, 2) = indice(4);
621 new_elems(i + 12 * oldsz, 3) = indice(0);
622 mettre_a_jour_sous_domaine(domaine, i, i + 12 * oldsz, 1);
623
624 new_elems(i + 13 * oldsz, 0) = indice(17);
625 new_elems(i + 13 * oldsz, 1) = les_elems(i, 6);
626 new_elems(i + 13 * oldsz, 2) = indice(4);
627 new_elems(i + 13 * oldsz, 3) = indice(0);
628 mettre_a_jour_sous_domaine(domaine, i, i + 13 * oldsz, 1);
629
630 new_elems(i + 14 * oldsz, 0) = les_elems(i, 6);
631 new_elems(i + 14 * oldsz, 1) = indice(18);
632 new_elems(i + 14 * oldsz, 2) = indice(4);
633 new_elems(i + 14 * oldsz, 3) = indice(0);
634 mettre_a_jour_sous_domaine(domaine, i, i + 14 * oldsz, 1);
635
636 new_elems(i + 15 * oldsz, 0) = indice(18);
637 new_elems(i + 15 * oldsz, 1) = les_elems(i, 4);
638 new_elems(i + 15 * oldsz, 2) = indice(4);
639 new_elems(i + 15 * oldsz, 3) = indice(0);
640 mettre_a_jour_sous_domaine(domaine, i, i + 15 * oldsz, 1);
641
642 // pyramide avant (base : 0-1-5-4)
643
644 new_elems(i + 16 * oldsz, 0) = les_elems(i, 0);
645 new_elems(i + 16 * oldsz, 1) = indice(7);
646 new_elems(i + 16 * oldsz, 2) = indice(3);
647 new_elems(i + 16 * oldsz, 3) = indice(0);
648 mettre_a_jour_sous_domaine(domaine, i, i + 16 * oldsz, 1);
649
650 new_elems(i + 17 * oldsz, 0) = indice(7);
651 new_elems(i + 17 * oldsz, 1) = les_elems(i, 1);
652 new_elems(i + 17 * oldsz, 2) = indice(3);
653 new_elems(i + 17 * oldsz, 3) = indice(0);
654 mettre_a_jour_sous_domaine(domaine, i, i + 17 * oldsz, 1);
655
656 new_elems(i + 18 * oldsz, 0) = les_elems(i, 1);
657 new_elems(i + 18 * oldsz, 1) = indice(11);
658 new_elems(i + 18 * oldsz, 2) = indice(3);
659 new_elems(i + 18 * oldsz, 3) = indice(0);
660 mettre_a_jour_sous_domaine(domaine, i, i + 18 * oldsz, 1);
661
662 new_elems(i + 19 * oldsz, 0) = indice(11);
663 new_elems(i + 19 * oldsz, 1) = les_elems(i, 5);
664 new_elems(i + 19 * oldsz, 2) = indice(3);
665 new_elems(i + 19 * oldsz, 3) = indice(0);
666 mettre_a_jour_sous_domaine(domaine, i, i + 19 * oldsz, 1);
667
668 new_elems(i + 20 * oldsz, 0) = les_elems(i, 5);
669 new_elems(i + 20 * oldsz, 1) = indice(15);
670 new_elems(i + 20 * oldsz, 2) = indice(3);
671 new_elems(i + 20 * oldsz, 3) = indice(0);
672 mettre_a_jour_sous_domaine(domaine, i, i + 20 * oldsz, 1);
673
674 new_elems(i + 21 * oldsz, 0) = indice(15);
675 new_elems(i + 21 * oldsz, 1) = les_elems(i, 4);
676 new_elems(i + 21 * oldsz, 2) = indice(3);
677 new_elems(i + 21 * oldsz, 3) = indice(0);
678 mettre_a_jour_sous_domaine(domaine, i, i + 21 * oldsz, 1);
679
680 new_elems(i + 22 * oldsz, 0) = les_elems(i, 4);
681 new_elems(i + 22 * oldsz, 1) = indice(14);
682 new_elems(i + 22 * oldsz, 2) = indice(3);
683 new_elems(i + 22 * oldsz, 3) = indice(0);
684 mettre_a_jour_sous_domaine(domaine, i, i + 22 * oldsz, 1);
685
686 new_elems(i + 23 * oldsz, 0) = indice(14);
687 new_elems(i + 23 * oldsz, 1) = les_elems(i, 0);
688 new_elems(i + 23 * oldsz, 2) = indice(3);
689 new_elems(i + 23 * oldsz, 3) = indice(0);
690 mettre_a_jour_sous_domaine(domaine, i, i + 23 * oldsz, 1);
691
692 // pyramide arriere (base : 2-3-7-6)
693
694 new_elems(i + 24 * oldsz, 0) = les_elems(i, 2);
695 new_elems(i + 24 * oldsz, 1) = indice(9);
696 new_elems(i + 24 * oldsz, 2) = indice(6);
697 new_elems(i + 24 * oldsz, 3) = indice(0);
698 mettre_a_jour_sous_domaine(domaine, i, i + 24 * oldsz, 1);
699
700 new_elems(i + 25 * oldsz, 0) = indice(9);
701 new_elems(i + 25 * oldsz, 1) = les_elems(i, 3);
702 new_elems(i + 25 * oldsz, 2) = indice(6);
703 new_elems(i + 25 * oldsz, 3) = indice(0);
704 mettre_a_jour_sous_domaine(domaine, i, i + 25 * oldsz, 1);
705
706 new_elems(i + 26 * oldsz, 0) = les_elems(i, 3);
707 new_elems(i + 26 * oldsz, 1) = indice(12);
708 new_elems(i + 26 * oldsz, 2) = indice(6);
709 new_elems(i + 26 * oldsz, 3) = indice(0);
710 mettre_a_jour_sous_domaine(domaine, i, i + 26 * oldsz, 1);
711
712 new_elems(i + 27 * oldsz, 0) = indice(12);
713 new_elems(i + 27 * oldsz, 1) = les_elems(i, 7);
714 new_elems(i + 27 * oldsz, 2) = indice(6);
715 new_elems(i + 27 * oldsz, 3) = indice(0);
716 mettre_a_jour_sous_domaine(domaine, i, i + 27 * oldsz, 1);
717
718 new_elems(i + 28 * oldsz, 0) = les_elems(i, 7);
719 new_elems(i + 28 * oldsz, 1) = indice(17);
720 new_elems(i + 28 * oldsz, 2) = indice(6);
721 new_elems(i + 28 * oldsz, 3) = indice(0);
722 mettre_a_jour_sous_domaine(domaine, i, i + 28 * oldsz, 1);
723
724 new_elems(i + 29 * oldsz, 0) = indice(17);
725 new_elems(i + 29 * oldsz, 1) = les_elems(i, 6);
726 new_elems(i + 29 * oldsz, 2) = indice(6);
727 new_elems(i + 29 * oldsz, 3) = indice(0);
728 mettre_a_jour_sous_domaine(domaine, i, i + 29 * oldsz, 1);
729
730 new_elems(i + 30 * oldsz, 0) = les_elems(i, 6);
731 new_elems(i + 30 * oldsz, 1) = indice(13);
732 new_elems(i + 30 * oldsz, 2) = indice(6);
733 new_elems(i + 30 * oldsz, 3) = indice(0);
734 mettre_a_jour_sous_domaine(domaine, i, i + 30 * oldsz, 1);
735
736 new_elems(i + 31 * oldsz, 0) = indice(13);
737 new_elems(i + 31 * oldsz, 1) = les_elems(i, 2);
738 new_elems(i + 31 * oldsz, 2) = indice(6);
739 new_elems(i + 31 * oldsz, 3) = indice(0);
740 mettre_a_jour_sous_domaine(domaine, i, i + 31 * oldsz, 1);
741
742 // pyramide gauche (base : 0-2-6-4)
743
744 new_elems(i + 32 * oldsz, 0) = les_elems(i, 0);
745 new_elems(i + 32 * oldsz, 1) = indice(10);
746 new_elems(i + 32 * oldsz, 2) = indice(2);
747 new_elems(i + 32 * oldsz, 3) = indice(0);
748 mettre_a_jour_sous_domaine(domaine, i, i + 32 * oldsz, 1);
749
750 new_elems(i + 33 * oldsz, 0) = indice(10);
751 new_elems(i + 33 * oldsz, 1) = les_elems(i, 2);
752 new_elems(i + 33 * oldsz, 2) = indice(2);
753 new_elems(i + 33 * oldsz, 3) = indice(0);
754 mettre_a_jour_sous_domaine(domaine, i, i + 33 * oldsz, 1);
755
756 new_elems(i + 34 * oldsz, 0) = les_elems(i, 2);
757 new_elems(i + 34 * oldsz, 1) = indice(13);
758 new_elems(i + 34 * oldsz, 2) = indice(2);
759 new_elems(i + 34 * oldsz, 3) = indice(0);
760 mettre_a_jour_sous_domaine(domaine, i, i + 34 * oldsz, 1);
761
762 new_elems(i + 35 * oldsz, 0) = indice(13);
763 new_elems(i + 35 * oldsz, 1) = les_elems(i, 6);
764 new_elems(i + 35 * oldsz, 2) = indice(2);
765 new_elems(i + 35 * oldsz, 3) = indice(0);
766 mettre_a_jour_sous_domaine(domaine, i, i + 35 * oldsz, 1);
767
768 new_elems(i + 36 * oldsz, 0) = les_elems(i, 6);
769 new_elems(i + 36 * oldsz, 1) = indice(18);
770 new_elems(i + 36 * oldsz, 2) = indice(2);
771 new_elems(i + 36 * oldsz, 3) = indice(0);
772 mettre_a_jour_sous_domaine(domaine, i, i + 36 * oldsz, 1);
773
774 new_elems(i + 37 * oldsz, 0) = indice(18);
775 new_elems(i + 37 * oldsz, 1) = les_elems(i, 4);
776 new_elems(i + 37 * oldsz, 2) = indice(2);
777 new_elems(i + 37 * oldsz, 3) = indice(0);
778 mettre_a_jour_sous_domaine(domaine, i, i + 37 * oldsz, 1);
779
780 new_elems(i + 38 * oldsz, 0) = les_elems(i, 4);
781 new_elems(i + 38 * oldsz, 1) = indice(14);
782 new_elems(i + 38 * oldsz, 2) = indice(2);
783 new_elems(i + 38 * oldsz, 3) = indice(0);
784 mettre_a_jour_sous_domaine(domaine, i, i + 38 * oldsz, 1);
785
786 new_elems(i + 39 * oldsz, 0) = indice(14);
787 new_elems(i + 39 * oldsz, 1) = les_elems(i, 0);
788 new_elems(i + 39 * oldsz, 2) = indice(2);
789 new_elems(i + 39 * oldsz, 3) = indice(0);
790 mettre_a_jour_sous_domaine(domaine, i, i + 39 * oldsz, 1);
791
792 // pyramide droite (base : 1-3-7-5)
793
794 new_elems(i + 40 * oldsz, 0) = les_elems(i, 1);
795 new_elems(i + 40 * oldsz, 1) = indice(8);
796 new_elems(i + 40 * oldsz, 2) = indice(5);
797 new_elems(i + 40 * oldsz, 3) = indice(0);
798 mettre_a_jour_sous_domaine(domaine, i, i + 40 * oldsz, 1);
799
800 new_elems(i + 41 * oldsz, 0) = indice(8);
801 new_elems(i + 41 * oldsz, 1) = les_elems(i, 3);
802 new_elems(i + 41 * oldsz, 2) = indice(5);
803 new_elems(i + 41 * oldsz, 3) = indice(0);
804 mettre_a_jour_sous_domaine(domaine, i, i + 41 * oldsz, 1);
805
806 new_elems(i + 42 * oldsz, 0) = les_elems(i, 3);
807 new_elems(i + 42 * oldsz, 1) = indice(12);
808 new_elems(i + 42 * oldsz, 2) = indice(5);
809 new_elems(i + 42 * oldsz, 3) = indice(0);
810 mettre_a_jour_sous_domaine(domaine, i, i + 42 * oldsz, 1);
811
812 new_elems(i + 43 * oldsz, 0) = indice(12);
813 new_elems(i + 43 * oldsz, 1) = les_elems(i, 7);
814 new_elems(i + 43 * oldsz, 2) = indice(5);
815 new_elems(i + 43 * oldsz, 3) = indice(0);
816 mettre_a_jour_sous_domaine(domaine, i, i + 43 * oldsz, 1);
817
818 new_elems(i + 44 * oldsz, 0) = les_elems(i, 7);
819 new_elems(i + 44 * oldsz, 1) = indice(16);
820 new_elems(i + 44 * oldsz, 2) = indice(5);
821 new_elems(i + 44 * oldsz, 3) = indice(0);
822 mettre_a_jour_sous_domaine(domaine, i, i + 44 * oldsz, 1);
823
824 new_elems(i + 45 * oldsz, 0) = indice(16);
825 new_elems(i + 45 * oldsz, 1) = les_elems(i, 5);
826 new_elems(i + 45 * oldsz, 2) = indice(5);
827 new_elems(i + 45 * oldsz, 3) = indice(0);
828 mettre_a_jour_sous_domaine(domaine, i, i + 45 * oldsz, 1);
829
830 new_elems(i + 46 * oldsz, 0) = les_elems(i, 5);
831 new_elems(i + 46 * oldsz, 1) = indice(11);
832 new_elems(i + 46 * oldsz, 2) = indice(5);
833 new_elems(i + 46 * oldsz, 3) = indice(0);
834 mettre_a_jour_sous_domaine(domaine, i, i + 46 * oldsz, 1);
835
836 new_elems(i + 47 * oldsz, 0) = indice(11);
837 new_elems(i + 47 * oldsz, 1) = les_elems(i, 1);
838 new_elems(i + 47 * oldsz, 2) = indice(5);
839 new_elems(i + 47 * oldsz, 3) = indice(0);
840 mettre_a_jour_sous_domaine(domaine, i, i + 47 * oldsz, 1);
841
842 }
843 sommets_dom.resize(dim_som_old, 3);
844 les_elems.ref(new_elems);
845
846 // Reconstruction de l'octree
847 Cerr << "We have split the cubes..." << finl;
848 domaine.invalide_octree();
849 domaine.typer("Tetraedre");
850 Cerr << " Reconstruction of the Octree" << finl;
851 domaine.construit_octree();
852 Cerr << " Octree rebuilt" << finl;
853
854 Cerr << "Splitting of the boundaries" << finl;
855 for (auto &itr : domaine.faces_bord())
856 {
857 Faces& les_faces = itr.faces();
858 les_faces.typer(Type_Face::triangle_3D);
859 decoupe(domaine, les_faces, new_soms_old_elems, fait_sommet, nface, fait_sommet_arete, narete);
860 }
861
862 Cerr << "Splitting of the connectors" << finl;
863 for (auto &itr : domaine.faces_raccord())
864 {
865 Faces& les_faces = itr->faces();
866 les_faces.typer(Type_Face::triangle_3D);
867 decoupe(domaine, les_faces, new_soms_old_elems, fait_sommet, nface, fait_sommet_arete, narete);
868 }
869
870 Cerr << "Splitting of the internal boundary faces" << finl;
871 for (auto &itr : domaine.bords_int())
872 {
873 Faces& les_faces = itr.faces();
874 les_faces.typer(Type_Face::triangle_3D);
875 decoupe(domaine, les_faces, new_soms_old_elems, fait_sommet, nface, fait_sommet_arete, narete);
876 }
877
878 Cerr << "Splitting of the group of faces" << finl;
879 for (auto &itr : domaine.groupes_faces())
880 {
881 Faces& les_faces = itr.faces();
882 les_faces.typer(Type_Face::triangle_3D);
883 decoupe(domaine, les_faces, new_soms_old_elems, fait_sommet, nface, fait_sommet_arete, narete);
884 }
885 Cerr << "END of Tetraedriser_homogene_fin..." << finl;
886 Cerr << " 1 NbElem=" << domaine.les_elems().dimension(0) << " NbNod=" << domaine.nb_som() << finl;
887 }
888 else
889 Cerr << "We do not yet know how to Tetraedriser_homogene_fin the " << domaine.type_elem()->que_suis_je() << "s" << finl;
890}
Entree
Class defining operators and methods for all reading operation in an input flow (file,...
Definition Entree.h:42
Faces_32_64::typer
void typer(const Motcle &)
Type les faces.
Definition Faces.cpp:390
Faces_32_64::voisins
IntTab_t & voisins()
Renvoie le tableau des voisins (des faces).
Definition Faces.h:89
Faces_32_64::les_sommets
const IntTab_t & les_sommets() const
Renvoie le tableau des sommets de toutes les faces.
Definition Faces.h:74
Interprete_geometrique_base_32_64< int >::mettre_a_jour_sous_domaine
void mettre_a_jour_sous_domaine(Domaine_t &domaine, int_t &elem, int_t num_premier_elem, int_t nb_elem) const
Interprete_geometrique_base_32_64< int >::domaine
Domaine_t & domaine(int i=0)
Definition Interprete_geometrique_base.h:48
Objet_U::readOn
virtual Entree & readOn(Entree &)
Lecture d'un Objet_U sur un flot d'entree Methode a surcharger.
Definition Objet_U.cpp:293
Objet_U::precision_geom
static double precision_geom
Definition Objet_U.h:86
Objet_U::printOn
virtual Sortie & printOn(Sortie &) const
Ecriture de l'objet sur un flot de sortie Methode a surcharger.
Definition Objet_U.cpp:282
Process::exit
static void exit(int exit_code=-1)
Routine de sortie de TRUST dans une region Kokkos.
Definition Process.cpp:455
Sortie
Classe de base des flux de sortie.
Definition Sortie.h:52
TRUSTTab::ref
virtual void ref(const TRUSTTab &)
Definition TRUSTTab.tpp:308
TRUSTTab::resize
void resize(_SIZE_ n, RESIZE_OPTIONS opt=RESIZE_OPTIONS::COPY_INIT)
Definition TRUSTTab.tpp:469
TRUSTTab::dimension
_SIZE_ dimension(int d) const
Definition TRUSTTab.tpp:133
Tetraedriser_homogene_fin
class Tetra_homogene_fin Realise un maillage en decoupant chaque pave en 48 tetraedres
Definition Tetraedriser_homogene_fin.h:37
Tetraedriser_homogene_fin::trianguler
void trianguler(Domaine &) const override
Definition Tetraedriser_homogene_fin.cpp:338
Tetraedriser_homogene_fin::creer_sommet
int creer_sommet(Domaine &, DoubleTab &, IntTab &, IntTab &, int, IntTab &, int &, int, int &, IntTab &, int &, IntTab &, int &) const
Definition Tetraedriser_homogene_fin.cpp:190
Tetraedriser_homogene_fin::decoupe
void decoupe(Domaine &, Faces &, IntTab &, IntTab &, int, IntTab &, int) const
Definition Tetraedriser_homogene_fin.cpp:37
Triangulation_base
Triangulation_base Classe destinee a factoriser l'action de triangulation des interpretes.
Definition Triangulation_base.h:29