2 // C++ : Table d'hexaedres
4 // Copyright (C) 2009-2014 CEA/DEN, EDF R&D
6 // This library is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 2.1 of the License, or (at your option) any later version.
11 // This library is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14 // Lesser General Public License for more details.
16 // You should have received a copy of the GNU Lesser General Public
17 // License along with this library; if not, write to the Free Software
18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
23 #include "HexElements.hxx"
25 #include "HexDocument.hxx"
26 #include "HexVector.hxx"
27 #include "HexVertex.hxx"
28 #include "HexHexa.hxx"
29 #include "HexEdge.hxx"
31 #include "HexGlobale.hxx"
32 #include "HexNewShape.hxx"
38 // ====================================================== getHexaIJK
39 Hexa* Elements::getHexaIJK (int nx, int ny, int nz)
43 if (nx<0 || nx>=size_hx || ny<0 || ny>=size_hy || nz<0 || nz>=size_hz)
48 int nro = nx + size_hx*ny + size_hx*size_hy*nz;
50 DumpStart ("getHexaIJK", nx << ny << nz);
51 DumpReturn (tab_hexa [nro]);
52 return tab_hexa [nro];
54 // ====================================================== getQuadIJ
55 Quad* Elements::getQuadIJ (int nx, int ny, int nz)
59 if (nx<0 || nx>=size_qx || ny<0 || ny>=size_qy || nz<0 || nz>=size_qz)
64 int nro = nx + size_qx*ny + size_qx*size_qy*nz
65 + size_qx*size_qy*size_qz*dir_z;
67 DumpStart ("getQuadIJ", nx << ny << nz);
68 DumpReturn (tab_quad [nro]);
69 return tab_quad [nro];
71 // ====================================================== getQuadJK
72 Quad* Elements::getQuadJK (int nx, int ny, int nz)
76 if (nx<0 || nx>=size_qx || ny<0 || ny>=size_qy || nz<0 || nz>=size_qz)
81 int nro = nx + size_qx*ny + size_qx*size_qy*nz; // + dir_x*...
83 DumpStart ("getQuadJK", nx << ny << nz);
84 DumpReturn (tab_quad [nro]);
85 return tab_quad [nro];
87 // ====================================================== getQuadIK
88 Quad* Elements::getQuadIK (int nx, int ny, int nz)
92 if (nx<0 || nx>=size_qx || ny<0 || ny>=size_qy || nz<0 || nz>=size_qz)
97 int nro = nx + size_qx*ny + size_qx*size_qy*nz + size_qx*size_qy*size_qz;
99 DumpStart ("getQuadIK", nx << ny << nz);
100 DumpReturn (tab_quad [nro]);
101 return tab_quad [nro];
103 // ====================================================== getEdgeI
104 Edge* Elements::getEdgeI (int nx, int ny, int nz)
108 if (nx<0 || nx>=size_ex || ny<0 || ny>=size_ey || nz<0 || nz>=size_ez)
110 else if (grid_nocart)
113 int nro = nx + size_ex*ny + size_ex*size_ey*nz;
115 DumpStart ("getEdgeI", nx << ny << nz);
116 DumpReturn (tab_edge [nro]);
117 return tab_edge [nro];
119 // ====================================================== getEdgeJ
120 Edge* Elements::getEdgeJ (int nx, int ny, int nz)
124 if (nx<0 || nx>=size_ex || ny<0 || ny>=size_ey || nz<0 || nz>=size_ez)
126 else if (grid_nocart)
129 int nro = nx + size_ex*ny + size_ex*size_ey*nz + size_ex*size_ey*size_ez;
131 DumpStart ("getEdgeJ", nx << ny << nz);
132 DumpReturn (tab_edge [nro]);
133 return tab_edge [nro];
135 // ====================================================== getEdgeK
136 Edge* Elements::getEdgeK (int nx, int ny, int nz)
140 if (nx<0 || nx>=size_ex || ny<0 || ny>=size_ey || nz<0 || nz>=size_ez)
142 else if (grid_nocart)
145 int nro = nx + size_ex*ny + size_ex*size_ey*nz
146 + size_ex*size_ey*size_ez*dir_z;
148 DumpStart ("getEdgeK", nx << ny << nz);
149 DumpReturn (tab_edge [nro]);
150 return tab_edge [nro];
152 // ====================================================== getVertexIJK
153 Vertex* Elements::getVertexIJK (int nx, int ny, int nz)
157 if (nx<0 || nx>=size_vx || ny<0 || ny>=size_vy || nz<0 || nz>=size_vz)
159 else if (grid_nocart)
162 int nro = nx + size_vx*ny + size_vx*size_vy*nz;
164 DumpStart ("getVertexIJK", nx << ny << nz);
165 DumpReturn (tab_vertex [nro]);
166 return tab_vertex [nro];
168 // ====================================================== setVertex
169 void Elements::setVertex (Vertex* elt, int nx, int ny, int nz)
171 if ( nx < 0 || nx >= size_vx || ny < 0 || ny >= size_vy
172 || nz < 0 || nz >= size_vz) return;
174 int nro = nx + size_vx*ny + size_vx*size_vy*nz;
175 tab_vertex [nro] = elt;
177 // ====================================================== setVertex (2)
178 void Elements::setVertex (int nx, int ny, int nz, double px, double py,
181 if ( nx < 0 || nx >= size_vx || ny < 0 || ny >= size_vy
182 || nz < 0 || nz >= size_vz) return;
184 Vertex* node = el_root->addVertex (px, py, pz);
185 setVertex (node, nx, ny, nz);
187 // ====================================================== setEdge
188 void Elements::setEdge (Edge* elt, EnumCoord dir, int nx, int ny, int nz)
190 if (nx<0 || nx>=size_ex || ny<0 || ny>=size_ey || nz<0 || nz>=size_ez
191 || dir < dir_x || dir > dir_z )
194 int nro = nx + size_ex*ny + size_ex*size_ey*nz + size_ex*size_ey*size_ez*dir;
195 tab_edge [nro] = elt;
197 // ====================================================== setQuad
198 void Elements::setQuad (Quad* elt, EnumCoord dir, int nx, int ny, int nz)
200 if (nx<0 || nx>=size_ex || ny<0 || ny>=size_ey || nz<0 || nz>=size_ez
201 || dir < dir_x || dir > dir_z )
204 int nro = nx + size_ex*ny + size_ex*size_ey*nz + size_ex*size_ey*size_ez*dir;
205 tab_quad [nro] = elt;
207 // ====================================================== setHexa
208 void Elements::setHexa (Hexa* elt, int nx, int ny, int nz)
210 if ( nx < 0 || nx >= size_hx || ny < 0 || ny >= size_hy
211 || nz < 0 || nz >= size_hz) return;
213 int nro = nx + size_hx*ny + size_hx*size_hy*nz;
214 tab_hexa [nro] = elt;
216 // ====================================================== remove
217 void Elements::remove ()
219 int nbre=tab_hexa.size ();
221 for (int nh=0 ; nh<nbre ; nh++)
222 if (tab_hexa[nh] != NULL)
223 tab_hexa[nh]->remove();
225 // ====================================================== prismHexas
226 int Elements::prismHexas (int nro, Quad* qbase, int hauteur)
228 int ind_node [QUAD4];
230 Real3 koord, transfo;
232 // ----------------------------- Vertex + aretes verticales
233 for (int ns=0 ; ns<QUAD4 ; ns++)
235 Vertex* vbase = qbase->getVertex (ns);
236 int indx = vbase->getMark ();
239 bool asso = vbase->isAssociated();
240 vbase->getAssoCoord (koord);
243 vbase->setMark (indx);
249 Real3 centre, point, om;
250 revo_center->getPoint (centre);
251 vbase ->getPoint (point);
253 calc_vecteur (centre, point, om);
254 double oh = prod_scalaire (om, revo_axe);
257 for (int dd=dir_x; dd<=dir_z ; dd++)
259 ph [dd] = centre [dd] + oh*revo_axe[dd];
260 hm [dd] = point [dd] - ph[dd];
261 rayon += hm[dd] * hm[dd];
263 subid = grid_geom->addCircle (ph, sqrt(rayon), revo_axe, hm);
266 for (int nh=0 ; nh<hauteur ; nh++)
268 nd1 = el_root->addVertex (nd0->getX(), nd0->getY(), nd0->getZ());
270 gen_matrix.perform (nd1);
271 tab_vertex.push_back (nd1);
272 Edge* pilier = newEdge (nd0, nd1);
273 tab_pilier.push_back (pilier);
277 beta = alpha + gen_values[nh];
278 grid_geom->addAssociation (pilier, subid, alpha/360, beta/360);
282 cum_matrix.perform (koord, transfo);
283 nd1->setAssociation (transfo);
288 ind_node [ns] = indx;
290 // ----------------------------- Aretes horizontales
291 // ----------------------------- + face verticales
292 int ind_poutre [QUAD4];
293 for (int ns=0 ; ns<QUAD4 ; ns++)
295 Edge* ebase = qbase->getEdge (ns);
296 int indx = ebase->getMark ();
300 ebase->setMark (indx);
301 int nd1 = ind_node [ns];
302 int nd2 = ind_node [(ns+1) MODULO QUAD4];
304 Edge *ed1, *ed2, *ed3;
305 for (int nh=0 ; nh<hauteur ; nh++)
308 ed0 = newEdge (tab_vertex [nd1*hauteur + nh],
309 tab_vertex [nd2*hauteur + nh]);
310 ed1 = tab_pilier [nd1*hauteur + nh];
311 ed3 = tab_pilier [nd2*hauteur + nh];
313 Quad* mur = newQuad (ed0, ed1, ed2, ed3);
314 tab_edge.push_back (ed0);
315 tab_quad.push_back (mur);
316 if (ebase->isAssociated ())
317 prismAssociation (ebase, ed0, nh);
320 ind_poutre [ns] = indx;
322 // ----------------------------- Faces horizontales
323 // ----------------------------- + Hexas
325 Quad *qb, *qc, *qd, *qe, *qf;
326 int nv0 = hauteur*ind_poutre [0];
327 int nv1 = hauteur*ind_poutre [1];
328 int nv2 = hauteur*ind_poutre [2];
329 int nv3 = hauteur*ind_poutre [3];
330 for (int nh=0 ; nh<hauteur ; nh++)
332 qb = newQuad (tab_edge [nh+nv0], tab_edge [nh+nv1],
333 tab_edge [nh+nv2], tab_edge [nh+nv3]);
334 qc = tab_quad [nh + nv0];
335 qd = tab_quad [nh + nv2];
336 qe = tab_quad [nh + nv1];
337 qf = tab_quad [nh + nv3];
339 // *** tab_hexa [nh*hauteur + nro] = newHexa (qa, qb, qc, qd, qe, qf); Abu
340 tab_hexa [nro*hauteur + nh] = newHexa (qa, qb, qc, qd, qe, qf);
341 ker_hquad.push_back (qb);
346 // ====================================================== updateMatrix
347 void Elements::updateMatrix (int nh)
351 gen_matrix.defRotation (revo_center, revo_axe, gen_values[nh]);
352 cum_matrix.defRotation (revo_center, revo_axe, cum_values[nh]);
356 double hauteur = (nh+1)*prism_len;
357 double dh = prism_len;
362 hauteur = gen_values[nh];
363 dh = nh>0 ? hauteur-gen_values[nh-1] : hauteur;
367 hauteur = cum_values[nh];
372 for (int nc=dir_x ; nc<=dir_z ; nc++)
374 decal [nc] = prism_dir [nc]*dh;
375 trans [nc] = prism_dir [nc]*hauteur;
378 gen_matrix.defTranslation (decal);
379 cum_matrix.defTranslation (trans);
382 // ====================================================== endPrism
383 void Elements::endPrism ()
387 int nbelts = ker_hquad.size();
388 for (int nro=0 ; nro<nbelts ; nro++)
389 tab_quad.push_back (ker_hquad[nro]);
391 nbelts = tab_pilier.size();
392 for (int nro=0 ; nro<nbelts ; nro++)
393 tab_edge.push_back (tab_pilier[nro]);
396 nbr_hexas = tab_hexa.size ();
397 nbr_edges = tab_edge.size ();
398 nbr_quads = tab_quad.size ();
399 nbr_vertex = tab_vertex.size ();
401 // ====================================================== getShape
402 NewShape* Elements::getShape()
406 string name = "extrud_" + el_name;
407 grid_geom = el_root->addShape (name.c_str(), SH_EXTRUD);
408 grid_geom -> openShape();
413 // ====================================================== closeShape
414 void Elements::closeShape()
419 grid_geom -> closeShape();