1 // Copyright (C) 2009-2014 CEA/DEN, EDF R&D
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License, or (at your option) any later version.
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 // Lesser General Public License for more details.
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
20 // C++ : Table des noeuds
22 #include "HexElements.hxx"
24 #include "HexDocument.hxx"
25 #include "HexVector.hxx"
26 #include "HexVertex.hxx"
27 #include "HexHexa.hxx"
28 #include "HexEdge.hxx"
30 #include "HexGlobale.hxx"
35 // ====================================================== fillGrid
36 int Elements::fillGrid ()
39 for (int nx=0 ; nx<size_vx ; nx++)
40 for (int nz=0 ; nz<size_vz ; nz++)
41 setVertex (getVertexIJK (nx, 0, nz), nx, size_vy-1, nz);
43 for (int nz=0 ; nz<size_vz ; nz++)
44 for (int ny=0 ; ny<size_vy ; ny++)
45 for (int nx=0 ; nx<size_vx ; nx++)
47 Vertex* v0 = getVertexIJK (nx, ny, nz );
48 Vertex* vx = getVertexIJK (nx+1, ny, nz);
49 Vertex* vy = getVertexIJK (nx, ny+1, nz);
50 Vertex* vz = getVertexIJK (nx, ny, nz+1);
53 Edge* e2 = newEdge (v0, vy);
56 if (cyl_closed && ny==size_vy-1)
58 e1 = getEdgeI (nx, 0, nz);
59 e3 = getEdgeK (nx, 0, nz);
63 e3 = newEdge (v0, vz);
64 e1 = newEdge (v0, vx);
67 setEdge (e1, dir_x, nx, ny, nz);
68 setEdge (e2, dir_y, nx, ny, nz);
69 setEdge (e3, dir_z, nx, ny, nz);
72 for (int nz=0 ; nz<size_vz ; nz++)
73 for (int ny=0 ; ny<size_vy ; ny++)
74 for (int nx=0 ; nx<size_vx ; nx++)
76 Edge* eae = getEdgeI (nx, ny, nz);
77 Edge* ebe = getEdgeI (nx, ny, nz+1);
78 Edge* eaf = getEdgeI (nx, ny+1, nz);
80 Edge* eac = getEdgeJ (nx, ny, nz);
81 Edge* ead = getEdgeJ (nx+1, ny, nz);
82 Edge* ebc = getEdgeJ (nx, ny, nz+1);
84 Edge* ece = getEdgeK (nx, ny, nz);
85 Edge* ede = getEdgeK (nx+1, ny, nz);
86 Edge* ecf = getEdgeK (nx, ny+1, nz);
88 Quad* q1 = newQuad (eac, eaf, ead, eae);
89 Quad* q2 = newQuad (eac, ecf, ebc, ece);
91 Quad* q30 = getQuadIK (nx, 0 ,nz);
93 if (cyl_closed && ny==size_vy-1 && q30!= NULL)
100 q3 = newQuad (eae, ede, ebe, ece);
103 setQuad (q1, dir_z, nx,ny,nz);
104 setQuad (q2, dir_x, nx,ny,nz);
105 setQuad (q3, dir_y, nx,ny,nz);
108 for (int nz=0 ; nz<size_hz ; nz++)
109 for (int ny=0 ; ny<size_hy ; ny++)
110 for (int nx=0 ; nx<size_hx ; nx++)
112 Quad* qa = getQuadIJ (nx, ny, nz);
113 Quad* qb = getQuadIJ (nx, ny, nz+1);
115 Quad* qc = getQuadIK (nx, ny, nz);
116 Quad* qd = getQuadIK (nx, ny+1, nz);
118 Quad* qe = getQuadJK (nx, ny, nz);
119 Quad* qf = getQuadJK (nx+1, ny, nz);
121 setHexa (newHexa (qa, qb, qc, qd, qe, qf), nx, ny, nz);
127 case CYL_PEER : fillCenter ();
129 case CYL_CL4 : fillCenter4 ();
131 case CYL_CL6 : fillCenter6 ();
133 case CYL_ODD : fillCenterOdd ();
140 // ====================================================== fillCenter
141 // === Remplissage radial
142 #define IndElt(nc,nz) (nbsecteurs*(nz) + nc)
143 #define IndRedge(nc,nz) (nbrayons *(nz) + nc)
144 #define IndVquad(nc,nz) (nbrayons *(nz-1) + nc)
145 void Elements::fillCenter ()
148 int nbsecteurs = size_hy / 2;
149 int nbrayons = cyl_closed ? nbsecteurs : nbsecteurs + 1;
151 size_hplus = nbsecteurs * size_hz;
152 size_qhplus = nbsecteurs * size_vz;
153 size_qvplus = nbrayons * size_hz;
154 size_ehplus = nbrayons * size_vz;
155 size_evplus = size_hz;
157 ker_hexa .resize (size_hplus, NULL);
158 ker_hquad.resize (size_qhplus, NULL);
159 ker_vquad.resize (size_qvplus, NULL);
160 ker_hedge.resize (size_ehplus, NULL);
161 ker_vedge.resize (size_evplus, NULL);
163 Vertex* pcenter = NULL;
165 for (int nz=0 ; nz<size_vz ; nz++)
168 Vertex* center = getVertex (ker_vertex+nz);
169 // Edges horizontaux radiaux
170 for (int nc=0 ; nc<nbrayons ; nc++)
172 Vertex* vv = getVertexIJK (nx0, 2*nc, nz);
173 Edge* edge = newEdge (center, vv);
174 ker_hedge [IndRedge(nc,nz)] = edge;
177 for (int nc=0; nc<nbsecteurs ; nc++)
179 int nc1 = (nc + 1) MODULO nbrayons;
180 Edge* e1 = ker_hedge [IndRedge (nc, nz)];
181 Edge* e2 = getEdgeJ (nx0, 2*nc, nz);
182 Edge* e3 = getEdgeJ (nx0, 2*nc+1,nz);
183 Edge* e4 = ker_hedge [IndRedge (nc1, nz)];
185 ker_hquad [IndElt (nc,nz)] = newQuad (e1, e2, e3, e4);
188 printf ("hquad (%d,%d) = ", nc, nz);
189 ker_hquad [IndElt (nc,nz)]->dumpPlus();
195 // Edges verticaux + cloisons interieures
196 Edge* pilier = ker_vedge [nz-1] = newEdge (pcenter, center);
198 for (int nc=0 ; nc<nbrayons ; nc++)
200 Edge* e1 = ker_hedge [IndRedge (nc, nz)];
201 Edge* e2 = getEdgeK (nx0, 2*nc, nz-1);
202 Edge* e3 = ker_hedge [IndRedge (nc, nz-1)];
203 ker_vquad [IndVquad (nc,nz)] = newQuad (e1, e2, e3, pilier);
206 printf ("vquad (%d,%d) = ", nc, nz);
207 ker_vquad [IndVquad (nc,nz)]->dumpPlus();
211 for (int nc=0 ; nc < nbsecteurs ; nc++)
215 nc1 = nc1 MODULO nbsecteurs;
216 Quad* qa = ker_hquad [IndElt (nc, nz-1)];
217 Quad* qb = ker_hquad [IndElt (nc, nz)];
219 Quad* qc = ker_vquad [IndVquad (nc, nz)];
220 Quad* qd = getQuadJK (nx0, 2*nc+1, nz-1);
222 Quad* qe = getQuadJK (nx0, 2*nc, nz-1);
223 Quad* qf = ker_vquad [IndVquad (nc1, nz)];
227 printf (" --------------- Hexa : nc=%d, nz=%d\n", nc, nz);
235 Hexa* cell = newHexa (qa, qb, qc, qd, qe, qf);
236 tab_hexa.push_back (cell);
237 ker_hexa [IndElt (nc, nz-1)] = cell;
243 // ====================================================== fillCenter4
244 // === Remplissage radial
245 void Elements::fillCenter4 ()
250 for (int nz=0 ; nz<size_vz ; nz++)
254 Quad* plafond = newQuad (getEdgeJ (nx0, 0, nz), getEdgeJ (nx0, 1, nz),
255 getEdgeJ (nx0, 2, nz), getEdgeJ (nx0, 3, nz));
259 Quad* qc = getQuadJK (nx0, 0, nz-1);
260 Quad* qd = getQuadJK (nx0, 2, nz-1);
261 Quad* qe = getQuadJK (nx0, 1, nz-1);
262 Quad* qf = getQuadJK (nx0, 3, nz-1);
266 printf (" --------------- Hexa grille4 : nz=%d\n", nz);
274 Hexa* cell = newHexa (plafond, sol, qc, qd, qe, qf);
275 tab_hexa.push_back (cell);
280 // ====================================================== fillCenter6
281 void Elements::fillCenter6 ()
284 int nydemi = size_hy / 2;
286 Edge* s_barre = NULL;
287 Quad* sr_quad = NULL;
288 Quad* sl_quad = NULL;
290 for (int nz=0 ; nz<size_vz ; nz++)
292 // Edges horizontal radial
293 Edge* p_barre = newEdge (getVertexIJK (nx0, 0, nz),
294 getVertexIJK (nx0, nydemi, nz));
296 Edge* e0 = getEdgeJ (nx0, 0, nz);
297 Edge* e1 = getEdgeJ (nx0, 1, nz);
298 Edge* e2 = getEdgeJ (nx0, 2, nz);
299 Edge* e3 = getEdgeJ (nx0, 3, nz);
300 Edge* e4 = getEdgeJ (nx0, 4, nz);
301 Edge* e5 = getEdgeJ (nx0, 5, nz);
303 Quad* pl_quad = newQuad (e0, e1, e2, p_barre);
304 Quad* pr_quad = newQuad (e3, e4, e5, p_barre);
308 // Cloison interieure
309 Quad* cloison = newQuad (p_barre, getEdgeK (nx0, 0, nz-1),
310 s_barre, getEdgeK (nx0, nydemi, nz-1));
312 Quad* q0 = getQuadJK (nx0, 0, nz-1);
313 Quad* q1 = getQuadJK (nx0, 1, nz-1);
314 Quad* q2 = getQuadJK (nx0, 2, nz-1);
315 Quad* q3 = getQuadJK (nx0, 3, nz-1);
316 Quad* q4 = getQuadJK (nx0, 4, nz-1);
317 Quad* q5 = getQuadJK (nx0, 5, nz-1);
319 Hexa* left = newHexa (sl_quad, pl_quad, q0, q2, q1, cloison);
320 Hexa* right = newHexa (sr_quad, pr_quad, q3, q5, q4, cloison);
321 tab_hexa.push_back (left);
322 tab_hexa.push_back (right);
329 // ====================================================== fillCenterOdd
331 #define IndElt(nc,nz) (nbsecteurs*(nz) + nc)
332 void Elements::fillCenterOdd ()
335 int nbsecteurs = size_hy / 2;
337 vector <Edge*> ker_hedge (nbsecteurs*size_vz);
338 vector <Quad*> ker_hquad (nbsecteurs*size_vz);
339 vector <Quad*> ker_vquad (nbsecteurs*size_vz);
341 for (int nz=0 ; nz<size_vz ; nz++)
343 // Edges horizontaux radiaux
344 for (int nc=0 ; nc<nbsecteurs ; nc++)
346 int nc1 = size_hy - nc;
347 Edge* edge = newEdge (getVertexIJK (nx0, nc, nz),
348 getVertexIJK (nx0, nc1, nz));
349 ker_hedge [IndElt(nc,nz)] = edge;
352 for (int nc=0; nc<nbsecteurs ; nc++)
354 Edge* e1 = ker_hedge [IndElt (nc, nz)];
355 Edge* e2 = getEdgeJ (nx0, nc, nz);
356 Edge* e4 = getEdgeJ (nx0, size_hy-nc-1, nz);
358 Edge* e3 = nc<nbsecteurs-1 ? ker_hedge [IndElt (nc+1, nz)]
359 : getEdgeJ (nx0, nbsecteurs, nz);
361 ker_hquad [IndElt (nc,nz)] = newQuad (e1, e2, e3, e4);
364 printf ("hquad (%d,%d) = ", nc, nz);
365 ker_hquad [IndElt (nc,nz)]->dumpPlus();
371 // Edges verticaux + cloisons interieures
372 for (int nc=0 ; nc<nbsecteurs ; nc++)
374 int nc1 = size_hy - nc;
375 Edge* e1 = ker_hedge [IndElt (nc, nz)];
376 Edge* e3 = ker_hedge [IndElt (nc, nz-1)];
377 Edge* e2 = getEdgeK (nx0, nc, nz-1);
378 Edge* e4 = getEdgeK (nx0, nc1, nz-1);
379 ker_vquad [IndElt (nc,nz)] = newQuad (e1, e2, e3, e4);
382 printf ("vquad (%d,%d) = ", nc, nz);
383 ker_vquad [IndElt (nc,nz)]->dumpPlus();
387 for (int nc=0 ; nc < nbsecteurs ; nc++)
389 int nc1 = size_hy - nc-1;
390 Quad* qa = ker_hquad [IndElt (nc, nz-1)];
391 Quad* qb = ker_hquad [IndElt (nc, nz)];
393 Quad* qc = getQuadJK (nx0, nc, nz-1);
394 Quad* qd = getQuadJK (nx0, nc1, nz-1);
396 Quad* qe = ker_vquad [IndElt (nc, nz)];
397 Quad* qf = nc<nbsecteurs-1 ? ker_vquad [IndElt (nc+1, nz)]
398 : getQuadJK (nx0, nbsecteurs, nz-1);
401 printf (" --------------- Hexa : nc=%d, nz=%d\n", nc, nz);
409 Hexa* cell = newHexa (qa, qb, qc, qd, qe, qf);
410 tab_hexa.push_back (cell);
415 // --------------------------------------------------------------------------
416 // ----------------------------------------- Evols Hexa 3
417 // --------------------------------------------------------------------------
419 // ====================================================== makeBasicCylinder
420 // ==== Version avec vecteurs
421 int Elements::makeBasicCylinder (RealVector& tdr, RealVector& tda,
422 RealVector& tdh, bool fill)
424 int na = tda.size()-1;
425 cyl_dispo = CYL_NOFILL;
434 else if (na MODULO 2 == 0)
435 cyl_dispo = CYL_CLOSED;
437 else if ((na MODULO 2)==0)
438 cyl_dispo = CYL_PEER;
443 cyl_fill = cyl_dispo != CYL_NOFILL;
446 int nb_secteurs = cyl_closed ? size_vy-1 : size_vy;
448 for (int ny=0 ; ny<nb_secteurs ; ny++)
451 double theta = M_PI*alpha/180;
452 double cos_theta = cos (theta);
453 double sin_theta = sin (theta);
455 for (int nx=0 ; nx<size_vx ; nx++)
457 double rayon = tdr [nx];
458 double px = rayon*cos_theta;
459 double py = rayon*sin_theta;
461 for (int nz=0 ; nz<size_vz ; nz++)
463 double pz = tdh [nz];
464 Vertex* node = el_root->addVertex (px, py, pz);
465 setVertex (node, nx, ny, nz);
472 for (int nx=0 ; nx<size_vx ; nx++)
473 for (int nz=0 ; nz<size_vz ; nz++)
475 Vertex* node = getVertexIJK (nx, 0, nz);
476 setVertex (node, nx, size_vy-1, nz);
479 // Les vertex centraux
482 ker_vertex = nbr_vertex;
483 for (int nz=0 ; nz<size_vz ; nz++)
485 double pz = tdh [nz];
486 Vertex* node = el_root->addVertex (0, 0, pz);
487 tab_vertex.push_back (node);