// C++ : Gestion des cylindres croises
-// Copyright (C) 2009-2013 CEA/DEN, EDF R&D
+// Copyright (C) 2009-2024 CEA, EDF
//
// 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 "HexNewShape.hxx"
+#include "HexBiCylinderShape.hxx"
-static bool db = false;
+/// static const double UnSur2pi = DEMI/M_PI;
BEGIN_NAMESPACE_HEXA
+static bool db = on_debug ();
static const double cos45 = cos (M_PI/4);
-
-
-void geom_define_line (string& brep);
-void geom_asso_point (double angle, Vertex* node);
-
-void geom_create_circle (double* milieu, double rayon, double* normale,
- double* base, string& brep);
-int geom_create_cylcyl (double* borig, double* bnorm, double* bbase,
- double bray, double bhaut,
- double* sorig, double* snorm, double* sbase,
- double sray, double shaut);
-int geom_asso_cylcyl (Edge* edge);
+#define NameOf(x) ((x)?x->getName():"Null")
// ====================================================== Constructeur
BiCylinder::BiCylinder (Document* doc)
: Elements (doc)
{
cyl_fill = false;
- cross_cyl1 = NULL;
- cross_cyl1 = NULL;
- cross_cyl2 = NULL;
- grid_type = GR_BICYL;
+ grid_type = GR_BIPIPE;
at_right = at_left = true;
tab_hexa.push_back (NULL);
tab_edge.push_back (NULL);
tab_vertex.push_back (NULL);
nbr_vertex = nbr_edges = nbr_quads = nbr_hexas = 1;
+
+ grid_geom = NULL;
}
// ====================================================== getHexaIJK
Hexa* BiCylinder::getHexaIJK (int cyl, int nx, int ny, int nz)
if (tab_vertex[nro] == NULL)
printf ("NULL\n");
else
- printf ("%s = (%g, %g, %g)\n", tab_vertex[nro]->getName(),
+ printf ("%s = (%g, %g, %g)\n", NameOf (tab_vertex[nro]),
tab_vertex[nro]->getX(), tab_vertex[nro]->getY(),
tab_vertex[nro]->getZ());
}
if (node == NULL)
printf ("NULL\n");
else
- printf ("%s = (%g, %g, %g)\n", node->getName(), node->getX(),
+ printf ("%s = (%g, %g, %g)\n", NameOf (node), node->getX(),
node->getY(), node->getZ());
}
nbr_vertex ++;
printf ("pres_edge [%d,%d, %d,%d,%d] = E%d = ", cyl, dir, nx, ny, nz,
nbr_edges);
- printf ("%s ((%s, %s)) = (%s,%s)\n", edge->getName(),
- edge->getVertex(0)->getName(), edge->getVertex(1)->getName(),
- v1->getName(), v2->getName());
- if (NOT edge->definedBy (v1,v2))
+ printf ("%s ((%s, %s)) = (%s,%s)\n", NameOf(edge),
+ NameOf(edge->getVertex(0)), NameOf(edge->getVertex(1)),
+ NameOf(v1), NameOf(v2));
+ if (edge->definedBy (v1,v2))
+ return tab_edge [nro];
+ else
printf (" **** Incoherence !!\n");
}
- return tab_edge [nro];
}
if (v1==NULL || v2==NULL)
return NULL;
Edge* edge = findEdge (v1, v2);
- if (edge==NULL)
- edge = newEdge (v1, v2);
+ if (edge!=NULL)
+ return edge;
+ edge = newEdge (v1, v2);
map_edge [key] = nbr_edges;
tab_edge.push_back (edge);
+ nbr_edges ++;
if (db)
{
printf (" tab_edge [%d,%d, %d,%d,%d] = E%d = ", cyl, dir, nx, ny, nz,
- nbr_edges);
+ nbr_edges-1);
if (edge == NULL)
printf ("NULL\n");
else
- printf ("%s = (%s, %s)\n", edge->getName(),
- edge->getVertex(0)->getName(),
- edge->getVertex(1)->getName());
+ printf ("%s = (%s, %s)\n", NameOf(edge),
+ NameOf(edge->getVertex(0)),
+ NameOf(edge->getVertex(1)));
}
- nbr_edges ++;
return edge;
}
// ====================================================== addQuad
int key = getKey (cyl, dir, nx, ny, nz);
int nro = map_quad [key];
if (nro>0)
- return tab_quad [nro];
+ {
+ if (db)
+ {
+ Quad* quad = tab_quad [nro];
+
+ printf ("pres_quad [%d,%d, %d,%d,%d] = Q%d = ", cyl, dir, nx, ny, nz,
+ nbr_quads);
+ printf ("%s (%s, %s, %s, %s)\n", NameOf(quad),
+ NameOf(e1), NameOf(e2), NameOf(e3), NameOf(e4));
+ if (quad->definedBy (e1,e3))
+ return tab_quad [nro];
+ else
+ {
+ printf (" **** Incoherence !!\n");
+ printf ("%s = (%s, %s, %s, %s)\n", NameOf(quad),
+ NameOf(quad->getEdge(0)), NameOf(quad->getEdge(1)),
+ NameOf(quad->getEdge(1)), NameOf(quad->getEdge(2)));
+ printf ("%s = (%s, %s, %s, %s)\n", NameOf(quad),
+ NameOf(quad->getVertex(0)), NameOf(quad->getVertex(1)),
+ NameOf(quad->getVertex(1)), NameOf(quad->getVertex(2)));
+ }
+ }
+ }
+
+ Quad* quad = Elements::findQuad (e1, e2);
+ if (quad!=NULL)
+ return quad;
- Quad* quad = newQuad (e1, e2, e3, e4);
+ quad = newQuad (e1, e2, e3, e4);
map_quad [key] = nbr_quads;
tab_quad.push_back (quad);
+ nbr_quads ++;
if (db)
{
printf (" tab_quad [%d,%d, %d,%d,%d] = Q%d = ", cyl, dir, nx, ny, nz,
- nbr_quads);
+ nbr_quads-1);
if (quad == NULL)
printf ("NULL\n");
else
- printf ("%s = (%s, %s, %s, %s)\n", quad->getName(),
- quad->getEdge(0)->getName(), quad->getEdge(1)->getName(),
- quad->getEdge(2)->getName(), quad->getEdge(3)->getName());
+ printf ("%s = (%s, %s, %s, %s)\n", NameOf(quad),
+ NameOf(quad->getEdge(0)), NameOf(quad->getEdge(1)),
+ NameOf(quad->getEdge(2)), NameOf(quad->getEdge(3)));
}
- nbr_quads ++;
return quad;
}
// ====================================================== addHexa
int nro = map_hexa [key];
if (nro>0)
{
- printf (" tab_hexa [%d, %d,%d,%d] = H%d est deja la\n ",
+ if (tab_hexa [nro]->definedBy (q1,q2))
+ {
+ printf (" tab_hexa [%d, %d,%d,%d] = H%d est deja la\n ",
cyl, nx, ny, nz, nbr_hexas);
- return tab_hexa [nro];
+ return tab_hexa [nro];
+ }
+ printf (" Incoherence : tab_hexa [%d, %d,%d,%d] = H%d = %s\n ",
+ cyl, nx, ny, nz, nbr_hexas, tab_hexa [nro]->getName ());
}
Hexa* hexa = newHexa (q1, q2, q3, q4, q5, q6);
if (hexa == NULL)
printf ("NULL\n");
else
- printf ("%s = (%s, %s, %s, %s, %s, %s)\n", hexa->getName(),
- hexa->getQuad(0)->getName(), hexa->getQuad(1)->getName(),
- hexa->getQuad(2)->getName(), hexa->getQuad(3)->getName(),
- hexa->getQuad(4)->getName(), hexa->getQuad(5)->getName());
+ printf ("%s = (%s, %s, %s, %s, %s, %s)\n", NameOf(hexa),
+ NameOf(hexa->getQuad(0)), NameOf(hexa->getQuad(1)),
+ NameOf(hexa->getQuad(2)), NameOf(hexa->getQuad(3)),
+ NameOf(hexa->getQuad(4)), NameOf(hexa->getQuad(5)));
}
nbr_hexas ++;
return hexa;
if (db)
printf ("findVertex [Big,%d,%d,%d] = V%d = '%d' = %s\n",
nx, ny, nz, nro, it_map->first,
- tab_vertex[nro]->getName());
+ NameOf(tab_vertex[nro]));
return elt;
}
}
if (db)
printf ("findEdge [%d, %d,%d,%d] = E%d = '%d' = %s\n",
dir, nx, ny, nz, nro, it_map->first,
- elt->getName());
+ NameOf(elt));
return elt;
}
}
if (db)
printf ("findQuad [%d, %d,%d,%d] = E%d = '%d' = %s\n",
dir, nx, ny, nz, nro, it_map->first,
- elt->getName());
+ NameOf(elt));
return elt;
}
}
map_quad [key] = nro;
if (db)
printf ("findQuad [%d, %d,%d,%d] = E%d = '%d' = %s\n",
- dir, nx, ny, nz, nro, it_map->first, elt->getName());
+ dir, nx, ny, nz, nro, it_map->first, NameOf(elt));
+ return elt;
+ }
+ }
+ fatal_error ("HexBiCylinder : Quad non trouve");
+ return NULL;
+}
+/*******************************************
+// ====================================================== findQuad
+Quad* BiCylinder::findQuad (Edge* e1, Edge* e2)
+{
+ for (it_map=map_quad.begin() ; it_map!=map_quad.end() ; ++it_map)
+ {
+ int nro = it_map->second;
+ Quad* elt = tab_quad[nro];
+ if (elt != NULL && elt->definedBy (e1, e2))
+ {
+ int key = getKey (CylBig, dir, nx, ny, nz);
+ map_quad [key] = nro;
+ if (db)
+ printf ("findQuad [%d, %d,%d,%d] = E%d = '%d' = %s\n",
+ dir, nx, ny, nz, nro, it_map->first, NameOf(elt));
return elt;
}
}
fatal_error ("HexBiCylinder : Quad non trouve");
return NULL;
}
+ ********************************************/
// ====================================================== findHexa
Hexa* BiCylinder::findHexa (Quad* q1, Quad* q2, int nx, int ny, int nz)
{
fatal_error ("HexBiCylinder : Hexa non trouve");
return NULL;
}
-// ====================================================== crossCylinders
-int BiCylinder::crossCylinders (Cylinder* lun, Cylinder* lautre)
+// ====================================================== createLittlePipe
+void BiCylinder::createLittlePipe ()
{
- if (lun->getRadius() < lautre->getRadius())
+ Real lg = cross_hauteur[CylSmall];
+ Real h1 = calc_distance (cross_center, cross_orismall);
+ Real h2 = cross_rayext[CylBig]*cos45;
+ Real h3 = (cross_rayext[CylBig] - cross_rayint[CylBig])*cos45;
+ Real h5 = lg-h1;
+
+// double t_haut [NbrVslices] = { -h1, -h2, -h2+h3, h2-h3, h2, lg-h1 };
+
+ if (at_left)
{
- cross_cyl1 = lun;
- cross_cyl2 = lautre;
+ addSlice (CylSmall, 0, 0, -h1, cross_rayint [CylSmall]);
+ addSlice (CylSmall, 1, 0, -h1, cross_rayext [CylSmall]);
+ addSlice (CylSmall, 1, 1, -h2, cross_rayext [CylSmall]);
+ addSlice (CylSmall, 0, 2, h3-h2, cross_rayint [CylSmall]);
}
- else
+ else if (cyl_fill)
{
- cross_cyl1 = lautre;
- cross_cyl2 = lun;
- }
+ addSlice (CylSmall, 0, 1, -h2, cross_rayint [CylSmall]);
- int ier = cross_cyl2->interCylinder (cross_cyl1, at_left, at_right,
- cross_center);
- if (ier != HOK)
- return ier;
-
- cyl_fill = false;
- cross_rayext [CylSmall] = cross_cyl1->getRadius ();
- cross_rayext [CylBig] = cross_cyl2->getRadius ();
- cross_rayint [CylSmall] = cross_rayext [CylSmall] / 2;
- cross_rayint [CylBig ] = cross_rayext [CylBig ] / 2;
- cross_hauteur [CylSmall] = cross_cyl1->getHeight ();
- cross_hauteur [CylBig] = cross_cyl2->getHeight ();
-
- createLittleCyl ();
- createBigCyl ();
- adjustLittleSlice (1, 1);
- adjustLittleSlice (0, 2);
- adjustLittleSlice (0, 3);
- adjustLittleSlice (1, 4);
-
- Vector* iprim = new Vector (cross_cyl1->getDirection());
- Vector* kprim = new Vector (cross_cyl2->getDirection());
- Vector* jprim = new Vector (kprim);
-
- iprim->renormer ();
- kprim->renormer ();
- jprim->vectoriel (kprim, iprim);
+ addSlice (CylSmall, 1, 1, -h2, cross_rayext [CylSmall]);
+ addSlice (CylSmall, 0, 2, h3-h2, cross_rayint [CylSmall]);
+ }
- // transfoVertices (cross_center, iprim, jprim, kprim);
- // Real3 snorm, bnorm;
- // iprim->getCoord (snorm);
- // kprim->getCoord (bnorm);
+ addSlice (CylSmall, 2, 1, -h2, cross_rayext [CylBig]);
+ addSlice (CylSmall, 1, 2, h3-h2, cross_rayint [CylBig]);
+ addSlice (CylSmall, 1, 3, h2-h3, cross_rayint [CylBig]);
+ addSlice (CylSmall, 2, 4, h2, cross_rayext [CylBig]);
- Real3 bnorm = {0, 0, 1};
- Real3 snorm = {1, 0, 0};
- assoCylinders (snorm, bnorm);
+ if (at_right)
+ {
+ addSlice (CylSmall, 0, 3, h2-h3, cross_rayint [CylSmall]);
+ addSlice (CylSmall, 1, 4, h2, cross_rayext [CylSmall]);
+ addSlice (CylSmall, 0, 5, h5, cross_rayint [CylSmall]);
+ addSlice (CylSmall, 1, 5, h5, cross_rayext [CylSmall]);
+ }
+ else if (cyl_fill)
+ {
+ addSlice (CylSmall, 0, 3, h2-h3, cross_rayint [CylSmall]);
+ addSlice (CylSmall, 1, 4, h2, cross_rayext [CylSmall]);
-/*********************************************
+ addSlice (CylSmall, 0, 4, h2, cross_rayint [CylSmall]);
+ }
+ //--------------------------------- Remplissage
+ // ka kb kc kd
if (at_left)
{
- assoIntersection (NxExt, 1, snorm, bnorm);
- if (grid_type == GR_BIPIPE)
+ fillSlice (CylSmall, 0,0, 0,2, 1,1, 1,0);
+ fillSlice (CylSmall, 0,2, 1,2, 2,1, 1,1);
+ if (cyl_fill)
{
- assoIntersection (NxInt, 2, snorm, bnorm);
+ addCube (CylSmall, 0,0,2);
}
}
-
- if (at_right)
+ else if (cyl_fill)
{
- assoIntersection (NxExt, NbrSlices1-1, snorm, bnorm);
- if (grid_type == GR_BIPIPE)
- {
- assoIntersection (NxInt, NbrSlices1-2, snorm, bnorm);
- }
+ fillSlice (CylSmall, 0,1, 0,2, 1,2, 2,1);
+ addCube (CylSmall, 1,0);
}
+ else
+ addCube ( CylSmall, 1, 2);
- ******************************************* */
- assoResiduelle ();
- return HOK;
-}
-// ====================================================== createLittleCyl
-void BiCylinder::createLittleCyl ()
-{
- Real3 base;
- Vertex* vbase = cross_cyl1->getBase();
-
- Real lg = cross_hauteur[CylSmall];
- Real h1 = calc_distance (cross_center, vbase->getPoint (base));
- Real h2 = cross_rayext[CylBig]*cos45;
- Real h3 = (cross_rayext[CylBig] - cross_rayint[CylBig])*cos45;
-
- double t_haut [NbrVslices] = { -h1, -h2, -h2+h3, h2-h3, h2, lg-h1 };
+ fillSlice (CylSmall, 1,2, 1,3, 2,4, 2,1, true);
+ if (cyl_fill)
+ {
+ addCube (CylSmall, 2,0, 3);
+ fillSlice (CylSmall, 0,2, 0,3, 1,3, 1,2);
+ }
- for (int nk=0; nk<NbrVslices ; nk++)
- {
- switch (nk)
- {
- case 0 : case 5 :
- addSlice (CylSmall, 0, nk, t_haut[nk], cross_rayint [CylSmall]);
- addSlice (CylSmall, 1, nk, t_haut[nk], cross_rayext [CylSmall]);
- break;
- case 1 : case 4 :
- addSlice (CylSmall, 1, nk, t_haut[nk], cross_rayext [CylSmall]);
- addSlice (CylSmall, 2, nk, t_haut[nk], cross_rayext [CylBig]);
- break;
- case 2 : case 3 :
- addSlice (CylSmall, 0, nk, t_haut[nk], cross_rayint [CylSmall]);
- addSlice (CylSmall, 1, nk, t_haut[nk], cross_rayint [CylBig]);
- break;
- }
- }
- // ka kb kc kd
- fillSlice (CylSmall, 0,0, 0,2, 1,1, 1,0);
- fillSlice (CylSmall, 0,2, 1,2, 2,1, 1,1); // OK
- // fillSlice (CylSmall, 1,1, 0,2, 1,2, 2,1); // Test
- fillSlice (CylSmall, 1,2, 1,3, 2,4, 2,1, true);
- fillSlice (CylSmall, 1,3, 0,3, 1,4, 2,4);
- fillSlice (CylSmall, 0,3, 0,5, 1,5, 1,4);
- // fillSmallCyl ();
+ if (at_right)
+ {
+ fillSlice (CylSmall, 1,3, 0,3, 1,4, 2,4);
+ fillSlice (CylSmall, 0,3, 0,5, 1,5, 1,4);
+ if (cyl_fill)
+ addCube (CylSmall, 3, 0, 5);
+ }
+ else if (cyl_fill)
+ {
+ fillSlice (CylSmall, 0,3, 0,4, 2,4, 1,3);
+ addCube (CylSmall, 3, 0);
+ }
+ else
+ addCube ( CylSmall, 3, 1);
}
-// ====================================================== createBigCyl
-void BiCylinder::createBigCyl ()
+// ====================================================== createBigPipe
+void BiCylinder::createBigPipe ()
{
const int cyl = CylBig;
- Real3 base;
- Vertex* vbase = cross_cyl2->getBase();
Real lg = cross_hauteur[cyl];
Real rext = cross_rayext [cyl];
Real rint = cross_rayint [cyl];
- Real h1 = calc_distance (cross_center, vbase->getPoint (base));
+
+ Real h1 = calc_distance (cross_center, cross_oribig);
Real h2 = rext * cos45;
Real h3 = lg - h1;
Real dh = (rext - rint)*cos45;
addSlice (CylBig, 1, 3, h3, rext);
// A B C D
+ if (cyl_fill)
+ addCube (CylBig, 0);
fillSlice (CylBig, 0,0, 0,1, 1,1, 1,0);
- fillSlice (CylBig, 0,2, 0,3, 1,3, 1,2);
+ if (cyl_fill)
+ addCube (CylBig, 2);
+ fillSlice (CylBig, 0,2, 0,3, 1,3, 1,2);
}
// ====================================================== adjustLittleSlice
void BiCylinder::adjustLittleSlice (int ni, int nk)
double prayon = cross_rayext[CylSmall];
if (ni==0)
{
- grayon2 = cross_rayint[CylBig] * cross_rayint[CylBig];
+ if (NOT cyl_fill)
+ grayon2 = cross_rayint[CylBig] * cross_rayint[CylBig];
prayon = cross_rayint[CylSmall];
}
addVertex (hauteur, px, py, CylSmall, ni, nj, nk);
}
}
+// ====================================================== addCarre
+void BiCylinder::addCarre (int cyl, int nk, double hauteur, double rayon,
+ bool find)
+{
+ for (int nj=0 ; nj<NbrCotes ; nj++)
+ {
+ double theta = getAngle (nj);
+ double px = rayon*cos(theta);
+ double py = rayon*sin(theta);
+ if (find)
+ findVertex (px, py, hauteur, 0, nj, nk);
+ else if (cyl==CylBig)
+ addVertex (px, py, hauteur, CylBig, 0, nj, nk);
+ else // CylSmall
+ addVertex (hauteur, px, py, cyl, 0, nj, nk);
+ }
+}
+// ====================================================== adjustLittleSquare
+void BiCylinder::adjustLittleSquare (int nk)
+{
+ double grayon2 = cross_rayext[CylBig] * cross_rayext[CylBig];
+
+ for (int nj=0 ; nj<NbrCotes ; nj++)
+ {
+ Vertex* node = getVertexIJK (CylSmall, 0, nj, nk);
+ double py = node->getY();
+ double px = sqrt (grayon2 - py*py);
+
+ if (node->getX() > 0) node->setX ( px);
+ else node->setX (-px);
+ }
+}
+// ====================================================== addCube
+void BiCylinder::addCube (int cyl, int nk0, int nvi0, int k1)
+{
+/* *****************************************************************
+ H=bed +----bd-----+ bdf=G......0
+ /| /|
+ be | B bf |
+ / | / |
+ E=bce +----bc-----+...|...bcf=F
+ | de D | df
+ | E | | F | J
+ ce | C cf | ^
+ D=ade...|...+----ad-|---+ adf=C.....3 | I
+ | / | / | /
+ | ae A | af | /
+ |/ |/ |/
+ A=ace +----ac-----+ acf=B.....2 +-----> K
+
+ ***************************************************************** */
+ enum { nj0, nj1, nj2, nj3, ni0=0 };
+ const int nk1 = k1>0 ? k1 : nk0+1;
+ const int nvi1 = nvi0==1 ? 2
+ : nvi0==2 ? 1 : 0;
+
+ // La face F
+ Vertex* vacf = getVertexIJK (cyl, nvi1, nj2, nk1);
+ Vertex* vadf = getVertexIJK (cyl, nvi1, nj3, nk1);
+ Vertex* vbdf = getVertexIJK (cyl, nvi1, nj0, nk1);
+ Vertex* vbcf = getVertexIJK (cyl, nvi1, nj1, nk1);
+
+ Edge* eaf = addEdge (vacf, vadf, cyl, dir_y, ni0, nj2, nk1);
+ Edge* edf = addEdge (vadf, vbdf, cyl, dir_y, ni0, nj3, nk1);
+ Edge* ebf = addEdge (vbcf, vbdf, cyl, dir_y, ni0, nj0, nk1);
+ Edge* ecf = addEdge (vacf, vbcf, cyl, dir_y, ni0, nj1, nk1);
+
+ Quad* qf = addQuad (eaf, edf, ebf, ecf, cyl, dir_z, ni0 , nj0, nk1);
+
+ Vertex* vace = getVertexIJK (cyl, nvi0, nj2, nk0);
+ Vertex* vade = getVertexIJK (cyl, nvi0, nj3, nk0);
+ Vertex* vbde = getVertexIJK (cyl, nvi0, nj0, nk0);
+ Vertex* vbce = getVertexIJK (cyl, nvi0, nj1, nk0);
+
+ Edge* eac = addEdge (vace, vacf, cyl, dir_z, ni0, nj2, nk0);
+ Edge* ead = addEdge (vade, vadf, cyl, dir_z, ni0, nj3, nk0);
+ Edge* ebd = addEdge (vbde, vbdf, cyl, dir_z, ni0, nj0, nk0);
+ Edge* ebc = addEdge (vbce, vbcf, cyl, dir_z, ni0, nj1, nk0);
+
+ Edge* eae = addEdge (vace, vade, cyl, dir_y, ni0, nj2, nk0);
+ Edge* ede = addEdge (vade, vbde, cyl, dir_y, ni0, nj3, nk0);
+ Edge* ebe = addEdge (vbce, vbde, cyl, dir_y, ni0, nj0, nk0);
+ Edge* ece = addEdge (vace, vbce, cyl, dir_y, ni0, nj1, nk0);
+
+// Quad* qe = addQuad (eae, ede, ebe, ece, cyl, dir_z, ni0 , nj0, nk0);
+ Quad* qe = addQuad (eae, ece, ebe, ede, cyl, dir_z, ni0 , nj0, nk0);
+
+// Quad* qa = addQuad (eac, eaf, ead, eae, cyl, dir_y, ni0 , nj2, nk0);
+ Quad* qa = addQuad (eac, eae, ead, eaf, cyl, dir_y, ni0 , nj2, nk0);
+// Quad* qd = addQuad (ead, edf, ebd, ede, cyl, dir_x, ni0 , nj3, nk0);
+ Quad* qd = addQuad (ead, ede, ebd, edf, cyl, dir_x, ni0 , nj3, nk0);
+ Quad* qb = addQuad (ebc, ebf, ebd, ebe, cyl, dir_y, ni0 , nj0, nk0);
+ Quad* qc = addQuad (eac, ecf, ebc, ece, cyl, dir_x, ni0 , nj1, nk0);
+
+ addHexa (qa, qb, qc, qd, qe, qf, cyl, ni0, nj0, nk0);
+}
// ====================================================== fillSlice
/* *****************************************************************
- H=bed +----bd-----+ bdf=G
+ H=bde +----bd-----+ bdf=G
/| /|
be | B bf |
/ | / |
Vertex* vbdf = getVertexIJK (cyl, nic, nj1, nkc);
Vertex* vbde = getVertexIJK (cyl, nid, nj1, nkd);
-/* *******************
- PutName (vace);
- PutName (vacf);
- PutName (vadf);
- PutName (vade);
- PutName (vbce);
- PutName (vbcf);
- PutName (vbdf);
- PutName (vbde);
- ******************* */
-
Edge* eac = addEdge (vace, vacf, cyl, dir_z, nia, nj0, nka);
- Edge* ead = addEdge (vade, vadf, cyl, dir_z, nid, nj0, nkd);
- Edge* eae = addEdge (vace, vade, cyl, dir_x, nia, nj0, nka);
+ Edge* ead = addEdge (vadf, vade, cyl, dir_z, nid, nj0, nkd);
+ Edge* eae = addEdge (vade, vace, cyl, dir_x, nia, nj0, nka);
Edge* eaf = addEdge (vacf, vadf, cyl, dir_x, nib, nj0, nkb);
Edge* ebc = addEdge (vbce, vbcf, cyl, dir_z, nia, nj1, nka);
- Edge* ebd = addEdge (vbde, vbdf, cyl, dir_z, nid, nj1, nkd);
- Edge* ebe = addEdge (vbce, vbde, cyl, dir_x, nia, nj1, nka);
+ Edge* ebd = addEdge (vbdf, vbde, cyl, dir_z, nid, nj1, nkd);
+ Edge* ebe = addEdge (vbde, vbce, cyl, dir_x, nia, nj1, nka);
Edge* ebf = addEdge (vbcf, vbdf, cyl, dir_x, nib, nj1, nkb);
Edge* ece = addEdge (vace, vbce, cyl, dir_y, nia, nj0, nka);
Edge* edf = addEdge (vadf, vbdf, cyl, dir_y, nic, nj0, nkc);
Edge* ede = addEdge (vade, vbde, cyl, dir_y, nid, nj0, nkd);
- Quad* qa = addQuad (eac, eaf, ead, eae, cyl, dir_y, nia , nj0, nka);
- Quad* qb = addQuad (ebc, ebf, ebd, ebe, cyl, dir_y, nia , nj1, nka);
- Quad* qc = addQuad (eac, ecf, ebc, ece, cyl, dir_x, nia , nj0, nka);
- Quad* qd = addQuad (ead, edf, ebd, ede, cyl, dir_x, nid , nj0, nkd);
- Quad* qe = addQuad (eae, ede, ebe, ece, cyl, dir_z, nia , nj0, nka);
- Quad* qf = addQuad (eaf, edf, ebf, ecf, cyl, dir_z, nib , nj0, nkb);
+ // Quad* qa = addQuad (eac, eaf, ead, eae, cyl, dir_y, nia+1, nj0, nka);
+ Quad* qa = addQuad (eac, eae, ead, eaf, cyl, dir_y, nia+1, nj0, nka);
+ Quad* qb = addQuad (ebc, ebf, ebd, ebe, cyl, dir_y, nia+1, nj1, nka);
+
+ Quad* qc = addQuad (eac, ecf, ebc, ece, cyl, dir_x, nia+1, nj0, nka);
+ // Quad* qd = addQuad (ead, edf, ebd, ede, cyl, dir_x, nid+1, nj0, nkd);
+ Quad* qd = addQuad (ead, ede, ebd, edf, cyl, dir_x, nid+1, nj0, nkd);
+
+ // Quad* qe = addQuad (eae, ede, ebe, ece, cyl, dir_z, nia+1, nj0, nka);
+ Quad* qe = addQuad (eae, ece, ebe, ede, cyl, dir_z, nia+1, nj0, nka);
+ Quad* qf = addQuad (eaf, edf, ebf, ecf, cyl, dir_z, nib+1, nj0, nkb);
- addHexa (qa, qb, qc, qd, qe, qf, cyl, nia, nj0, nka);
+ addHexa (qa, qb, qc, qd, qe, qf, cyl, nia+1, nj0, nka);
}
}
// ====================================================== assoCylinders
-void BiCylinder::assoCylinders (double* snormal, double* gnormal)
+void BiCylinder::assoCylinders ()
{
- assoSlice (CylSmall, 0, 0, snormal);
- assoSlice (CylSmall, 1, 0, snormal);
- assoSlice (CylSmall, 0, 5, snormal);
- assoSlice (CylSmall, 1, 5, snormal);
+ char name [12];
+ sprintf (name, "grid_%02d", el_id);
+ grid_geom = el_root->addShape (name, SH_INTER);
+ grid_geom -> openShape();
+
+ int s_kmax = 5;
+ int imin = cyl_fill ? 1 : 0;
+ int g_ksize = 4;
+
+ for (int ni=imin ; ni<2 ; ni++)
+ {
+ assoSlice (CylSmall, ni, 0, cross_dirsmall);
+ assoSlice (CylSmall, ni, s_kmax, cross_dirsmall);
+
+ for (int nk=0 ; nk<g_ksize ; nk++)
+ assoSlice (CylBig, ni, nk, cross_dirbig);
+ }
- for (int nk=0 ; nk<4 ; nk++)
- for (int ni=0 ; ni<2 ; ni++)
- assoSlice (CylBig, ni, nk, gnormal);
- assoIntersection (1, 1, snormal, gnormal);
- assoIntersection (0, 2, snormal, gnormal);
- assoIntersection (0, 3, snormal, gnormal);
- assoIntersection (1, 4, snormal, gnormal);
+ if (at_left)
+ {
+ assoIntersection (1, 1);
+ if (NOT cyl_fill)
+ assoIntersection (0, 2);
+ }
+
+ if (at_right)
+ {
+ assoIntersection (1, 4);
+ if (NOT cyl_fill)
+ assoIntersection (0, 3);
+ }
+ grid_geom -> closeShape();
}
// ====================================================== assoSlice
void BiCylinder::assoSlice (int cyl, int nx, int nzs, double* normal)
{
- Real3 center, pnt1, pnt2, vbase;
- string brep;
+ Real3 center, pnt0, pnt1, pnt2, vbase;
int ny0 = 0;
- int nyd = NbrCotes/2;
+ int ny1 = 1;
+ int ny2 = NbrCotes/2;
Vertex* v0 = getVertexIJK (cyl, nx, ny0 , nzs);
- Vertex* vd = getVertexIJK (cyl, nx, nyd , nzs);
+ Vertex* v1 = getVertexIJK (cyl, nx, ny1 , nzs);
+ Vertex* v2 = getVertexIJK (cyl, nx, ny2 , nzs);
- if (vd==NULL || v0==NULL)
+ if (v0==NULL || v1==NULL || v2==NULL)
return;
- v0->getPoint (pnt1);
- vd->getPoint (pnt2);
+ v0->getPoint (pnt0);
+ v1->getPoint (pnt1);
+ v2->getPoint (pnt2);
double rayon = 0;
for (int nro=0 ; nro<DIM3 ; nro++)
{
- center[nro] = (pnt1[nro] + pnt2[nro])/2;
- vbase [nro] = (pnt1[nro] - center[nro]);
- rayon += vbase [nro]*vbase[nro];
+ center[nro] = (pnt0[nro] + pnt2[nro])/2;
+ vbase [nro] = (pnt0[nro] - pnt1[nro]);
+ rayon += carre (center [nro] - pnt0 [nro]);
}
- rayon = sqrt (rayon);
-
- PutCoord (pnt1);
- PutCoord (pnt2);
- PutCoord (vbase);
- PutCoord (center);
- PutData (rayon);
-
- geom_create_circle (center, rayon, normal, vbase, brep);
- geom_define_line (brep); // pour geom_asso_point
+ int subid = grid_geom->addCircle (center, sqrt(rayon), normal, vbase);
for (int ny=0 ; ny<NbrCotes ; ny++)
{
- assoArc (cyl, nx, ny, nzs, brep, rayon);
+ int ny1 = (ny+1) MODULO NbrCotes;
+ Vertex* va = getVertexIJK (cyl, nx, ny , nzs);
+ Vertex* vb = getVertexIJK (cyl, nx, ny1, nzs);
+
+ assoArc (cyl, ny, va, vb, subid);
}
}
// ===================================================== assoArc
-void BiCylinder::assoArc (int cyl, int nx, int ny, int nz, string& brep,
- double rayon)
+void BiCylinder::assoArc (int cyl, int ny, Vertex* v1, Vertex* v2, int subid)
{
- static const double Alpha = 1.0 / NbrCotes;
- static const double Theta = 2*M_PI/ NbrCotes;
- int ny1 = (ny+1) MODULO NbrCotes;
-
- Edge* edge = getEdgeJ (cyl, nx, ny, nz);
- Vertex* node0 = getVertexIJK (cyl, nx, ny, nz);
- Vertex* node1 = getVertexIJK (cyl, nx, ny1, nz);
-
- if (db)
- printf ("AssoArc : Edge = %s = (%s,%s) -> (%s,%s)\n", edge->getName(),
- edge->getVertex(0)->getName(), edge->getVertex(1)->getName(),
- node0->getName(), node1->getName());
+ const double Decal = 1.0 / NbrCotes;
+ const double Start = Decal / 2;
+
+ int ny2 = ny+1;
+ double posit = Start + ny*Decal;
+ Edge* edge = findEdge (v1, v2);
+ if (edge==NULL)
+ return;
+ // Vertex
+ grid_geom->addAssociation (v1, subid, posit);
+
+ if (ny2 < NbrCotes)
+ {
+ grid_geom->addAssociation (edge, subid, posit, posit+Decal);
+ }
+ else
+ {
+ grid_geom->addAssociation (edge, subid, posit, 1.0);
+ grid_geom->addAssociation (edge, subid, 0, Start);
+ }
+}
+// ===================================================== assoArc
+void BiCylinder::assoArc (int cyl, int nx, int ny, int nz, int subid)
+{
+ const double Decal = 1.0 / NbrCotes;
+ const double Start = Decal / 2;
- // Shape* shape = new Shape (brep);
+ int ny2 = ny+1;
+ double posit = Start + ny*Decal;
+ Edge* edge = getEdgeJ (cyl, nx, ny, nz);
+ if (edge==NULL)
+ return;
- // shape->setBounds (ny*Alpha, (ny+1)*Alpha);
- // edge ->addAssociation (shape);
+ // Vertex
+ Vertex* node = getVertexIJK (cyl, nx, ny, nz);
+ grid_geom->addAssociation (node, subid, posit);
- geom_asso_point ( ny *Theta*rayon, node0);
- geom_asso_point ((ny+1)*Theta*rayon, node1);
+ if (ny2 < NbrCotes)
+ {
+ grid_geom->addAssociation (edge, subid, posit, posit+Decal);
+ // node = getVertexIJK (cyl, nx, ny2, nz);
+ // grid_geom->addAssociation (node, subid, angle2*UnSur2pi);
+ }
+ else
+ {
+ grid_geom->addAssociation (edge, subid, posit, 1.0);
+ grid_geom->addAssociation (edge, subid, 0, Start);
+ }
}
// ===================================================== assoIntersection
-int BiCylinder::assoIntersection (int nxs, int nzs, double* snorm,
- double* bnorm)
+int BiCylinder::assoIntersection (int nxs, int nzs)
{
- Real3 X_center = {0, 0, 0}; // provisoire, a la place de cross_center
+ int ier = HOK;
Real3 pse, psw, sorig, sbase;
Real3 pbe, pbw, borig, bbase;
- string brep;
+ std::string brep;
int ny0 = 0;
int nyd = NbrCotes/2;
int MiddleSlice1 = 3;
int nz = nzs < MiddleSlice1 ? 0 : 5;
- getVertexIJK (CylSmall, nxs, ny0 , nz)->getPoint (pse);
- getVertexIJK (CylSmall, nxs, nyd , nz)->getPoint (psw);
- getVertexIJK (CylBig, nxs, ny0 , 0) ->getPoint (pbe);
- getVertexIJK (CylBig, nxs, nyd , 0) ->getPoint (pbw);
+ Vertex* vse = getVertexIJK (CylSmall, nxs, ny0 , nz);
+ Vertex* vsw = getVertexIJK (CylSmall, nxs, nyd , nz);
+ Vertex* vbe = getVertexIJK (CylBig, nxs, ny0 , 0);
+ Vertex* vbw = getVertexIJK (CylBig, nxs, nyd , 0);
+
+ if (vse==NULL || vsw==NULL || vbe==NULL || vbw==NULL)
+ return HERR;
+
+ vse->getPoint (pse);
+ vsw->getPoint (psw);
+ vbe->getPoint (pbe);
+ vbw->getPoint (pbw);
double srayon = calc_distance (pse, psw)/2;
double brayon = calc_distance (pbe, pbw)/2;
calc_vecteur (psw, pse, sbase);
calc_vecteur (pbw, pbe, bbase);
- double shaut = calc_distance (X_center, sorig);
- double bhaut = calc_distance (X_center, borig)*2;
- double* orig = nzs < MiddleSlice1 ? sorig : X_center; // Pb orientation
+ double shaut = calc_distance (cross_center, sorig);
+ double bhaut = calc_distance (cross_center, borig)*2;
+ double* orig = nzs < MiddleSlice1 ? sorig : cross_center; // Pb orientation
- if (db)
- {
- PutCoord (borig);
- PutCoord (sorig);
- PutCoord (orig);
- PutData (nz);
- PutData (srayon);
- PutData (brayon);
- }
-
- int ier = geom_create_cylcyl (borig, bnorm, bbase, brayon, bhaut,
- orig, snorm, sbase, srayon, shaut);
+ BiCylinderShape bicyl_shape (el_root);
+ ier = bicyl_shape.defineCyls (borig, cross_dirbig, bbase, brayon, bhaut,
+ orig, cross_dirsmall, sbase, srayon, shaut);
if (ier != HOK)
return ier;
Vertex* node = getVertexIJK (CylSmall, nxs, ny, nzs);
if (node!=NULL)
node->clearAssociation ();
+ // Edge* edge = getEdgeJ (CylSmall, nxs, ny, nzs);
+ // if (edge!=NULL) edge->clearAssociation ();
}
for (int ny=0 ; ny<NbrCotes ; ny++)
{
- Edge* edge = getEdgeJ (CylSmall, nxs, ny, nzs);
- geom_asso_cylcyl (edge);
+ int ny1 = (ny+1) MODULO QUAD4;
+ Vertex* node0 = getVertexIJK (CylSmall, nxs, ny, nzs);
+ Vertex* node1 = getVertexIJK (CylSmall, nxs, ny1, nzs);
+ // Edge* edge = getEdgeJ (CylSmall, nxs, ny, nzs);
+ Edge* edge = findEdge (node0, node1);
+ bicyl_shape.associate (edge);
}
return HOK;
}
-// ======================================================== test_bicylinder
-BiCylinder* test_bicylinder (Document* docu, int option)
+// ====================================================== makeCylinders
+int BiCylinder::makeCylinders (Vertex* ori1, double* vz1, double r1, double h1,
+ Vertex* ori2, double* vz2, double r2, double h2)
+{
+ grid_type = GR_BICYL;
+ int ier = makePipes (ori1,vz1, r1/2, r1,h1, ori2,vz2, r2/2, r2,h2, true);
+ return ier;
+}
+// ====================================================== makePipes
+int BiCylinder::makePipes (Vertex* ori1, double* vz1, double rint1,
+ double rext1, double h1, Vertex* ori2, double* vz2,
+ double rint2, double rext2, double h2, bool fill)
{
- Vertex* ori1 = docu->addVertex ( 0,0,0);
- Vertex* ori2 = docu->addVertex (-5,0,5);
- Vector* vz = docu->addVector ( 0,0,1);
- Vector* vx = docu->addVector ( 1,0,0);
+ cyl_fill = fill;
+ cross_hauteur [CylBig] = h2;
+ cross_rayext [CylBig ] = rext2;
+ cross_rayint [CylBig ] = rint2;
+
+ cross_hauteur [CylSmall] = h1;
+ cross_rayext [CylSmall] = rext1;
+ cross_rayint [CylSmall] = rint1;
+
+ ori1->getPoint (cross_orismall);
+ ori2->getPoint (cross_oribig );
+ copy_vecteur (vz2, cross_dirbig);
+ copy_vecteur (vz1, cross_dirsmall);
+
+// Intersection des 2 cylindres :
+ int ier = checkInter (cross_oribig, vz2, rext2, h2,
+ cross_orismall, vz1, rext1, h1,
+ cross_center, at_left, at_right);
+ if (ier != HOK)
+ return ier;
- double r1 = 2*sqrt (2.0);
- double r2 = 3*sqrt (2.0);
- double l2 = 10;
- double l1 = 10;
+ createLittlePipe ();
+ createBigPipe ();
- Cylinder* cyl1 = docu->addCylinder (ori1, vz, r1, l1);
- Cylinder* cyl2 = docu->addCylinder (ori2, vx, r2, l2);
+ if (at_left)
+ {
+ adjustLittleSlice (1, 1);
+ adjustLittleSlice (0, 2);
+ }
+ else if (cyl_fill)
+ {
+ adjustLittleSlice (0, 1);
+ }
- BiCylinder* grid = new BiCylinder (docu);
- grid->crossCylinders (cyl1, cyl2);
- return grid;
+ if (at_right)
+ {
+ adjustLittleSlice (0, 3);
+ adjustLittleSlice (1, 4);
+ }
+ else if (cyl_fill)
+ {
+ adjustLittleSlice (0, 4);
+ }
+
+
+ transfoVertices (cross_center, cross_dirsmall, cross_dirbig);
+ assoCylinders ();
+
+ // assoResiduelle ();
+ // el_root->reorderQuads ();
+ return HOK;
}
END_NAMESPACE_HEXA
