src/TEST_CPP/test_hexa1.cxx

   1
   2 // C++ : Tests unitaires
   3
   4 // Copyright (C) 2009-2012  CEA/DEN, EDF R&D
   5 //
   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.
  10 //
  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.
  15 //
  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
  19 //
  20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  21 //
  22
  23 #include "test_unit.hxx"
  24
  25 #include "test_unit.hxx"
  26
  27 #include "Hex.hxx"
  28 #include "HexDocument.hxx"
  29 #include "HexElements.hxx"
  30 #include "HexCrossElements.hxx"
  31
  32 #include "HexHexa.hxx"
  33 #include "HexQuad.hxx"
  34 #include "HexEdge.hxx"
  35 #include "HexVertex.hxx"
  36
  37 #include "HexPropagation.hxx"
  38 #include "HexLaw.hxx"
  39 #include "HexMatrix.hxx"
  40 #include "HexCramer.hxx"
  41 #include "HexGroup.hxx"
  42
  43
  44 #include <cstdlib>
  45
  46 static int nbr_vtk = 0;
  47 static cpchar case_name = "hexa";
  48 static Hex::Document*   docu = NULL;
  49
  50 // ======================================================== print_propagations
  51 void print_propagations (Hex::Document* doc)
  52 {
  53    int nb = doc->countPropagation ();
  54    HexDisplay (nb);
  55    for (int nro=0 ; nro<nb ; nro++)
  56        {
  57        Hex::Propagation*  prop  = doc ->getPropagation (nro);
  58        const Hex::Edges&  table = prop->getEdges ();
  59        printf (" ____________________________________ Prop nro %d\n", nro);
  60        for (int ned=0 ; ned<(int)table.size() ; ned++)
  61            {
  62            bool way = table [ned]->getWay ();
  63            Hex::Edge*   edge = table [ned];
  64            Hex::Vertex* v0   = edge->getVertex (0);
  65            Hex::Vertex* v1   = edge->getVertex (1);
  66
  67            if (way)
  68               {
  69               printf ("     (");
  70               v0->printName (", ");
  71               v1->printName (")\n");
  72               }
  73           else
  74               {
  75               v1->printName (", ");
  76               v0->printName (")\n");
  77               printf ("     (");
  78               }
  79            }
  80        }
  81 }
  82 // ======================================================== save_vtk
  83 void save_vtk ()
  84 {
  85    if (docu==NULL)
  86       return;
  87
  88    docu->saveVtk (case_name, nbr_vtk);
  89 }
  90 // ======================================================== remove_hexa
  91 void remove_hexa (Hex::Hexa* hexa)
  92 {
  93    if (hexa==NULL)
  94       return;
  95
  96    hexa->remove();
  97
  98    if (docu==NULL)
  99       return;
 100
 101    docu->saveVtk (case_name, nbr_vtk);
 102 }
 103 // ======================================================== test_sphere
 104 int test_sphere (int nbargs, cpchar tabargs[])
 105 {
 106    Hex::Hex mon_ex;
 107    Hex::Document* doc = mon_ex.addDocument ();
 108    Hex::Vertex*  orig = doc->addVertex (0,0,0);
 109
 110    int    ncouches = 0;
 111    double k      = 1;
 112    double rayon  = 1;
 113    Hex::Elements* sphere = doc->makeSpherical (orig, rayon, ncouches, k);
 114
 115    if (ncouches>0)
 116       {
 117    for (int nc=0 ; nc <= ncouches ; nc++)
 118        {
 119        Hex::Hexa* cell = sphere->getStrate (nc, Hex::Q_A);
 120        if (nc==0)
 121           cell->dumpFull();
 122        cell->remove ();
 123        // sphere->getStrate (nc, Hex::Q_A)->remove ();
 124        // sphere->getStrate (nc, Hex::Q_B)->remove ();
 125        }
 126        }
 127
 128    sphere->saveVtk ("sphere.vtk");
 129    return HOK;
 130 }
 131 // ======================================================== test_cartesi1
 132 int test_cartesi1 ()
 133 {
 134    const int size_x = 15;
 135    const int size_y = 12;
 136    const int size_z = 8;
 137
 138    Hex::Hex mon_ex;
 139    Hex::Document* doc = mon_ex.addDocument ();
 140    Hex::Vertex* orig = doc->addVertex (0,0,0);
 141
 142    Hex::Vector*   dir  = doc->addVector (1,1,1);
 143    Hex::Elements* grid = doc->makeCartesian (orig, dir, size_x,size_y,size_z);
 144
 145    //  print_propagations (doc);
 146
 147    for (int nz=0; nz<size_z ; nz++)
 148        for (int ny=nz+1; ny<size_y-nz-1 ; ny++)
 149            for (int nx=nz+1; nx<size_x-nz-1 ; nx++)
 150                {
 151                Hex::Hexa* cell = grid->getHexaIJK (nx, ny, nz);
 152                cell->remove ();
 153                }
 154
 155    doc->setLevel (1);
 156    print_propagations (doc);
 157    grid->saveVtk ("grid_cart.vtk");
 158
 159    // doc->dump ();
 160    return HOK;
 161 }
 162 // ======================================================== afficher
 163 #define Afficher(elt) afficher (#elt, elt)
 164 int afficher (cpchar nom, Hex::EltBase* elt)
 165 {
 166    if (elt==NULL)
 167       {
 168       printf (" .... %s = 0x0\n", nom);
 169       return HOK;
 170       }
 171
 172    printf (" .... %s = 0x%08lx = %03d\n", nom, (unsigned long) elt, elt->getId());
 173    return HOK;
 174 }
 175 // ======================================================== test_find
 176 int test_find ()
 177 {
 178    const int size_x = 2;
 179    const int size_y = 2;
 180    const int size_z = 2;
 181
 182    Hex::Hex mon_ex;
 183    Hex::Document* doc = mon_ex.addDocument ();
 184
 185    Hex::Vertex* orig = doc->addVertex (0,0,0);
 186    Hex::Vector* dir  = doc->addVector (1,1,1);
 187    Hex::Elements*  grid = doc->makeCartesian (orig, dir, size_x,size_y,size_z);
 188
 189    grid->saveVtk ("mini1.vtk");
 190    doc->dump ();
 191
 192    Hex::Vertex *v00, *v02, *v06, *v08, *v10, *v22, *v26;
 193
 194    Afficher ( v00 = doc->findVertex (0, 0, 0));
 195    Afficher ( v02 = doc->findVertex (1, 1, 0));
 196    Afficher ( v06 = doc->findVertex (1, 1, 1));
 197    Afficher ( v08 = doc->findVertex (2, 1, 0));
 198    Afficher ( v10 = doc->findVertex (2, 1, 1));
 199    Afficher ( v22 = doc->findVertex (2, 1, 2));
 200    Afficher ( v26 = doc->findVertex (2, 2, 2));
 201
 202    printf ("\n");
 203
 204    Afficher (doc->findEdge (v06, v10));
 205    Afficher (doc->findEdge (v10, v06));
 206    printf ("\n");
 207
 208    Afficher (doc->findQuad (v02, v10));
 209    Afficher (doc->findQuad (v06, v08));
 210    Afficher (doc->findQuad (v02, v06));
 211
 212    printf ("\n");
 213    Afficher (doc->findHexa (v00, v06));
 214    Afficher (doc->findHexa (v06, v26));
 215    Afficher (doc->findHexa (v26, v06));
 216
 217    return HOK;
 218 }
 219 // ======================================================== test_joint
 220 int test_joint (int nbargs, cpchar tabargs[])
 221 {
 222    const int dimx = 11;
 223    const int dimy = 11;
 224    const int dimz = 2;
 225
 226    Hex::Hex mon_ex;
 227    Hex::Document* doc = mon_ex.addDocument ();
 228
 229    Hex::Vertex* orig1 = doc->addVertex (0,0,0);
 230    Hex::Vector* dir   = doc->addVector (1,1,1);
 231
 232    Hex::Elements* grid1 = doc->makeCartesian   (orig1, dir, dimx,dimy,dimz);
 233
 234    Hex::Vertex*   orig2 = doc->addVertex  (dimx/2.0,0,8);
 235    Hex::Vector*   vectj = doc->addVector (0,1,0);
 236    Hex::Vector*   vecti = doc->addVector (1,0,0);
 237    Hex::Elements* grid2 = doc->makeCylindrical (orig2, vecti, vectj,
 238                    1, 180, 1,        dimz,dimy,dimx);
 239
 240    int mx = dimx/2;
 241    int my = dimy/2;
 242    Hex::Quad* prems = grid1->getQuadIJ (mx, my, dimz);
 243    Hex::Quad* cible = grid2->getQuadJK (dimz, mx, my);
 244
 245    Hex::Vertex* v1 = prems->getVertex (0);
 246    Hex::Vertex* v3 = prems->getVertex (1);
 247
 248    Hex::Vertex* v2 = cible->getVertex (1);
 249    Hex::Vertex* v4 = cible->getVertex (2);
 250
 251    // v1->setScalar (3);
 252    // v2->setScalar (3);
 253    // v3->setScalar (6);
 254    // v4->setScalar (6);
 255
 256    Hex::Quads liste;
 257
 258    liste.push_back (prems);
 259    for (int nx=0; nx<dimx; nx++)
 260        if (nx!=mx)
 261           liste.push_back (grid1->getQuadIJ (nx, my, dimz));
 262
 263    for (int ny=0; ny<dimy; ny++)
 264        if (ny!=my)
 265           liste.push_back (grid1->getQuadIJ (mx, ny, dimz));
 266
 267    doc->saveVtk ("joint1.vtk");
 268    const int hauteur = 5;
 269    Hex::Elements* joint = doc->joinQuads(liste, cible, v1,v2,v3,v4, hauteur);
 270    // for (int nh=0 ; nh<hauteur ; nh++) joint->getHexa(nh)->setScalar (5);
 271
 272    int nbr_joint_vertex =  joint->countVertex ();
 273    int nbr_surf_vertex  =  nbr_joint_vertex/(hauteur+1);
 274
 275    HexDisplay (nbr_joint_vertex);
 276    HexDisplay (nbr_surf_vertex);
 277
 278    int indice0 = joint->findVertex (v1);
 279    HexDisplay (indice0);
 280
 281    for (int nh=0 ; nh<nbr_surf_vertex ; nh++)
 282        joint->getVertex(nh)->setScalar (5);
 283
 284    for (int nh=0 ; nh<=hauteur ; nh++)
 285        joint->getVertex(nh*nbr_surf_vertex)->setScalar (3);
 286
 287    doc->saveVtk ("joint2.vtk");
 288    return HOK;
 289 }
 290 // ======================================================== test_prism
 291 int test_prism (int nbargs, cpchar tabargs[])
 292 {
 293    const int dimx = 11;
 294    const int dimy = 11;
 295    const int dimz = 2;
 296
 297    Hex::Hex mon_ex;
 298    Hex::Document* doc = mon_ex.addDocument ();
 299
 300    Hex::Vertex* orig1 = doc->addVertex ( 0,0,0);
 301    Hex::Vector* dir1  = doc->addVector ( 1,1,1);
 302    Hex::Vector* dir2  = doc->addVector ( 1,1,-1);
 303
 304    Hex::Elements* grid1 = doc->makeCartesian (orig1, dir1, dimx,dimy,dimz);
 305
 306    int mx = dimx/2;
 307    int my = dimy/2;
 308    Hex::Quads liste1, liste2;
 309
 310    liste1.push_back (grid1->getQuadIJ (mx, my, dimz));
 311    liste2.push_back (grid1->getQuadIJ (mx, my, 0));
 312    for (int nx=0; nx<dimx; nx++)
 313        if (nx!=mx)
 314           {
 315           liste1.push_back (grid1->getQuadIJ (nx, my, dimz));
 316           liste2.push_back (grid1->getQuadIJ (nx, my, 0));
 317           }
 318
 319    for (int ny=0; ny<dimy; ny++)
 320        if (ny!=my)
 321           {
 322           liste1.push_back (grid1->getQuadIJ (mx, ny, dimz));
 323           liste2.push_back (grid1->getQuadIJ (mx, ny, 0));
 324           }
 325
 326    Hex::RealVector tlen;
 327    double dh = 2;
 328    for (int nro=0; nro<5; nro++)
 329        {
 330        dh = 2*dh + 1;
 331        tlen.push_back (dh);
 332        }
 333
 334    const int nbiter = 5;
 335    doc->saveVtk ("prisme1.vtk");
 336    Hex::Elements* prisme2 = doc->prismQuads    (liste2, dir2, nbiter);
 337    doc->saveVtk ("prisme2.vtk");
 338
 339    Hex::Elements* prisme1 = doc->prismQuadsVec (liste1, dir1, tlen, 0);
 340
 341    PutData (liste1.size());
 342    PutData (tlen.size());
 343    PutData (prisme1->countHexa());
 344    PutData (prisme1->countQuad());
 345    PutData (prisme1->countEdge());
 346    PutData (prisme1->countVertex());
 347
 348    for (int nro=0 ; nro <nbiter  ; nro++)
 349        {
 350        Hex::Hexa* cell = prisme2-> getHexa (nbiter+nro);
 351        cell->setScalar (5);
 352        }
 353
 354    doc->saveVtk ("prisme3.vtk");
 355    return HOK;
 356 }
 357 // ======================================================== test_revolution9
 358 int test_revolution9 (int nbargs, cpchar tabargs[])
 359 {
 360    const int dimx = 11;
 361    const int dimy = 11;
 362    const int dimz = 2;
 363
 364    Hex::Hex mon_ex;
 365    Hex::Document* doc = mon_ex.addDocument ();
 366
 367    Hex::Vertex* orig1 = doc->addVertex (0,0,0);
 368    Hex::Vector* dir   = doc->addVector (1,1,1);
 369
 370    Hex::Elements* grid1 = doc->makeCartesian   (orig1, dir, dimx,dimy,dimz);
 371
 372    int mx = dimx/2;
 373    int my = dimy/2;
 374    Hex::Quad* prems = grid1->getQuadIJ (mx, my, dimz);
 375    Hex::Quads liste;
 376
 377    liste.push_back (prems);
 378    prems -> setScalar (5);
 379    for (int nx=0; nx<dimx; nx++)
 380        if (nx!=mx)
 381           {
 382           Hex::Quad* cell = grid1->getQuadIJ (nx, my, dimz);
 383           liste.push_back (cell);
 384           cell -> setScalar (5);
 385           }
 386
 387    for (int ny=0; ny<dimy; ny++)
 388        if (ny!=my)
 389           {
 390           Hex::Quad* cell = grid1->getQuadIJ (mx, ny, dimz);
 391           liste.push_back (cell);
 392           cell -> setScalar (5);
 393           }
 394
 395    Hex::Vertex* center = doc->addVertex (0, -10, 0);
 396    Hex::Vector* axis   = doc->addVector (1,   0, 0);
 397    Hex::RealVector  angles;
 398
 399    Hex::Vector*   dir1  = doc->addVector (10,0.3,0.3);
 400    Hex::Elements* grid2 = doc->makeCartesian (center, dir1, 1,1,1);
 401    Hex::Hexa*     cell  = grid2->getHexaIJK (0,0,0);
 402    cell->setScalar (5);
 403
 404    doc->saveVtk ("revolution1.vtk");
 405
 406    double alpha = 5;
 407    int niter    = 5;
 408    double coeff = 1.5;
 409    for (int na=0 ; na<niter ; na++)
 410        {
 411        angles.push_back (alpha);
 412        alpha *= coeff;
 413        }
 414    for (int na=1 ; na<niter ; na++)
 415        {
 416        alpha /= coeff;
 417        angles.push_back (alpha);
 418        }
 419
 420    Hex::Elements* bloc = doc->revolutionQuads  (liste, center, axis, angles);
 421    if (bloc != NULL)
 422        doc->saveVtk ("revolution2.vtk");
 423
 424    return HOK;
 425 }
 426 // ======================================================== test_revolution
 427 int test_revolution (int nbargs, cpchar tabargs[])
 428 {
 429    Hex::Hex mon_ex;
 430    Hex::Document* doc = mon_ex.addDocument ();
 431
 432    Hex::Vertex* ori = doc->addVertex (0,0,0);
 433    Hex::Vector* vx  = doc->addVector (1,0,0);
 434    Hex::Vector* vz  = doc->addVector (0,0,1);
 435
 436    int dr = 1;
 437    int da = 360;
 438    int dl = 1;
 439
 440    int nr = 1;
 441    int na = 6;
 442    int nl = 1;
 443
 444    Hex::Elements* grid = doc->makeCylindrical (ori, vx,vz, dr,da,dl,
 445                                                nr,na,nl, false);
 446
 447    Hex::Quads liste;
 448    for (int nx=0; nx<nr; nx++)
 449        for (int ny=0; ny<na; ny++)
 450            {
 451            Hex::Quad* cell = grid->getQuadIJ (nx, ny, nl);
 452            liste.push_back (cell);
 453            cell -> setScalar (5);
 454            }
 455
 456    Hex::Vertex* center = doc->addVertex (0, -10, 0);
 457    Hex::Vector* axis   = doc->addVector (1, 0, 0);
 458    Hex::RealVector  angles;
 459
 460    Hex::Vector*   dir1  = doc->addVector (10,0.3,0.3);
 461    Hex::Elements* grid2 = doc->makeCartesian (center, dir1, 1,1,1);
 462    Hex::Hexa*     cell  = grid2->getHexaIJK (0,0,0);
 463    cell->setScalar (5);
 464
 465    doc->saveVtk ("revolution1.vtk");
 466
 467    double alpha = 5;
 468    int niter    = 5;
 469    double coeff = 1.5;
 470    for (int na=0 ; na<niter ; na++)
 471        {
 472        angles.push_back (alpha);
 473        alpha *= coeff;
 474        }
 475    for (int na=1 ; na<niter ; na++)
 476        {
 477        alpha /= coeff;
 478        angles.push_back (alpha);
 479        }
 480
 481    Hex::Elements* bloc = doc->revolutionQuads  (liste, center, axis, angles);
 482    if (bloc != NULL)
 483        doc->saveVtk ("revolution2.vtk");
 484
 485    return HOK;
 486 }
 487 // ======================================================== test_coude
 488 int test_coude (int nbargs, cpchar tabargs[])
 489 {
 490 #if 0
 491    const int dimx = 11;
 492    const int dimy = 11;
 493    const int dimz = 2;
 494
 495    Hex::Hex mon_ex;
 496    Hex::Document* doc = mon_ex.addDocument ();
 497
 498    Hex::Vertex* orig = doc->addVertex (0,0,0);
 499    Hex::Vector* vx   = doc->addVector (1,0,0);
 500    Hex::Vector* vz   = doc->addVector (0,0,1);
 501
 502    //   grid1 = doc->makeCartesian   (orig1, dir, dimx,dimy,dimz);
 503    double dr = 1;
 504    int    dl = 5;
 505    int    nr = 4;
 506    int    na = 8;
 507
 508    Hex::Elements* grid1 = doc->makeCylindrical (orig1, vx,vz,dr,360, dl,
 509                                                 nr, 10, nl, false);
 510    int mx = dimx/2;
 511    int my = dimy/2;
 512    Hex::Quad* prems = grid1->getQuadIJ (mx, my, dimz);
 513    Hex::Quads liste;
 514
 515    liste.push_back (prems);
 516    prems -> setScalar (5);
 517    for (int nx=0; nx<dimx; nx++)
 518        if (nx!=mx)
 519           {
 520           Hex::Quad* cell = grid1->getQuadIJ (nx, my, dimz);
 521           liste.push_back (cell);
 522           cell -> setScalar (5);
 523           }
 524
 525    for (int ny=0; ny<dimy; ny++)
 526        if (ny!=my)
 527           {
 528           Hex::Quad* cell = grid1->getQuadIJ (mx, ny, dimz);
 529           liste.push_back (cell);
 530           cell -> setScalar (5);
 531           }
 532
 533
 534    Hex::Vertex* center = doc->addVertex (0, -10, 0);
 535    Hex::Vector* axis   = doc->addVector (1, 0, 0);
 536    Hex::RealVector  angles;
 537
 538    Hex::Vector*   dir1  = doc->addVector (10,0.3,0.3);
 539    Hex::Elements* grid2 = doc->makeCartesian (center, dir1, 1,1,1);
 540    Hex::Hexa*     cell  = grid2->getHexaIJK (0,0,0);
 541    cell->setScalar (5);
 542
 543    doc->saveVtk ("revolution1.vtk");
 544
 545    double alpha = 5;
 546    int niter    = 5;
 547    double coeff = 1.5;
 548    for (int na=0 ; na<niter ; na++)
 549        {
 550        angles.push_back (alpha);
 551        alpha *= coeff;
 552        }
 553    for (int na=1 ; na<niter ; na++)
 554        {
 555        alpha /= coeff;
 556        angles.push_back (alpha);
 557        }
 558
 559    Hex::Elements* bloc = doc->revolutionQuads  (liste, center, axis, angles);
 560    if (bloc != NULL)
 561        doc->saveVtk ("revolution2.vtk");
 562 #endif
 563    return HOK;
 564 }
 565 // ======================================================== test_count
 566 int test_count (int nbargs, cpchar tabargs[])
 567 {
 568    Hex::Hex mon_ex;
 569    Hex::Document* doc = mon_ex.addDocument ();
 570
 571    Hex::Vertex* orig1 = doc->addVertex ( 0, 0,0);
 572    Hex::Vector* vx    = doc->addVector (1,0,0);
 573    Hex::Vector* vz    = doc->addVector (0,0,1);
 574
 575    double dr = 1;
 576    double dl = 1;
 577    int    nr = 2;
 578    int    nl = 3;
 579
 580    Hex::Elements* c1;
 581
 582    c1 = doc->makeCylindrical (orig1, vx,vz,dr, 360, dl,nr, 10, nl, false);
 583
 584    HexDisplay (doc->countVertex ());
 585    HexDisplay (doc->countUsedVertex ());
 586    doc ->saveVtk ("hexa1.vtk");
 587
 588    return HOK;
 589 }
 590 // ======================================================== test_decoupage
 591 int test_decoupage (int nbargs, cpchar tabargs[])
 592 {
 593    const int size_x = 2;
 594    const int size_y = 1;
 595    const int size_z = 1;
 596
 597    Hex::Hex mon_ex;
 598    Hex::Document* doc = mon_ex.addDocument ();
 599
 600    Hex::Vertex* orig = doc->addVertex (0,0,0);
 601    Hex::Vector* dir  = doc->addVector (1,1,1);
 602
 603    Hex::Elements* grid  = doc->makeCartesian (orig, dir, size_x,size_y,size_z);
 604    Hex::Edge*     arete = grid->getEdgeK (0, 0, 0);
 605
 606    //  doc ->dump ();
 607    int nvtk=0;
 608    doc ->saveVtk ("decoupe", nvtk);
 609 /* Hex::Elements* grid2 = */  doc->cut (arete, 1);
 610
 611    ///  doc ->dump ();
 612    doc ->saveVtk ("decoupe", nvtk);
 613
 614    return HOK;
 615 }
 616 // ======================================================== test_gen_xml
 617 int test_gen_xml (int nbargs, cpchar tabargs[])
 618 {
 619    const int size_x = 2;
 620    const int size_y = 2;
 621    const int size_z = 2;
 622
 623    Hex::Hex mon_ex;
 624    Hex::Document* doc = mon_ex.addDocument ();
 625
 626    Hex::Vertex* orig = doc->addVertex (0,0,0);
 627    Hex::Vector* dir  = doc->addVector (1,1,1);
 628    Hex::Elements*  grid = doc->makeCartesian (orig, dir, size_x,size_y,size_z);
 629
 630    Hex::Hexa*   cell    = grid->getHexa (0);
 631    Hex::Quad*   face    = cell->getQuad (0);
 632    Hex::Edge*   arete   = cell->getEdge (0);
 633    Hex::Vertex* noeud   = cell->getVertex (0);
 634
 635    Hex::Law* law1 = doc->addLaw("loi1", 1);
 636    Hex::Law* law2 = doc->addLaw("loi2", 2);
 637    Hex::Law* law3 = doc->addLaw("loi3", 3);
 638
 639    law1->setKind (Hex::Uniform);
 640    law2->setKind (Hex::Arithmetic);
 641    law3->setKind (Hex::Geometric);
 642
 643    Hex::Propagation* prop1 = doc->getPropagation (0);
 644    Hex::Propagation* prop2 = doc->getPropagation (1);
 645    Hex::Propagation* prop3 = doc->getPropagation (2);
 646
 647    prop1->setLaw (law1);
 648    prop2->setLaw (law2);
 649    prop3->setLaw (law3);
 650
 651    prop1->setWay (true);
 652    prop2->setWay (false);
 653    prop3->setWay (true);
 654
 655    doc ->saveVtk ("mini.vtk");
 656    doc ->save ("Essai");
 657
 658    return HOK;
 659 }
 660 // ======================================================== test_string_xml
 661 int test_string_xml (int nbargs, cpchar tabargs[])
 662 {
 663    const int size_x = 2;
 664    const int size_y = 2;
 665    const int size_z = 2;
 666
 667    Hex::Hex mon_ex;
 668    Hex::Document* doc = mon_ex.addDocument ();
 669
 670    Hex::Vertex* orig = doc->addVertex (0,0,0);
 671    Hex::Vector* dir  = doc->addVector (1,1,1);
 672    Hex::Elements*  grid = doc->makeCartesian (orig, dir, size_x,size_y,size_z);
 673
 674    Hex::Hexa*   cell    = grid->getHexa (0);
 675    Hex::Quad*   face    = cell->getQuad (0);
 676    Hex::Edge*   arete   = cell->getEdge (0);
 677    Hex::Vertex* noeud   = cell->getVertex (0);
 678
 679    Hex::Law* law1 = doc->addLaw("loi1", 1);
 680    Hex::Law* law2 = doc->addLaw("loi2", 2);
 681    Hex::Law* law3 = doc->addLaw("loi3", 3);
 682
 683    law1->setKind (Hex::Uniform);
 684    law2->setKind (Hex::Arithmetic);
 685    law3->setKind (Hex::Geometric);
 686
 687    Hex::Propagation* prop1 = doc->getPropagation (0);
 688    Hex::Propagation* prop2 = doc->getPropagation (1);
 689    Hex::Propagation* prop3 = doc->getPropagation (2);
 690
 691    prop1->setLaw (law1);
 692    prop2->setLaw (law2);
 693    prop3->setLaw (law3);
 694
 695    prop1->setWay (true);
 696    prop2->setWay (false);
 697    prop3->setWay (true);
 698
 699    doc ->saveVtk ("mini.vtk");
 700    doc ->save ("Essai");
 701
 702    cpchar flux = doc ->getXml ();
 703    Hex::Document* docbis = mon_ex.addDocument ();
 704    docbis->setXml  (flux);
 705    docbis->saveVtk ("clone.vtk");
 706
 707    return HOK;
 708 }
 709 // ======================================================== test_relecture
 710 int test_relecture (int nbargs, cpchar tabargs[])
 711 {
 712    Hex::Hex mon_ex;
 713    cpchar nomdoc = "Essai";
 714    if (nbargs>1)
 715       nomdoc = tabargs[1];
 716    Hex::Document* doc = mon_ex.loadDocument (nomdoc);
 717
 718 /*********************
 719    Hex::Vertex* v4 = doc->findVertex (80.0, 0.0,  0.0);
 720    Hex::Vertex* v5 = doc->findVertex (80.0, 0.0, 40.0);
 721    Hex::Edge*   e4 = doc->findEdge   (v4, v5);
 722
 723    HexDump (v4);
 724    HexDump (v5);
 725    HexDump (e4);
 726
 727    e4->setScalar (5);
 728 ***********************/
 729
 730    doc ->dump ();
 731    doc ->saveVtk ("restore.vtk");
 732    doc ->save ("restore");
 733
 734    // doc ->reorderFaces ();
 735    // doc ->dump ();
 736
 737    // Hex::Elements* grid2 = doc->cut (e4, 2);
 738    return HOK;
 739 }
 740 // ======================================================== test_clone
 741 int test_clone ()
 742 {
 743    const int size_x = 2;
 744    const int size_y = 2;
 745    const int size_z = 2;
 746
 747    Hex::Hex mon_ex;
 748    Hex::Document* doc = mon_ex.addDocument ();
 749
 750    Hex::Vertex* orig = doc->addVertex (0,0,0);
 751    Hex::Vector* dir  = doc->addVector (1,1,1);
 752    Hex::Elements* grid = doc->makeCartesian (orig, dir, size_x,size_y,size_z);
 753
 754    Hex::Vector*   bond  = doc->addVector (0, 0, 7);
 755    Hex::Elements* grid2 = doc->makeTranslation (grid, bond);
 756
 757    doc ->saveVtk ("clonage.vtk");
 758    doc ->dump();
 759
 760    HexDump (grid2->getHexa      (1));
 761    HexDump (grid2->getHexaIJK   (1,1,1));
 762    HexDump (grid2->getVertexIJK (1,1,1));
 763
 764    return HOK;
 765 }
 766 // ======================================================== test_separ
 767 int test_separ ()
 768 {
 769    const int size_x = 2;
 770    const int size_y = 2;
 771    const int size_z = 2;
 772
 773    Hex::Hex mon_ex;
 774    Hex::Document* doc = mon_ex.addDocument ();
 775
 776    Hex::Vertex* orig = doc->addVertex (0,0,0);
 777    Hex::Vector* dir  = doc->addVector (1,1,1);
 778    //  Hex::Elements*  grid =
 779    doc->makeCartesian (orig, dir, size_x,size_y,size_z);
 780
 781    doc ->saveVtk ("separ.vtk");
 782    doc ->dump();
 783
 784    return HOK;
 785 }
 786 // ======================================================== test_shperical
 787 int test_spherical (int nbargs, const char* tabargs[])
 788 {
 789    Hex::Hex mon_ex;
 790    Hex::Document* doc = mon_ex.addDocument ();
 791
 792    Hex::Vertex* orig = doc->addVertex (0,0,0);
 793    double       rayon  = 1;
 794    int          nbr    = 3;
 795
 796    Hex::Elements* grid = doc->makeSpherical (orig, rayon, nbr);
 797
 798    int nbhexas = grid->countHexa ();
 799    HexDisplay (nbhexas);
 800    for (int nro=0 ; nro<nbhexas ; nro +=3)
 801        grid->getHexa(nro)->remove();
 802    HexDisplay (doc->countHexa ());
 803    doc->saveVtk ("shperical.vtk");
 804    // doc->dump ();
 805
 806    return HOK;
 807 }
 808 // ================================================== test_grille_cyl
 809 int test_grille_cyl (int nbargs, cpchar tabargs[])
 810 {
 811    Hex::Hex mon_ex;
 812    Hex::Document* doc = mon_ex.addDocument ();
 813
 814    Hex::Vertex* orig1 = doc->addVertex ( 0, 0,0);
 815    Hex::Vertex* orig2 = doc->addVertex (10, 0,0);
 816
 817    Hex::Vertex* orig3 = doc->addVertex ( 0,10,0);
 818    Hex::Vertex* orig4 = doc->addVertex (10,10,0);
 819
 820    Hex::Vertex* orig5 = doc->addVertex ( 0,20,0);
 821    Hex::Vertex* orig6 = doc->addVertex (10,20,0);
 822
 823    Hex::Vector* vz = doc->addVector (0,0,1);
 824    Hex::Vector* vx = doc->addVector (1,0,0);
 825
 826    double dr = 1;
 827    double dl = 1;
 828    int    nr = 2;
 829    int    nl = 3;
 830
 831    Hex::Elements *c1, *c2, *c3, *c4, *c5, *c6;
 832
 833    c1 = doc->makeCylindrical (orig1, vx,vz,dr, 360, dl,nr, 4, nl, true);
 834    c2 = doc->makeCylindrical (orig2, vx,vz,dr, 360, dl,nr, 8, nl, true);
 835    c3 = doc->makeCylindrical (orig3, vx,vz,dr, 270, dl,nr, 8, nl, true);
 836    c4 = doc->makeCylindrical (orig4, vx,vz,dr, 270, dl,nr, 7, nl, true);
 837    c5 = doc->makeCylindrical (orig5, vx,vz,dr, 360, dl,nr, 5, nl, true);
 838    c6 = doc->makeCylindrical (orig6, vx,vz,dr, 360, dl,nr, 6, nl, true);
 839
 840    int base2 = nr*nl*8;
 841    c2->getHexa(base2 + 0)->setScalar (5);
 842    c2->getHexa(base2 + 1)->setScalar (5);
 843    c2->getHexa(base2 + 2)->setScalar (5);
 844
 845    doc->saveVtk ("cylindres.vtk");
 846    // doc->dump ();
 847
 848    return HOK;
 849 }
 850 // ===================================================== test_cylindrical
 851 int test_cylindrical (int nbargs, cpchar tabargs[])
 852 {
 853    cpchar fic_vtk = "cylindre";
 854
 855    Hex::Hex mon_ex;
 856    Hex::Document* doc = mon_ex.addDocument ();
 857
 858    Hex::Vertex* orig = doc->addVertex (0, 0, 0);
 859    Hex::Vector* vx   = doc->addVector (1, 0, 0);
 860    Hex::Vector* vz   = doc->addVector (0, 0, 1);
 861
 862    double dr = 1;
 863    double da = 360;
 864    double dl = 1;
 865
 866    int nr = 1;
 867    int na = 8;
 868    int nl = 2;
 869
 870    if (nbargs>1)
 871       {
 872       na = atoi (tabargs[1]);
 873       HexDisplay (na);
 874       if (na <= 2)
 875           na = 8;
 876       }
 877
 878    if (nbargs>2)
 879       {
 880       da = atof (tabargs[2]);
 881       HexDisplay (da);
 882       }
 883
 884    // Hex::Cylinder* cyl  = doc->addCylinder   (orig, vz, nr, nl);
 885    // Hex::Elements* grid = doc->makeCylinder (cyl, vx, nr, na, nl);
 886    doc->makeCylindrical (orig,vx, vz, dr,da,dl, nr,na, nl, true);
 887    doc ->saveVtk (fic_vtk, na);
 888    return HOK;
 889 }
 890 // ===================================================== test_cylinder
 891 int test_cylinder (int nbargs, cpchar tabargs[])
 892 {
 893    int    nvtk    = 1;
 894    cpchar fic_vtk = "cylindre";
 895
 896    Hex::Hex mon_ex;
 897    Hex::Document* doc = mon_ex.addDocument ();
 898
 899    Hex::Vertex* orig = doc->addVertex (0, 0, 0);
 900    Hex::Vector* vx   = doc->addVector (1, 0, 0);
 901    Hex::Vector* vz   = doc->addVector (0, 0, 1);
 902
 903    double rayon   = 10;
 904    double hauteur = 6;
 905
 906    int nr = 2;
 907    int na = 8;
 908    int nl = 5;
 909
 910    Hex::Cylinder* cyl  = doc->addCylinder   (orig, vz, rayon, hauteur);
 911    doc->makeCylinder (cyl, vx, nr, na, nl);
 912    doc ->saveVtk (fic_vtk, nvtk);
 913    return HOK;
 914 }
 915 // ===================================================== test_xml_cylinder
 916 int test_xml_cylinder (int nbargs, cpchar tabargs[])
 917 {
 918    int    nvtk    = 0;
 919    cpchar fic_vtk = "cylindre";
 920
 921    Hex::Hex mon_ex;
 922    Hex::Document* doc = mon_ex.addDocument ();
 923
 924    Hex::Vertex* orig1 = doc->addVertex (0, 0,0);
 925    Hex::Vertex* orig2 = doc->addVertex (50,0,0);
 926    Hex::Vector* vz    = doc->addVector (0,0,1);
 927    Hex::Vector* vx    = doc->addVector (1,0,0);
 928
 929    vx->setName ("vx");
 930    vz->setName ("vz");
 931    orig1->setName ("orig1");
 932    orig2->setName ("orig2");
 933
 934    int nr  = 4;
 935    int nri = 3;
 936    int nre = nr;
 937    int na = 9;
 938    int nl = 5;
 939
 940    Hex::Cylinder* cyl  = doc->addCylinder   (orig1, vz, nr, nl);
 941    Hex::Pipe*     pipe = doc->addPipe       (orig2, vz, nri, nre, nl);
 942
 943    Hex::Elements* grid = doc->makeCylinder (cyl,  vx, nr, na, nl);
 944    doc ->saveVtk (fic_vtk, nvtk);
 945
 946    Hex::Group* groupe = doc->addGroup ("GroupeAMA", Hex::HexaCell);
 947    groupe->addElement (grid->getHexaIJK (0,0,0));
 948    groupe->addElement (grid->getHexaIJK (1,0,0));
 949    groupe->addElement (grid->getHexaIJK (0,1,0));
 950    groupe->addElement (grid->getHexaIJK (1,1,0));
 951    groupe->addElement (grid->getHexaIJK (2,1,0));
 952
 953    grid->getHexaIJK  (0,0,0)->setName ("Hexa0");
 954    grid->getQuadIJ   (0,0,0)->setName ("QuadIJ0");
 955    grid->getEdgeK    (0,0,0)->setName ("EdgeK0");
 956
 957    doc->makePipe     (pipe, vx, nr, na, nl);
 958    doc ->saveVtk (fic_vtk, nvtk);
 959    doc->save ("cylindre");
 960
 961    return HOK;
 962 }
 963 // ===================================================== test_pipe
 964 int test_pipe (int nbargs, cpchar tabargs[])
 965 {
 966    int    nvtk    = 0;
 967    cpchar fic_vtk = "cylindre";
 968
 969    Hex::Hex mon_ex;
 970    Hex::Document* doc = mon_ex.addDocument ();
 971
 972    Hex::Vertex* orig1 = doc->addVertex (0, 0,0);
 973    Hex::Vector* vx    = doc->addVector (1,0,0);
 974    Hex::Vector* vy    = doc->addVector (0,1,0);
 975
 976    int nr  = 1;
 977    int nri = 1;
 978    int nre = 2;
 979    int na = 2;
 980    int nl = 1;
 981
 982    Hex::Pipe*  pipe = doc->addPipe  (orig1, vx, nri, nre, nl);
 983    doc->makePipe     (pipe, vy, nr, na, nl);
 984    doc ->saveVtk (fic_vtk, nvtk);
 985
 986    return HOK;
 987 }
 988 // ======================================================== del_hexa
 989 void del_hexa (Hex::CrossElements* gr, int cyl, int ni, int nj, int nk, int dr)
 990 {
 991    Hex::Hexa* hexa = gr->getHexaIJK (cyl, ni, nj, nk);
 992    if (hexa!=NULL)
 993       {
 994       hexa->remove ();
 995       if (dr>1)
 996           save_vtk ();
 997       }
 998 }
 999 // ======================================================== del_tranche
1000 int del_tranche (Hex::CrossElements* grid, int cyl, int ni, int nk, int dr=1)
1001 {
1002    for (int nj = 0 ; nj < 8 ; nj++)
1003         del_hexa (grid, cyl, ni, nj, nk, dr);
1004
1005    if (dr==1)
1006       save_vtk ();
1007    printf ("del_tranche (g=%d, i=%d, k=%d) : fic = %d\n",
1008                          cyl, ni, nk, nbr_vtk-1);
1009    return nbr_vtk;
1010 }
1011 // ======================================================== test_joint2
1012 int test_joint2 (int nbargs, cpchar tabargs[])
1013 {
1014    Hex::Hex mon_ex;
1015    docu = mon_ex.addDocument ();
1016    case_name = "pb_joint";
1017
1018    Hex::Vector* vx    = docu->addVector (1, 0, 0);
1019    Hex::Vector* vz    = docu->addVector (0, 0, 1);
1020    Hex::Vertex* hori  = docu->addVertex (0, 0, 0);
1021
1022    double da = 360;
1023    double dr = 2;
1024    double dl = 1;
1025    int    nr = 1;
1026    int    na = 8;
1027    int    nl = 1;
1028    bool   fill = true;
1029
1030    Hex::Elements *bgrid=NULL, *hgrid=NULL;
1031
1032    hgrid = docu->makeCylindrical (hori, vx,vz, dr,da,dl, nr,na,nl, fill);
1033    docu->dump ();
1034    save_vtk ();
1035
1036    Hex::Vertex* bori  = docu->addVertex (0, 0, -5);
1037    bgrid = docu->makeCylindrical (bori, vx,vz, dr,da,dl, nr,na,nl, fill);
1038    save_vtk ();
1039
1040    Hex::Quads qsource, qdest;
1041    printf (" Source = ");
1042    for (int ny=0 ; ny< na ; ny++)
1043        {
1044        Hex::Quad* quad = hgrid->getQuadIJ (0, ny, 0);
1045        PrintName (quad);
1046        qsource.push_back (quad);
1047        }
1048    printf ("\n          ");
1049    for (int ny=0 ; ny<4 ; ny++)
1050        {
1051        Hex::Quad* quad = hgrid->getKerHQuad (ny);
1052        PrintName (quad);
1053        qsource.push_back (quad);
1054        }
1055    printf ("\n");
1056
1057
1058    printf (" Cible  = ");
1059    for (int ny=0 ; ny< na ; ny++)
1060        {
1061        Hex::Quad* quad = bgrid->getQuadIJ (0, ny, 1);
1062        PrintName (quad);
1063        qdest.push_back (quad);
1064        }
1065    printf ("\n          ");
1066    for (int ny=0 ; ny<4 ; ny++)
1067        {
1068        Hex::Quad* quad = bgrid->getKerHQuad (ny+4);
1069        PrintName (quad);
1070        qdest.push_back (quad);
1071        }
1072    printf ("\n");
1073    docu ->setLevel (1);
1074    Hex::Quad*   cible = bgrid->getQuadIJ    (0, 1, 1);
1075    Hex::Vertex* vc1   = bgrid->getVertexIJK (0, 1, 1);
1076    Hex::Vertex* vc2   = bgrid->getVertexIJK (1, 1, 1);
1077
1078    Hex::Vertex* vs1 = hgrid->getVertexIJK (0, 0, 0);
1079    Hex::Vertex* vs2 = hgrid->getVertexIJK (1, 0, 0);
1080
1081    docu->joinQuads (qsource, cible, vs1, vc1, vs2, vc2, 1);
1082    save_vtk ();
1083
1084    return HOK;
1085 }
1086 // ======================================================== test_croix
1087 int test_croix (int nbargs, cpchar tabargs[])
1088 {
1089    Hex::Hex mon_ex;
1090    docu = mon_ex.addDocument ();
1091
1092    Hex::Vertex* ori1 = docu->addVertex ( 0,0,0);
1093    Hex::Vertex* ori2 = docu->addVertex (-5,0,5);
1094    Hex::Vector* vz   = docu->addVector ( 0,0,1);
1095    Hex::Vector* vx   = docu->addVector ( 1,0,0);
1096
1097    double r1 = 2;
1098    double r2 = 3;
1099    double l2 = 10;
1100    double l1 = 10;
1101
1102    Hex::Cylinder*      cyl1 = docu->addCylinder (ori1, vz, r1, l1);
1103    Hex::Cylinder*      cyl2 = docu->addCylinder (ori2, vx, r2, l2);
1104    Hex::CrossElements* grid = docu->makeCylinders (cyl1, cyl2);
1105
1106    case_name = "croix";
1107    save_vtk ();
1108
1109    del_tranche (grid, 0, 1, 0);
1110    del_tranche (grid, 0, 1, 5);
1111
1112    del_tranche (grid, 1, 1, 0);
1113    del_tranche (grid, 1, 1, 3);
1114
1115    del_tranche (grid, 1, 0, 0);
1116    del_tranche (grid, 1, 0, 3);
1117                                   // Le trognon
1118    del_tranche (grid, 0, 0, 0);
1119    del_tranche (grid, 0, 0, 5);
1120    del_tranche (grid, 0, 0, 1);
1121    del_tranche (grid, 0, 0, 2);
1122    del_tranche (grid, 0, 0, 3);
1123    del_tranche (grid, 0, 0, 4);
1124                                   // Partie critique
1125
1126    del_tranche (grid, 1, 1, 1, 2);
1127
1128    del_tranche (grid, 0, 1, 1, 2);
1129    del_tranche (grid, 0, 1, 4, 2);
1130    del_tranche (grid, 0, 1, 3, 2);
1131    del_tranche (grid, 0, 1, 2, 2);
1132
1133    del_tranche (grid, 1, 1, 2, 2);
1134    return HOK;
1135 }
1136 // ======================================================== test_pipes
1137 int test_pipes (int nbargs, cpchar tabargs[])
1138 {
1139    Hex::Hex mon_ex;
1140    Hex::Document* doc = mon_ex.addDocument ();
1141
1142    Hex::Vertex* ori1 = doc->addVertex ( 0,0,0);
1143    Hex::Vertex* ori2 = doc->addVertex (-5,0,5);
1144    Hex::Vector* vz   = doc->addVector ( 0,0,1);
1145    Hex::Vector* vx   = doc->addVector ( 1,0,0);
1146
1147 // double h1  =  5, ri1 = 1, re1 = 2;
1148    double h1  = 10, ri1 = 1, re1 = 2;
1149    double h2  = 10, ri2 = 2, re2 = 3;
1150
1151    Hex::Pipe* pipe1  = doc->addPipe (ori1, vz, ri1, re1, h1);
1152    Hex::Pipe* pipe2  = doc->addPipe (ori2, vx, ri2, re2, h2);
1153    Hex::CrossElements* grid = doc->makePipes (pipe1, pipe2);
1154
1155    case_name = "pipe";
1156    docu      = doc;
1157    save_vtk ();
1158
1159    del_tranche (grid, 0, 1, 0);
1160    del_tranche (grid, 0, 1, 5);
1161
1162    del_tranche (grid, 1, 1, 0);
1163    del_tranche (grid, 1, 1, 3);
1164                                   // Partie critique
1165
1166    del_tranche (grid, 1, 1, 1, 2);
1167
1168    del_tranche (grid, 0, 1, 1, 2);
1169    del_tranche (grid, 0, 1, 4, 2);
1170    del_tranche (grid, 0, 1, 3, 2);
1171    del_tranche (grid, 0, 1, 2, 2);
1172
1173    del_tranche (grid, 1, 1, 2, 2);
1174    /* ***************************************************
1175
1176    int nbz [2] =  { 8, 4 };
1177    int ni = 1;
1178    for (int cyl = 0 ; cyl < 2 ; cyl++)
1179        for (int nk = 0 ; nk < nbz[cyl] ; nk++)
1180        for (int nj = 0 ; nj < 4 ; nj++)
1181            {
1182            int jj = nj;
1183            if (cyl==0) jj = (jj+6) MODULO 8 ;
1184            Hex::Hexa* hexa = grid->getHexaIJK (cyl, ni, jj, nk);
1185            if (hexa!=NULL)
1186               {
1187               hexa->remove ();
1188               doc->saveVtk (case_name,  nbr_vtk);
1189               }
1190            }
1191       *************************************************** */
1192    return HOK;
1193 }
1194 // ======================================================== test_lorraine
1195 int test_lorraine(int nbargs, cpchar tabargs[])
1196 {
1197    Hex::Hex mon_ex;
1198    Hex::Document* doc = mon_ex.addDocument ();
1199
1200    Hex::Vertex* ori1 = doc->addVertex ( 0,0,0);
1201    Hex::Vertex* ori2 = doc->addVertex (-5,0,5);
1202    Hex::Vertex* ori3 = doc->addVertex ( 0,0,12);
1203    Hex::Vertex* ori4 = doc->addVertex (-5,0,17);
1204
1205    Hex::Vector* vz   = doc->addVector ( 0,0,1);
1206    Hex::Vector* vx   = doc->addVector ( 1,0,0);
1207
1208    int nl1 = 10;
1209    int nl2 = 10;
1210
1211    double rsmall = 1;
1212    double rmoy   = 2;
1213    double rbig   = 3;
1214
1215    Hex::Cylinder* cyl1  = doc->addCylinder (ori1, vz, rmoy,   nl1);
1216    Hex::Cylinder* cyl2  = doc->addCylinder (ori2, vx, rsmall, nl2);
1217
1218    Hex::Cylinder* cyl3  = doc->addCylinder (ori3, vz, rmoy, nl1);
1219    Hex::Cylinder* cyl4  = doc->addCylinder (ori4, vx, rbig, nl2);
1220
1221    Hex::CrossElements* grid1 = doc->makeCylinders (cyl1, cyl2);
1222    Hex::CrossElements* grid2 = doc->makeCylinders (cyl4, cyl3);
1223
1224 #define Imprimer(x) printf (#x " = ") ; if (x) x->dump() ; else printf ("NULL\n")
1225    const int nx_int = 0;
1226    const int nx_ext = 1;
1227
1228    //           vc2 = grid1->getVertexIJK (Hex::CylBig, 0,0,0);
1229    //           vc3 = grid2->getVertexIJK (Hex::CylSmall, 0,0,0);
1230                                      //    Cyl     i     j     k
1231    Hex::Quad* qb = grid1-> getQuadIJ (Hex::CylBig, nx_ext, Hex::S_E, 4);
1232    Hex::Quad* qh = grid2-> getQuadIJ (Hex::CylSmall, nx_ext, Hex::S_N, 0);
1233
1234    Hex::Vertex* vb0 = qb->getVertex (3);
1235    Hex::Vertex* vb1 = qb->getVertex (2);
1236    Hex::Vertex* vh0 = qh->getVertex (0);
1237    Hex::Vertex* vh1 = qh->getVertex (1);
1238
1239    vb0 = grid1->getVertexIJK (Hex::CylBig, 2, Hex::S_E,  4);  // cible
1240    vb1 = grid1->getVertexIJK (Hex::CylBig, 2, Hex::S_NE, 4);
1241    vh0 = grid2->getVertexIJK (Hex::CylSmall, 2, Hex::S_N,  0);   // depart
1242    vh1 = grid2->getVertexIJK (Hex::CylSmall, 2, Hex::S_NW, 0);
1243
1244    Imprimer (vh0);
1245    Imprimer (vh1);
1246    Imprimer (vb0);
1247    Imprimer (vb1);
1248
1249    // qb->remove ();
1250    // qh->remove ();
1251    Hex::Quads hliste;
1252
1253    hliste.push_back (qh);
1254    for (int ny=1; ny<Hex::S_MAXI; ny++)
1255        {
1256        int ns = (ny + Hex::S_N) MODULO Hex::S_MAXI;
1257        hliste.push_back (grid2->getQuadIJ (Hex::CylSmall, nx_ext, ns, 0));
1258        }
1259
1260    for (int ny=0; ny<4 ;  ny++)
1261        hliste.push_back (grid2->getQuadIJ (Hex::CylSmall, nx_int, ny, 0));
1262
1263    int hauteur = 3;
1264    doc->joinQuads  (hliste, qb, vh0, vb0, vh1, vb1, hauteur);
1265    doc->saveVtk ("lorraine.vtk");
1266
1267    // doc->dump ();
1268    return HOK;
1269 }
1270 // ======================================================== test_disconnect2
1271 // === Disconnect Edge seul
1272 int test_disconnect2 (int nbargs, cpchar tabargs[])
1273 {
1274    const int size_x = 2;
1275    const int size_y = 2;
1276    const int size_z = 1;
1277
1278    Hex::Hex mon_ex;
1279    Hex::Document* doc = mon_ex.addDocument ();
1280
1281    Hex::Vertex*   orig2 = doc->addVertex (0,0,0);
1282    Hex::Vector*   dir   = doc->addVector (1,1,1);
1283    Hex::Elements* grid2 = doc->makeCartesian (orig2, dir, size_x,size_y,size_z);
1284
1285    doc->dump ();
1286
1287    int nvtk = 0;
1288    doc->setLevel (1);
1289    Hex::Matrix  matrice;
1290    Hex::Vector* ecart  = doc->addVector (0.5,0.5,0);
1291    matrice.defTranslation (ecart);
1292
1293    Hex::Hexa* hexa2 = grid2->getHexaIJK (1,1,0);
1294    Hex::Edge* edge  = grid2->getEdgeK   (1,2,0);
1295
1296    hexa2->setScalar  (2);
1297    edge->setScalar   (5);
1298
1299    doc->saveVtk ("test_disco", nvtk);
1300
1301    doc->setLevel (4);
1302
1303    Hex::Elements* disco_edges =  doc->disconnectEdge (hexa2, edge);
1304    HexDisplay (disco_edges->countVertex());
1305    HexDisplay (disco_edges->countEdge());
1306    HexDisplay (disco_edges->countQuad());
1307    HexDisplay (disco_edges->countHexa());
1308
1309    // hexa2->transform (&matrice);
1310  /**********************************
1311    for (int ns=0; ns<disco_edges->countVertex(); ns++)
1312        {
1313        Hex::Vertex* sommet = disco_edges->getVertex(ns);
1314        sommet->setX (sommet->getX()+0.5);
1315        sommet->setY (sommet->getY()+0.5);
1316        }
1317    ********************************* */
1318
1319    doc->saveVtk ("test_disco", nvtk);
1320    doc->save ("test_disco");
1321    doc->dump ();
1322    hexa2->dumpFull ();
1323
1324    doc->setLevel (4);
1325    return HOK;
1326 }
1327 // ======================================================== test_disconnect4
1328 // === Disconnect Edges
1329 int test_disconnect4 (int nbargs, cpchar tabargs[])
1330 {
1331    const int size_x = 2;
1332    const int size_y = 2;
1333    const int size_z = 5;
1334
1335    Hex::Hex mon_ex;
1336    Hex::Document* doc = mon_ex.addDocument ();
1337
1338    Hex::Vertex*   orig2 = doc->addVertex (0,0,0);
1339    Hex::Vector*   dir   = doc->addVector (1,1,1);
1340    Hex::Elements* grid2 = doc->makeCartesian (orig2, dir, size_x,size_y,size_z);
1341
1342    // doc->dump ();
1343
1344    int nvtk = 0;
1345    doc->setLevel (1);
1346
1347    Hex::Hexas t_hexas;
1348    Hex::Edges t_edges;
1349    for (int nk=0 ; nk< size_z; nk++)
1350        {
1351        Hex::Hexa* hexa2 = grid2->getHexaIJK (1,1,nk);
1352        Hex::Edge* edge  = grid2->getEdgeK   (1,2,nk);
1353
1354        hexa2->setScalar  (2);
1355        edge->setScalar   (5);
1356        t_hexas.push_back (hexa2);
1357        t_edges.push_back (edge);
1358
1359        doc->setLevel (4);
1360
1361        }
1362
1363    doc->saveVtk ("test_disco", nvtk);
1364    Hex::Elements* disco_edges =  doc->disconnectEdges (t_hexas, t_edges);
1365    doc->saveVtk ("test_disco", nvtk);
1366    // doc->dump ();
1367    // hexa2->dumpFull ();
1368
1369    doc->setLevel (4);
1370    return HOK;
1371 }
1372 // ======================================================== test_disconnect
1373 // ==== Disconnect Quad
1374 int test_disconnect1 (int nbargs, cpchar tabargs[])
1375 {
1376    const int size_x = 2;
1377    const int size_y = 2;
1378    const int size_z = 1;
1379
1380    Hex::Hex mon_ex;
1381    Hex::Document* doc = mon_ex.addDocument ();
1382
1383    Hex::Vertex*   orig1 = doc->addVertex (0,0,0);
1384    Hex::Vector*   dir   = doc->addVector (1,1,1);
1385    Hex::Elements* grid1 = doc->makeCartesian (orig1, dir, size_x,size_y,size_z);
1386
1387    int nvtk = 0;
1388    doc->setLevel (1);
1389    Hex::Matrix  matrice;
1390    Hex::Vector* ecart  = doc->addVector (0.5,0.5,0);
1391    matrice.defTranslation (ecart);
1392
1393    Hex::Hexa* hexa1 = grid1->getHexaIJK (1,1,0);
1394    Hex::Quad* quad  = grid1->getQuadJK  (1,1,0);
1395
1396    quad->setScalar   (5);
1397
1398    doc->saveVtk ("test_disco", nvtk);
1399    doc->disconnectQuad (hexa1, quad);
1400    // hexa1 ->transform (&matrice);
1401    doc->saveVtk ("test_disco", nvtk);
1402
1403    // doc->dumpPropagation ();
1404    // doc->dump  ();
1405
1406    doc->save  ("disco_all");
1407    return HOK;
1408 }
1409 // ======================================================== test_disconnect3
1410 // ==== disconnectVertex
1411 int test_disconnect3 (int nbargs, cpchar tabargs[])
1412 {
1413    const int size_x = 2;
1414    const int size_y = 2;
1415    const int size_z = 1;
1416
1417    Hex::Hex mon_ex;
1418    Hex::Document* doc = mon_ex.addDocument ();
1419
1420    Hex::Vertex*   orig1 = doc->addVertex (0,0,0);
1421
1422    Hex::Vector*   dir   = doc->addVector (1,1,1);
1423    Hex::Elements* grid1 = doc->makeCartesian (orig1, dir, size_x,size_y,size_z);
1424
1425    int nvtk = 0;
1426    doc->setLevel (1);
1427    Hex::Matrix  matrice;
1428    Hex::Vector* ecart  = doc->addVector (0.5,0.5,0);
1429    matrice.defTranslation (ecart);
1430
1431    Hex::Hexa* hexa1 = grid1->getHexaIJK (1,1,0);
1432    Hex::Vertex* vertex = grid1->getVertexIJK (1,1,1);
1433
1434    vertex->setScalar (5);
1435
1436    doc->saveVtk ("test_disco", nvtk);
1437
1438    doc->disconnectVertex (hexa1, vertex);
1439    // hexa1->transform (&matrice);
1440    doc->saveVtk ("test_disco", nvtk);
1441
1442    // doc->dumpPropagation ();
1443    // doc->dump  ();
1444
1445    doc->save  ("disco_all");
1446    return HOK;
1447 }
1448 // ======================================================== contraction
1449 void contraction (Hex::Hexa* hexa, Hex::Elements* grid)
1450 {
1451    return;
1452    Hex::Real3 cg = { 0, 0, 0 };
1453
1454    for (int nro=0; nro<Hex::HV_MAXI ; nro++)
1455        {
1456        cg [0] += hexa->getVertex(nro)->getX()/Hex::HV_MAXI;
1457        cg [1] += hexa->getVertex(nro)->getY()/Hex::HV_MAXI;
1458        cg [2] += hexa->getVertex(nro)->getZ()/Hex::HV_MAXI;
1459        }
1460
1461    int nbvertex = grid->countVertex();
1462    const double coeff = 0.5;
1463    for (int nro=0; nro<nbvertex ; nro++)
1464        {
1465        Hex::Vertex* pv = grid->getVertex(nro);
1466        Hex::Real3 pold = { pv->getX(), pv->getY(), pv->getZ() };
1467        Hex::Real3 pnew;
1468        for (int dd=0; dd<3 ; dd++)
1469            pnew [dd] = cg[dd] + coeff * (pold[dd]-cg[dd]);
1470
1471        pv->setX (pnew[0]);
1472        pv->setY (pnew[1]);
1473        pv->setZ (pnew[2]);
1474        }
1475 }
1476 // ======================================================== test_disconnect
1477 // ==== Les 3 disconnect
1478 int test_disconnect (int nbargs, cpchar tabargs[])
1479 {
1480    const int size_x = 2;
1481    const int size_y = 2;
1482    const int size_z = 1;
1483
1484    Hex::Hex mon_ex;
1485    Hex::Document* doc = mon_ex.addDocument ();
1486
1487    Hex::Vertex*   orig1 = doc->addVertex (0,0,0);
1488    Hex::Vertex*   orig2 = doc->addVertex (4,0,0);
1489    Hex::Vertex*   orig3 = doc->addVertex (8,0,0);
1490
1491    Hex::Vector*   dir   = doc->addVector (1,1,1);
1492    Hex::Elements* grid1 = doc->makeCartesian (orig1, dir, size_x,size_y,size_z);
1493    Hex::Elements* grid2 = doc->makeCartesian (orig2, dir, size_x,size_y,size_z);
1494    Hex::Elements* grid3 = doc->makeCartesian (orig3, dir, size_x,size_y,size_z);
1495
1496    int nvtk = 0;
1497    doc->setLevel (1);
1498    Hex::Matrix  matrice;
1499    Hex::Vector* ecart  = doc->addVector (0.5,0.5,0);
1500    matrice.defTranslation (ecart);
1501
1502    Hex::Hexa* hexa1 = grid1->getHexaIJK (1,1,0);
1503    Hex::Hexa* hexa2 = grid2->getHexaIJK (1,1,0);
1504    Hex::Hexa* hexa3 = grid3->getHexaIJK (1,1,0);
1505
1506    Hex::Quad* quad  = grid1->getQuadJK  (1,1,0);
1507    Hex::Edge* edge  = grid2->getEdgeK   (1,2,0);
1508    Hex::Vertex* vertex = grid3->getVertexIJK (1,1,1);
1509
1510    quad->setScalar   (5);
1511    edge->setScalar   (5);
1512    vertex->setScalar (5);
1513
1514    doc->saveVtk ("test_disco", nvtk);
1515    doc->disconnectQuad (hexa1, quad);
1516    doc->saveVtk ("test_disco", nvtk);
1517
1518    doc->disconnectEdge (hexa2, edge);
1519    doc->saveVtk ("test_disco", nvtk);
1520
1521    doc->disconnectVertex (hexa3, vertex);
1522    doc->saveVtk ("test_disco", nvtk);
1523
1524    // doc->dumpPropagation ();
1525    // doc->dump  ();
1526
1527    doc->save  ("disco_all");
1528    return HOK;
1529 }
1530 // ======================================================== test_propagation
1531 int test_propagation ()
1532 {
1533    const int size_x = 2;
1534    const int size_y = 1;
1535    const int size_z = 1;
1536
1537    Hex::Hex mon_ex;
1538    Hex::Document* doc = mon_ex.addDocument ();
1539
1540    Hex::Vertex*   orig = doc->addVertex (0,0,0);
1541    Hex::Vector*   dir  = doc->addVector (1,1,1);
1542    //  Hex::Elements* grid =
1543    doc->makeCartesian (orig, dir, size_x,size_y,size_z);
1544
1545    int nb = doc->countPropagation ();
1546    for (int nro=0 ; nro<nb ; nro++)
1547        {
1548        Hex::Propagation*  prop  = doc ->getPropagation (nro);
1549        const Hex::Edges&  table = prop->getEdges ();
1550        printf (" ____________________________________ Prop nro %d\n", nro);
1551        for (int ned=0 ; ned<(int)table.size() ; ned++)
1552            {
1553            bool way = table [ned]->getWay ();
1554
1555            if (way)
1556               {
1557               printf ("     (");
1558               table [ned]->getVertex (0)->printName (", ");
1559               table [ned]->getVertex (1)->printName (")\n");
1560               }
1561           else
1562               {
1563               printf ("     (");
1564               table [ned]->getVertex (1)->printName (", ");
1565               table [ned]->getVertex (0)->printName (")\n");
1566               }
1567            }
1568        }
1569
1570    doc->dump  ();
1571    doc->saveVtk ("test_propagation.vtk");
1572    doc->save ("test_propagation");
1573
1574    return HOK;
1575 }
1576 // ======================================================== test_move
1577 int test_move ()
1578 {
1579    const int size_x = 1;
1580    const int size_y = 1;
1581    const int size_z = 2;
1582
1583    Hex::Hex mon_ex;
1584    Hex::Document* doc = mon_ex.addDocument ();
1585
1586    Hex::Vertex* orig = doc->addVertex (0,0,0);
1587    Hex::Vector* dir  = doc->addVector (1,1,1);
1588    Hex::Elements*  grid = doc->makeCartesian (orig, dir, size_x,size_y,size_z);
1589
1590    Hex::Vector*   enhaut  = doc->addVector (0, 0, 5);
1591    Hex::Vector*   devant  = doc->addVector (5, 0, 0);
1592    // Hex::Vector*   agauche = doc->addVector (0, 5, 0);
1593
1594    Hex::Matrix matrice;
1595    matrice.defTranslation (enhaut);
1596
1597    Hex::Hexa* cube    = grid->getHexa (1);
1598    Hex::Quad* dessous = cube->getQuad (Hex::Q_A);
1599    dessous->dump();
1600
1601    Hex::Elements* grid2 = doc->makeTranslation (grid, devant);
1602    /* Hex::Elements* grid3 = doc->makeTranslation (grid, agauche); */
1603    Hex::Hexa* cube2     = grid2->getHexa (1);
1604
1605    doc ->saveVtk ("move0.vtk");
1606
1607    cube ->disconnectQuad (dessous);
1608    cube ->transform (&matrice);
1609    cube2->transform (&matrice);
1610
1611    doc ->saveVtk ("move1.vtk");
1612    doc ->dump();
1613
1614    return HOK;
1615 }
1616 // ======================================================== test_transfo2
1617 int test_transfo2 (int nbargs, cpchar tabargs[])
1618 {
1619    const int size_x = 1;
1620    const int size_y = 1;
1621    const int size_z = 2;
1622
1623    int    nvtk    = 0;
1624    cpchar fic_vtk = "transfo";
1625
1626    Hex::Hex mon_ex;
1627    Hex::Document* doc = mon_ex.addDocument ();
1628    doc ->setLevel (1);
1629
1630    Hex::Vertex* orig = doc->addVertex (0,0,0);
1631    Hex::Vector* dir  = doc->addVector (1,1,1);
1632    Hex::Elements* grid = doc->makeCartesian (orig, dir, size_x, size_y,
1633                                                         size_z);
1634    if (grid==NULL)
1635       return HERR;
1636
1637    orig->setScalar(2);
1638
1639    doc ->saveVtk (fic_vtk, nvtk);
1640
1641    Hex::Vector*   devant  = doc->addVector (5, 0, 0);
1642
1643    Hex::Elements* grid2 = doc->makeTranslation (grid, devant);
1644    if (grid2==NULL)
1645       return HERR;
1646    doc ->saveVtk (fic_vtk, nvtk);
1647
1648    Hex::Elements* grid3  = doc->makeScale (grid2, orig, 2);
1649    if (grid3==NULL)
1650       return HERR;
1651    doc ->saveVtk (fic_vtk, nvtk);
1652
1653    Hex::Elements* grid4 = doc->makeRotation (grid2, orig, dir, 45);
1654    if (grid4==NULL)
1655       return HERR;
1656    doc ->saveVtk (fic_vtk, nvtk);
1657
1658    Hex::Elements* grid5 = doc->makeSymmetryPoint (grid4, orig);
1659    if (grid5==NULL)
1660       return HERR;
1661
1662    doc ->saveVtk (fic_vtk, nvtk);
1663
1664    Hex::Vector* dir1  = doc->addVector (1,0,0);
1665    Hex::Elements* grid6 = doc->makeSymmetryLine (grid4, orig, dir1);
1666    if (grid6==NULL)
1667       return HERR;
1668
1669    grid4->getHexa(0)->getVertex(0)->setScalar(3);
1670    grid6->getHexa(0)->getVertex(0)->setScalar(3);
1671    doc ->saveVtk (fic_vtk, nvtk);
1672
1673    grid4->getHexa(0)->getVertex(0)->setScalar(0);
1674    grid6->getHexa(0)->getVertex(0)->setScalar(0);
1675
1676    Hex::Elements* grid7 = doc->makeSymmetryLine (grid2, orig, dir1);
1677    if (grid7==NULL)
1678       return HERR;
1679
1680    grid2->getHexa(0)->getVertex(0)->setScalar(4);
1681    grid7->getHexa(0)->getVertex(0)->setScalar(4);
1682    doc ->saveVtk (fic_vtk, nvtk);
1683
1684    grid2->getHexa(0)->getVertex(0)->setScalar(0);
1685    grid7->getHexa(0)->getVertex(0)->setScalar(0);
1686
1687    Hex::Elements* grid8 = doc->makeSymmetryPlane (grid2, orig, dir1);
1688    if (grid8==NULL)
1689       return HERR;
1690
1691    grid2->getHexa(0)->getVertex(0)->setScalar(4);
1692    grid8->getHexa(0)->getVertex(0)->setScalar(4);
1693    doc ->saveVtk (fic_vtk, nvtk);
1694    grid2->getHexa(0)->getVertex(0)->setScalar(0);
1695    grid8->getHexa(0)->getVertex(0)->setScalar(0);
1696
1697    Hex::Elements* grid9 = doc->makeSymmetryPlane (grid3, orig, dir);
1698    if (grid9==NULL)
1699       return HERR;
1700
1701    grid3->getHexa(0)->getVertex(0)->setScalar(4);
1702    grid9->getHexa(0)->getVertex(0)->setScalar(4);
1703    doc ->saveVtk (fic_vtk, nvtk);
1704
1705    grid9->getHexa(0)->removeConnected ();
1706    doc ->saveVtk (fic_vtk, nvtk);
1707
1708    return HOK;
1709 }
1710 // ======================================================== test_transfo
1711 int test_transfo (int nbargs, cpchar tabargs[])
1712 {
1713    int    nvtk    = 0;
1714    cpchar fic_vtk = "transfo";
1715
1716    Hex::Hex mon_ex;
1717    Hex::Document* doc = mon_ex.addDocument ();
1718    doc ->setLevel (1);
1719
1720    Hex::Vertex* orig = doc->addVertex (0,0,0);
1721    Hex::Vector* vx   = doc->addVector (1,0,0);
1722    Hex::Vector* vz   = doc->addVector (0,0,1);
1723    double dr = 1;
1724    double da = 360;
1725    double dl = 1;
1726    int nr = 3;
1727    int na = 8;
1728    int nl = 3;
1729    Hex::Elements* grid = doc->makeCylindrical (orig, vx,vz, dr, da, dl,
1730                                                             nr, na, nl, false);
1731    if (grid==NULL)
1732       return HERR;
1733
1734    doc ->saveVtk (fic_vtk, nvtk);
1735    Hex::Vector*   devant  = doc->addVector (10, 0, 0);
1736
1737    Hex::Elements* grid2 = doc->makeTranslation (grid, devant);
1738    if (grid2==NULL)
1739       return HERR;
1740    doc ->saveVtk (fic_vtk, nvtk);
1741
1742    return HOK;
1743 }
1744 // ======================================================== test_copy_document
1745 int test_copy_document (int nbargs, cpchar tabargs[])
1746 {
1747    Hex::Hex mon_ex;
1748    Hex::Document* doc = mon_ex.loadDocument ("Essai");
1749    doc ->saveVtk ("restore1.vtk");
1750
1751    Hex::Document* clone = doc->copyDocument();
1752    clone->saveVtk ("restore2.vtk");
1753
1754    return HOK;
1755 }
1756 // ======================================================== test_remove
1757 int test_remove ()
1758 {
1759    const int size_x = 2;
1760    const int size_y = 2;
1761    const int size_z = 2;
1762
1763    Hex::Hex mon_ex;
1764    Hex::Document* doc = mon_ex.addDocument ();
1765
1766    Hex::Vertex* orig  = doc->addVertex (0,0,0);
1767    Hex::Vertex* orig1 = doc->addVertex (6,0,0);
1768    Hex::Vector* dir   = doc->addVector (1,1,1);
1769    Hex::Elements* grid  = doc->makeCartesian (orig, dir,  size_x,size_y,size_z);
1770    doc->makeCartesian (orig1, dir, 1,1,1);
1771    doc->saveVtk ("removeConn1.vtk");
1772
1773    Echo ("--------- Avant destruction");
1774    HexDisplay (doc->countVertex ());
1775    HexDisplay (doc->countEdge ());
1776    HexDisplay (doc->countQuad ());
1777    HexDisplay (doc->countHexa ());
1778    HexDisplay (doc->countUsedVertex ());
1779    HexDisplay (doc->countUsedEdge ());
1780    HexDisplay (doc->countUsedQuad ());
1781    HexDisplay (doc->countUsedHexa ());
1782
1783
1784    doc->removeConnectedHexa (grid->getHexaIJK (0,0,0));
1785
1786    Echo ("--------- Apres destruction");
1787    HexDisplay (doc->countVertex ());
1788    HexDisplay (doc->countEdge ());
1789    HexDisplay (doc->countQuad ());
1790    HexDisplay (doc->countHexa ());
1791
1792    HexDisplay (doc->countUsedVertex ());
1793    HexDisplay (doc->countUsedEdge ());
1794    HexDisplay (doc->countUsedQuad ());
1795    HexDisplay (doc->countUsedHexa ());
1796    doc->saveVtk ("removeConn2.vtk");
1797
1798    return HOK;
1799 }
1800 // ================================================== init_vec
1801 void init_vec (Hex::RealVector& tab, double n0=0, double n1=0, double n2=0,
1802                double n3=0, double n4=0, double n5=0, double n6=0,
1803                double n7=0, double n8=0, double n9=0, double n10=0,
1804                double n11=0, double n12=0, double n13=0, double n14=0,
1805                double n15=0, double n16=0)
1806 {
1807    if (n0>0.0) tab.push_back (n0);
1808    if (n1>0.0) tab.push_back (n1);
1809    if (n2>0.0) tab.push_back (n2);
1810    if (n3>0.0) tab.push_back (n3);
1811    if (n4>0.0) tab.push_back (n4);
1812    if (n5>0.0) tab.push_back (n5);
1813    if (n6>0.0) tab.push_back (n6);
1814    if (n7>0.0) tab.push_back (n7);
1815    if (n8>0.0) tab.push_back (n8);
1816    if (n9>0.0) tab.push_back (n9);
1817
1818    if (n10>0.0) tab.push_back (n10);
1819    if (n11>0.0) tab.push_back (n11);
1820    if (n12>0.0) tab.push_back (n12);
1821    if (n13>0.0) tab.push_back (n13);
1822    if (n14>0.0) tab.push_back (n14);
1823    if (n15>0.0) tab.push_back (n15);
1824    if (n16>0.0) tab.push_back (n16);
1825 }
1826 // ================================================== test_cylindricals
1827 int test_cylindricals (int nbargs, cpchar tabargs[])
1828 {
1829    Hex::Hex mon_ex;
1830    Hex::Document* doc = mon_ex.addDocument ();
1831
1832    Hex::Vertex* orig = doc->addVertex (0, 0, 0);
1833    Hex::Vector* vz   = doc->addVector (0, 0, 1);
1834    Hex::Vector* vx   = doc->addVector (1 ,0, 0);
1835
1836    Hex::RealVector tdr, tda, tdl;
1837
1838    /******************
1839    init_vec (tdr, 2, 1, 0.5);
1840    init_vec (tda, 40, 35, 30, 25, 20, 15, 10, 5,
1841                    5, 10, 15, 20, 25, 30, 35, 40);
1842    init_vec (tdl, 1, 2, 3 );
1843
1844    init_vec (tdr, 1, 1, 1, 1);
1845    init_vec (tda, 45,45, 45,45, 45,45, 45,45 );
1846    init_vec (tdl, 1, 1, 1 );
1847
1848     ****************** */
1849
1850
1851    init_vec (tdr, 1, 2, 1, 2);
1852    init_vec (tda, 20, 20, 20 );
1853    init_vec (tdl, 1 );
1854
1855    Hex::Elements* grid=doc->makeCylindricals (orig, vx,vz, tdr,tda,tdl, false);
1856
1857    doc->saveVtk ("cylindricals.vtk");
1858    doc->dump();
1859    grid->clearAssociation();
1860    doc->clearAssociation();
1861    return HOK;
1862 }
1863 // ======================================================== test_hexa
1864 int test_hexa (int nbargs, cpchar tabargs[])
1865 {
1866    goto_workspace ();
1867    int ier = test_cylindricals (nbargs, tabargs);
1868    ier = test_transfo (nbargs, tabargs);
1869    free_workspace ();
1870
1871    return ier;
1872 }