// C++ : Table d'hexaedres
-// Copyright (C) 2009-2012 CEA/DEN, EDF R&D
+// Copyright (C) 2009-2019 CEA/DEN, EDF R&D
//
// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
-// version 2.1 of the License.
+// version 2.1 of the License, or (at your option) any later version.
//
// This library is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
#include "HexEdge.hxx"
#include "HexGlobale.hxx"
-#include "HexCylinder.hxx"
-#include "HexShape.hxx"
+#include "HexNewShape.hxx"
#include <map>
BEGIN_NAMESPACE_HEXA
-void geom_dump_asso (Edge* edge);
-void geom_create_circle (double* milieu, double rayon, double* normale,
- double* base, string& brep);
-
// ====================================================== getHexaIJK
Hexa* Elements::getHexaIJK (int nx, int ny, int nz)
{
else if (grid_nocart)
return NULL;
- int nro = nx + size_hx*ny + size_hx*size_hy*nz;
+ int nro = nx + size_hx*ny + size_hx*size_hy*nz;
- return tab_hexa [nro];
+ DumpStart ("getHexaIJK", nx << ny << nz);
+ DumpReturn (tab_hexa [nro]);
+ return tab_hexa [nro];
}
// ====================================================== getQuadIJ
Quad* Elements::getQuadIJ (int nx, int ny, int nz)
else if (grid_nocart)
return NULL;
- int nro = nx + size_qx*ny + size_qx*size_qy*nz
+ int nro = nx + size_qx*ny + size_qx*size_qy*nz
+ size_qx*size_qy*size_qz*dir_z;
- return tab_quad [nro];
+
+ DumpStart ("getQuadIJ", nx << ny << nz);
+ DumpReturn (tab_quad [nro]);
+ return tab_quad [nro];
}
// ====================================================== getQuadJK
Quad* Elements::getQuadJK (int nx, int ny, int nz)
int nro = nx + size_qx*ny + size_qx*size_qy*nz; // + dir_x*...
- return tab_quad [nro];
+ DumpStart ("getQuadJK", nx << ny << nz);
+ DumpReturn (tab_quad [nro]);
+ return tab_quad [nro];
}
// ====================================================== getQuadIK
Quad* Elements::getQuadIK (int nx, int ny, int nz)
int nro = nx + size_qx*ny + size_qx*size_qy*nz + size_qx*size_qy*size_qz;
- return tab_quad [nro];
+ DumpStart ("getQuadIK", nx << ny << nz);
+ DumpReturn (tab_quad [nro]);
+ return tab_quad [nro];
}
// ====================================================== getEdgeI
Edge* Elements::getEdgeI (int nx, int ny, int nz)
int nro = nx + size_ex*ny + size_ex*size_ey*nz;
- return tab_edge [nro];
+ DumpStart ("getEdgeI", nx << ny << nz);
+ DumpReturn (tab_edge [nro]);
+ return tab_edge [nro];
}
// ====================================================== getEdgeJ
Edge* Elements::getEdgeJ (int nx, int ny, int nz)
int nro = nx + size_ex*ny + size_ex*size_ey*nz + size_ex*size_ey*size_ez;
- return tab_edge [nro];
+ DumpStart ("getEdgeJ", nx << ny << nz);
+ DumpReturn (tab_edge [nro]);
+ return tab_edge [nro];
}
// ====================================================== getEdgeK
Edge* Elements::getEdgeK (int nx, int ny, int nz)
else if (grid_nocart)
return NULL;
- int nro = nx + size_ex*ny + size_ex*size_ey*nz
+ int nro = nx + size_ex*ny + size_ex*size_ey*nz
+ size_ex*size_ey*size_ez*dir_z;
- return tab_edge [nro];
+
+ DumpStart ("getEdgeK", nx << ny << nz);
+ DumpReturn (tab_edge [nro]);
+ return tab_edge [nro];
}
// ====================================================== getVertexIJK
Vertex* Elements::getVertexIJK (int nx, int ny, int nz)
else if (grid_nocart)
return NULL;
- int nro = nx + size_vx*ny + size_vx*size_vy*nz;
+ int nro = nx + size_vx*ny + size_vx*size_vy*nz;
- return tab_vertex [nro];
+ DumpStart ("getVertexIJK", nx << ny << nz);
+ DumpReturn (tab_vertex [nro]);
+ return tab_vertex [nro];
}
// ====================================================== setVertex
void Elements::setVertex (Vertex* elt, int nx, int ny, int nz)
{
- if ( nx < 0 || nx >= size_vx || ny < 0 || ny >= size_vy
- || nz < 0 || nz >= size_vz) return;
+ if ( nx < 0 || nx >= size_vx || ny < 0 || ny >= size_vy
+ || nz < 0 || nz >= size_vz) return;
int nro = nx + size_vx*ny + size_vx*size_vy*nz;
tab_vertex [nro] = elt;
}
// ====================================================== setVertex (2)
-void Elements::setVertex (int nx, int ny, int nz, double px, double py,
+void Elements::setVertex (int nx, int ny, int nz, double px, double py,
double pz)
{
- if ( nx < 0 || nx >= size_vx || ny < 0 || ny >= size_vy
- || nz < 0 || nz >= size_vz) return;
+ if ( nx < 0 || nx >= size_vx || ny < 0 || ny >= size_vy
+ || nz < 0 || nz >= size_vz) return;
Vertex* node = el_root->addVertex (px, py, pz);
setVertex (node, nx, ny, nz);
return;
int nro = nx + size_ex*ny + size_ex*size_ey*nz + size_ex*size_ey*size_ez*dir;
- tab_edge [nro] = elt;
+ tab_edge [nro] = elt;
}
// ====================================================== setQuad
void Elements::setQuad (Quad* elt, EnumCoord dir, int nx, int ny, int nz)
return;
int nro = nx + size_ex*ny + size_ex*size_ey*nz + size_ex*size_ey*size_ez*dir;
- tab_quad [nro] = elt;
+ tab_quad [nro] = elt;
}
// ====================================================== setHexa
void Elements::setHexa (Hexa* elt, int nx, int ny, int nz)
{
- if ( nx < 0 || nx >= size_hx || ny < 0 || ny >= size_hy
- || nz < 0 || nz >= size_hz) return;
+ if ( nx < 0 || nx >= size_hx || ny < 0 || ny >= size_hy
+ || nz < 0 || nz >= size_hz) return;
int nro = nx + size_hx*ny + size_hx*size_hy*nz;
tab_hexa [nro] = elt;
if (tab_hexa[nh] != NULL)
tab_hexa[nh]->remove();
}
-// ====================================================== makeCylinder
-int Elements::makeCylinder (Cylinder* cyl, Vector* vx, int nr, int na, int nl)
-{
- if (BadElement (cyl) || BadElement (vx) || nr<=0 || na <=3 || nl <=0
- || vx->getNorm () <= Epsil)
- {
- setError ();
- return HERR;
- }
-
- Vertex* orig = cyl->getBase ();
- Vector* dir = cyl->getDirection ();
- double ray = cyl->getRadius ();
- double haut = cyl->getHeight ();
-
- resize (GR_CYLINDRIC, nr, na, nl);
- cyl_closed = true;
- makeCylindricalNodes (orig, vx, dir, ray/(nr+1), 360, haut/nl,
- nr, na, nl, true);
- fillGrid ();
- assoCylinder (orig, dir, 360);
- return HOK;
-}
-// ====================================================== makePipe
-int Elements::makePipe (Cylinder* cyl, Vector* vx, int nr, int na, int nl)
-{
- if (BadElement (cyl) || BadElement (vx) || nr<=0 || na <=3 || nl <=0
- || vx->getNorm () <= Epsil)
- {
- setError ();
- return HERR;
- }
-
- Vertex* orig = cyl->getBase ();
- Vector* dir = cyl->getDirection ();
- double ray = cyl->getRadius ();
- double haut = cyl->getHeight ();
-
- resize (GR_CYLINDRIC, nr, na, nl);
- cyl_closed = true;
- makeCylindricalNodes (orig, vx, dir, ray, 360, haut, nr, na, nl, false);
- fillGrid ();
- assoCylinder (orig, dir, 360);
- return HOK;
-}
-//
-// ---------------------------------------- prism Quads
-//
-// ====================================================== prismQuads
-int Elements::prismQuads (Quads& tstart, Vector* dir, int nbiter)
-{
- if (BadElement (dir) || dir->getNorm () <= Epsil || nbiter <= 0)
- {
- setError ();
- return HERR;
- }
-
- el_root->markAll (NO_USED);
- int nbcells = tstart.size ();
- nbr_vertex = 0;
- nbr_edges = 0;
-
- nbr_hexas = nbcells*nbiter;
-
- tab_hexa.resize (nbr_hexas);
- tab_quad.clear (); // verticaux
- ker_hquad.clear (); // Horizontaux
- tab_edge.clear ();
- tab_pilier.clear ();
- tab_vertex.clear ();
-
- revo_lution = false;
- prism_vec = false;
- gen_matrix.defTranslation (dir);
-
- for (int nro=0 ; nro<nbcells ; nro++)
- {
- prismHexas (nro, tstart[nro], nbiter);
- }
-
- endPrism ();
- return HOK;
-}
-// ====================================================== prismQuadsVec
-int Elements::prismQuadsVec (Quads& tstart, Vector* dir, RealVector& tlen,
- int mode)
-{
- int nbiter = tlen.size();
- if (BadElement (dir) || dir->getNorm () <= Epsil || nbiter <= 0)
- {
- setError ();
- return HERR;
- }
-
- el_root->markAll (NO_USED);
- int nbcells = tstart.size ();
- nbr_vertex = 0;
- nbr_edges = 0;
-
- nbr_hexas = nbcells*nbiter;
-
- tab_hexa.resize (nbr_hexas);
- tab_quad.clear (); // verticaux
- ker_hquad.clear (); // Horizontaux
- tab_edge.clear ();
- tab_pilier.clear ();
- tab_vertex.clear ();
-
- revo_lution = false;
- prism_vec = true;
- dir->getCoord (prism_dir);
- normer_vecteur (prism_dir);
- gen_values = tlen;
-
- for (int nro=0 ; nro<nbcells ; nro++)
- {
- prismHexas (nro, tstart[nro], nbiter);
- }
-
- endPrism ();
- return HOK;
-}
-// ======================================================== revolutionQuads
-int Elements::revolutionQuads (Quads& start, Vertex* center, Vector* axis,
- RealVector &angles)
-{
- int nbiter = angles.size();
- int nbcells = start.size ();
- if (BadElement (center) || BadElement(axis) || nbiter==0 || nbcells==0
- || axis->getNorm () <= Epsil)
- {
- setError ();
- return HERR;
- }
-
- el_root->markAll (NO_USED);
- nbr_vertex = 0;
- nbr_edges = 0;
-
- nbr_hexas = nbcells*nbiter;
-
- tab_hexa.resize (nbr_hexas);
- tab_quad.clear (); // verticaux
- ker_hquad.clear (); // Horizontaux
- tab_edge.clear ();
- tab_pilier.clear ();
- tab_vertex.clear ();
-
- revo_lution = true;
- prism_vec = false;
- revo_axis = axis;
- revo_center = center;
- gen_values = angles;
-
- for (int nro=0 ; nro<nbcells ; nro++)
- {
- prismHexas (nro, start[nro], nbiter);
- }
-
- endPrism ();
- return HOK;
-}
// ====================================================== prismHexas
int Elements::prismHexas (int nro, Quad* qbase, int hauteur)
{
- int ind_node [QUAD4];
- string c_rep;
+ int ind_node [QUAD4];
+ int subid = 0;
+ Real3 koord, transfo;
// ----------------------------- Vertex + aretes verticales
for (int ns=0 ; ns<QUAD4 ; ns++)
{
- Vertex* vbase = qbase ->getVertex (ns);
+ Vertex* vbase = qbase->getVertex (ns);
int indx = vbase->getMark ();
if (indx<0)
{
+ bool asso = vbase->isAssociated();
+ vbase->getAssoCoord (koord);
+
indx = nbr_vertex++;
vbase->setMark (indx);
Vertex* nd0 = vbase;
double beta = 0;
if (revo_lution)
{
- Real3 centre, vk, point, om;
+ Real3 centre, point, om;
revo_center->getPoint (centre);
vbase ->getPoint (point);
- revo_axis ->getCoord (vk);
- normer_vecteur (vk);
- calc_vecteur (centre, point, om);
- double oh = prod_scalaire (om, vk);
- double rayon = 0;
+ calc_vecteur (centre, point, om);
+ double oh = prod_scalaire (om, revo_axe);
+ double rayon = 0;
Real3 ph, hm;
for (int dd=dir_x; dd<=dir_z ; dd++)
{
- ph [dd] = centre [dd] + oh*vk[dd];
+ ph [dd] = centre [dd] + oh*revo_axe[dd];
hm [dd] = point [dd] - ph[dd];
rayon += hm[dd] * hm[dd];
}
- rayon = sqrt (rayon);
-/********************************
- PutCoord (centre);
- PutCoord (point);
- PutData (oh);
- PutCoord (ph);
- PutData (rayon);
- PutCoord (vk);
- PutCoord (hm);
-********************************/
- geom_create_circle (ph, rayon, vk, hm, c_rep);
+ subid = grid_geom->addCircle (ph, sqrt(rayon), revo_axe, hm);
}
for (int nh=0 ; nh<hauteur ; nh++)
if (revo_lution)
{
double alpha = beta;
- beta = alpha + gen_values[nh];
- Shape* shape = new Shape (c_rep);
- shape->setBounds (alpha/360, beta/360);
- pilier->addAssociation (shape);
- // geom_dump_asso (pilier);
+ beta = alpha + gen_values[nh];
+ grid_geom->addAssociation (pilier, subid, alpha/360, beta/360);
+ }
+ if (asso)
+ {
+ cum_matrix.perform (koord, transfo);
+ nd1->setAssociation (transfo);
}
nd0 = nd1;
}
for (int nh=0 ; nh<hauteur ; nh++)
{
ed2 = ed0;
- ed0 = newEdge (tab_vertex [nd1*hauteur + nh],
+ ed0 = newEdge (tab_vertex [nd1*hauteur + nh],
tab_vertex [nd2*hauteur + nh]);
ed1 = tab_pilier [nd1*hauteur + nh];
ed3 = tab_pilier [nd2*hauteur + nh];
Quad* mur = newQuad (ed0, ed1, ed2, ed3);
tab_edge.push_back (ed0);
tab_quad.push_back (mur);
- prismAssociation (ed2, ed0, nh, ed1);
+ if (ebase->isAssociated ())
+ prismAssociation (ebase, ed0, nh);
}
}
ind_poutre [ns] = indx;
int nv3 = hauteur*ind_poutre [3];
for (int nh=0 ; nh<hauteur ; nh++)
{
- qb = newQuad (tab_edge [nh+nv0], tab_edge [nh+nv1],
+ qb = newQuad (tab_edge [nh+nv0], tab_edge [nh+nv1],
tab_edge [nh+nv2], tab_edge [nh+nv3]);
qc = tab_quad [nh + nv0];
qd = tab_quad [nh + nv2];
qe = tab_quad [nh + nv1];
qf = tab_quad [nh + nv3];
-// *** tab_hexa [nh*hauteur + nro] = newHexa (qa, qb, qc, qd, qe, qf); Abu
+// *** tab_hexa [nh*hauteur + nro] = newHexa (qa, qb, qc, qd, qe, qf); Abu
tab_hexa [nro*hauteur + nh] = newHexa (qa, qb, qc, qd, qe, qf);
ker_hquad.push_back (qb);
qa = qb;
return HOK;
}
// ====================================================== updateMatrix
-void Elements::updateMatrix (int hauteur)
+void Elements::updateMatrix (int nh)
{
if (revo_lution)
{
- gen_matrix.defRotation (revo_center, revo_axis, gen_values[hauteur]);
+ gen_matrix.defRotation (revo_center, revo_axe, gen_values[nh]);
+ cum_matrix.defRotation (revo_center, revo_axe, cum_values[nh]);
}
- else if (prism_vec)
+ else
{
- double h0 = hauteur>0 ? gen_values[hauteur-1] : 0;
- double dh = gen_values[hauteur] - h0;
- Real3 decal;
+ double hauteur = (nh+1)*prism_len;
+ double dh = prism_len;
+ if (prism_vec)
+ {
+ if (under_v6)
+ {
+ hauteur = gen_values[nh];
+ dh = nh>0 ? hauteur-gen_values[nh-1] : hauteur;
+ }
+ else
+ {
+ hauteur = cum_values[nh];
+ dh = gen_values[nh];
+ }
+ }
+ Real3 trans, decal;
for (int nc=dir_x ; nc<=dir_z ; nc++)
- decal [nc] = prism_dir [nc]*dh;
+ {
+ decal [nc] = prism_dir [nc]*dh;
+ trans [nc] = prism_dir [nc]*hauteur;
+ }
+
gen_matrix.defTranslation (decal);
+ cum_matrix.defTranslation (trans);
}
}
// ====================================================== endPrism
void Elements::endPrism ()
{
+ closeShape();
+
int nbelts = ker_hquad.size();
for (int nro=0 ; nro<nbelts ; nro++)
tab_quad.push_back (ker_hquad[nro]);
nbr_quads = tab_quad.size ();
nbr_vertex = tab_vertex.size ();
}
+// ====================================================== getShape
+NewShape* Elements::getShape()
+{
+ if (grid_geom==NULL)
+ {
+ string name = "extrud_" + el_name;
+ grid_geom = el_root->addShape (name.c_str(), SH_EXTRUD);
+ grid_geom -> openShape();
+ }
+
+ return grid_geom;
+}
+// ====================================================== closeShape
+void Elements::closeShape()
+{
+ if (grid_geom==NULL)
+ return;
+
+ grid_geom -> closeShape();
+ grid_geom = NULL;
+}
END_NAMESPACE_HEXA