src/HEXABLOCK/test_hexa1.cxx

   1
   2 // C++ : Tests unitaires
   3
   4 // Copyright (C) 2009-2023  CEA, EDF
   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, or (at your option) any later version.
  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 "Hex.hxx"
  26 #include "HexDocument.hxx"
  27 #include "HexElements.hxx"
  28
  29 #include "HexHexa.hxx"
  30 #include "HexQuad.hxx"
  31 #include "HexEdge.hxx"
  32 #include "HexVertex.hxx"
  33
  34 #include "HexPropagation.hxx"
  35 #include "HexLaw.hxx"
  36 #include "HexMatrix.hxx"
  37 #include "HexCramer.hxx"
  38 #include "HexGroup.hxx"
  39
  40
  41 #include <cstdlib>
  42
  43 static int nbr_vtk = 0;
  44 static cpchar case_name = "hexa";
  45 static Hex::Document*   docu = NULL;
  46
  47 // ======================================================== print_propagations
  48 void print_propagations (Hex::Document* doc)
  49 {
  50    int nb = doc->countPropagation ();
  51    HexDisplay (nb);
  52    for (int nro=0 ; nro<nb ; nro++)
  53        {
  54        Hex::Propagation*  prop  = doc ->getPropagation (nro);
  55        const Hex::Edges&  table = prop->getEdges ();
  56        printf (" ____________________________________ Prop nro %d\n", nro);
  57        for (int ned=0 ; ned<(int)table.size() ; ned++)
  58            {
  59            bool way = table [ned]->getWay ();
  60            Hex::Edge*   edge = table [ned];
  61            Hex::Vertex* v0   = edge->getVertex (0);
  62            Hex::Vertex* v1   = edge->getVertex (1);
  63
  64            if (way)
  65               {
  66               printf ("     (");
  67               v0->printName (", ");
  68               v1->printName (")\n");
  69               }
  70           else
  71               {
  72               v1->printName (", ");
  73               v0->printName (")\n");
  74               printf ("     (");
  75               }
  76            }
  77        }
  78 }
  79 // ======================================================== save_vtk
  80 void save_vtk ()
  81 {
  82    if (docu==NULL)
  83       return;
  84
  85    docu->saveVtk (case_name, nbr_vtk);
  86 }
  87 // ======================================================== remove_hexa
  88 void remove_hexa (Hex::Hexa* hexa)
  89 {
  90    if (hexa==NULL)
  91       return;
  92
  93    hexa->remove();
  94
  95    if (docu==NULL)
  96       return;
  97
  98    docu->saveVtk (case_name, nbr_vtk);
  99 }
 100 // ======================================================== test_sphere
 101 int test_sphere (int nbargs, cpchar tabargs[])
 102 {
 103    Hex::Hex mon_ex;
 104    Hex::Document* doc = mon_ex.addDocument ();
 105
 106    int    ncouches = 2;
 107    double k      = 1;
 108    Hex::Elements* sphere = doc->makeSphericalTop (ncouches, k);
 109
 110    if (ncouches>0)
 111       {
 112    for (int nc=0 ; nc <= ncouches ; nc++)
 113        {
 114        Hex::Hexa* cell = sphere->getStrate (nc, Hex::Q_A);
 115        if (nc==0)
 116           cell->dumpFull();
 117        cell->remove ();
 118        // sphere->getStrate (nc, Hex::Q_A)->remove ();
 119        // sphere->getStrate (nc, Hex::Q_B)->remove ();
 120        }
 121        }
 122
 123    sphere->saveVtk ("sphere.vtk");
 124    return HOK;
 125 }
 126 // ======================================================== test_cartesi1
 127 int test_cartesi1 ()
 128 {
 129    const int size_x = 15;
 130    const int size_y = 12;
 131    const int size_z = 8;
 132
 133    Hex::Hex mon_ex;
 134    Hex::Document* doc = mon_ex.addDocument ();
 135
 136    Hex::Elements* grid = doc->makeCartesianTop (size_x,size_y,size_z);
 137
 138    //  print_propagations (doc);
 139
 140    for (int nz=0; nz<size_z ; nz++)
 141        for (int ny=nz+1; ny<size_y-nz-1 ; ny++)
 142            for (int nx=nz+1; nx<size_x-nz-1 ; nx++)
 143                {
 144                Hex::Hexa* cell = grid->getHexaIJK (nx, ny, nz);
 145                cell->remove ();
 146                }
 147
 148    doc->setLevel (1);
 149    print_propagations (doc);
 150    grid->saveVtk ("grid_cart.vtk");
 151
 152    // doc->dump ();
 153    return HOK;
 154 }
 155 // ======================================================== afficher
 156 #define Afficher(elt) afficher (#elt, elt)
 157 int afficher (cpchar nom, Hex::EltBase* elt)
 158 {
 159    if (elt==NULL)
 160       {
 161       printf (" .... %s = 0x0\n", nom);
 162       return HOK;
 163       }
 164
 165    printf (" .... %s = 0x%08lx = %03d\n", nom, (unsigned long) elt, elt->getId());
 166    return HOK;
 167 }
 168 // ======================================================== test_find
 169 int test_find ()
 170 {
 171    const int size_x = 2;
 172    const int size_y = 2;
 173    const int size_z = 2;
 174
 175    Hex::Hex mon_ex;
 176    Hex::Document* doc = mon_ex.addDocument ();
 177
 178    Hex::Elements*  grid = doc->makeCartesianTop (size_x,size_y,size_z);
 179
 180    grid->saveVtk ("mini1.vtk");
 181    doc->dump ();
 182
 183    Hex::Vertex *v00, *v02, *v06, *v08, *v10, *v22, *v26;
 184
 185    Afficher ( v00 = doc->findVertex (0, 0, 0));
 186    Afficher ( v02 = doc->findVertex (1, 1, 0));
 187    Afficher ( v06 = doc->findVertex (1, 1, 1));
 188    Afficher ( v08 = doc->findVertex (2, 1, 0));
 189    Afficher ( v10 = doc->findVertex (2, 1, 1));
 190    Afficher ( v22 = doc->findVertex (2, 1, 2));
 191    Afficher ( v26 = doc->findVertex (2, 2, 2));
 192
 193    printf ("\n");
 194
 195    Afficher (doc->findEdge (v06, v10));
 196    Afficher (doc->findEdge (v10, v06));
 197    printf ("\n");
 198
 199    Afficher (doc->findQuad (v02, v10));
 200    Afficher (doc->findQuad (v06, v08));
 201    Afficher (doc->findQuad (v02, v06));
 202
 203    printf ("\n");
 204    Afficher (doc->findHexa (v00, v06));
 205    Afficher (doc->findHexa (v06, v26));
 206    Afficher (doc->findHexa (v26, v06));
 207
 208    return HOK;
 209 }
 210 // ======================================================== test_joint
 211 int test_joint (int nbargs, cpchar tabargs[])
 212 {
 213    const int dimx = 11;
 214    const int dimy = 11;
 215    const int dimz = 2;
 216
 217    Hex::Hex mon_ex;
 218    Hex::Document* doc = mon_ex.addDocument ();
 219
 220    Hex::Elements* grid1 = doc->makeCartesianTop   (dimx,dimy,dimz);
 221
 222    Hex::Vertex*   orig2 = doc->addVertex  (dimx/2.0,0,8);
 223    Hex::Vector*   vectj = doc->addVector (0,1,0);
 224    Hex::Vector*   vecti = doc->addVector (1,0,0);
 225    Hex::Elements* grid2 = doc->makeCylinderUni (orig2, vecti, vectj,
 226                    1, 180, 1,        dimz,dimy,dimx, true);
 227
 228    int mx = dimx/2;
 229    int my = dimy/2;
 230    Hex::Quad* prems = grid1->getQuadIJ (mx, my, dimz);
 231    Hex::Quad* cible = grid2->getQuadJK (dimz, mx, my);
 232
 233    Hex::Vertex* v1 = prems->getVertex (0);
 234    Hex::Vertex* v3 = prems->getVertex (1);
 235
 236    Hex::Vertex* v2 = cible->getVertex (1);
 237    Hex::Vertex* v4 = cible->getVertex (2);
 238
 239    // v1->setScalar (3);
 240    // v2->setScalar (3);
 241    // v3->setScalar (6);
 242    // v4->setScalar (6);
 243
 244    Hex::Quads liste;
 245
 246    liste.push_back (prems);
 247    for (int nx=0; nx<dimx; nx++)
 248        if (nx!=mx)
 249           liste.push_back (grid1->getQuadIJ (nx, my, dimz));
 250
 251    for (int ny=0; ny<dimy; ny++)
 252        if (ny!=my)
 253           liste.push_back (grid1->getQuadIJ (mx, ny, dimz));
 254
 255    doc->saveVtk ("joint1.vtk");
 256    const int hauteur = 5;
 257    Hex::Elements* joint = doc->joinQuadsUni(liste,cible, v1,v2,v3,v4, hauteur);
 258    // for (int nh=0 ; nh<hauteur ; nh++) joint->getHexa(nh)->setScalar (5);
 259
 260    int nbr_joint_vertex =  joint->countVertex ();
 261    int nbr_surf_vertex  =  nbr_joint_vertex/(hauteur+1);
 262
 263    HexDisplay (nbr_joint_vertex);
 264    HexDisplay (nbr_surf_vertex);
 265
 266    int indice0 = joint->findVertex (v1);
 267    HexDisplay (indice0);
 268
 269    for (int nh=0 ; nh<nbr_surf_vertex ; nh++)
 270        joint->getVertex(nh)->setScalar (5);
 271
 272    for (int nh=0 ; nh<=hauteur ; nh++)
 273        joint->getVertex(nh*nbr_surf_vertex)->setScalar (3);
 274
 275    doc->saveVtk ("joint2.vtk");
 276    return HOK;
 277 }
 278 // ======================================================== test_count
 279 int test_count (int nbargs, cpchar tabargs[])
 280 {
 281    Hex::Hex mon_ex;
 282    Hex::Document* doc = mon_ex.addDocument ();
 283
 284    int    nr = 2;
 285    int    nl = 3;
 286
 287 // Hex::Elements* c1 =
 288    doc->makeCylinderTop (nr, 10, nl);
 289
 290    HexDisplay (doc->countVertex ());
 291    HexDisplay (doc->countUsedVertex ());
 292    doc ->saveVtk ("hexa1.vtk");
 293
 294    return HOK;
 295 }
 296 // ======================================================== test_decoupage
 297 int test_decoupage (int nbargs, cpchar tabargs[])
 298 {
 299    const int size_x = 2;
 300    const int size_y = 1;
 301    const int size_z = 1;
 302
 303    Hex::Hex mon_ex;
 304    Hex::Document* doc = mon_ex.addDocument ();
 305
 306    Hex::Elements* grid  = doc->makeCartesianTop (size_x,size_y,size_z);
 307    Hex::Edge*     arete = grid->getEdgeK (0, 0, 0);
 308
 309    //  doc ->dump ();
 310    int nvtk=0;
 311    doc ->saveVtk ("decoupe", nvtk);
 312 /* Hex::Elements* grid2 = */  doc->cutUni (arete, 1);
 313
 314    ///  doc ->dump ();
 315    doc ->saveVtk ("decoupe", nvtk);
 316
 317    return HOK;
 318 }
 319 // ======================================================== test_gen_xml
 320 int test_gen_xml (int nbargs, cpchar tabargs[])
 321 {
 322    const int size_x = 2;
 323    const int size_y = 2;
 324    const int size_z = 2;
 325
 326    Hex::Hex mon_ex;
 327    Hex::Document* doc = mon_ex.addDocument ();
 328
 329    // Hex::Elements*  grid =
 330    doc->makeCartesianTop (size_x,size_y,size_z);
 331
 332    // Hex::Hexa*   cell    = grid->getHexa (0);
 333    // Hex::Quad*   face    = cell->getQuad (0);
 334    // Hex::Edge*   arete   = cell->getEdge (0);
 335    // Hex::Vertex* noeud   = cell->getVertex (0);
 336
 337    // Hex::Shape* shape1 = new Hex::Shape("riri");
 338    // Hex::Shape* shape2 = new Hex::Shape("fifi");
 339    // Hex::Shape* shape3 = new Hex::Shape("loulou");
 340
 341    // noeud->setAssociation (shape1);
 342    // arete->setAssociation (shape2);
 343    // face ->setAssociation (shape3);
 344
 345    Hex::Law* law1 = doc->addLaw("loi1", 1);
 346    Hex::Law* law2 = doc->addLaw("loi2", 2);
 347    Hex::Law* law3 = doc->addLaw("loi3", 3);
 348
 349    law1->setKind (Hex::Uniform);
 350    law2->setKind (Hex::Arithmetic);
 351    law3->setKind (Hex::Geometric);
 352
 353    Hex::Propagation* prop1 = doc->getPropagation (0);
 354    Hex::Propagation* prop2 = doc->getPropagation (1);
 355    Hex::Propagation* prop3 = doc->getPropagation (2);
 356
 357    prop1->setLaw (law1);
 358    prop2->setLaw (law2);
 359    prop3->setLaw (law3);
 360
 361    prop1->setWay (true);
 362    prop2->setWay (false);
 363    prop3->setWay (true);
 364
 365    doc ->saveVtk ("mini.vtk");
 366    doc ->save ("Essai");
 367
 368    return HOK;
 369 }
 370 // ======================================================== test_string_xml
 371 int test_string_xml (int nbargs, cpchar tabargs[])
 372 {
 373    const int size_x = 2;
 374    const int size_y = 2;
 375    const int size_z = 2;
 376
 377    Hex::Hex mon_ex;
 378    Hex::Document* doc = mon_ex.addDocument ();
 379
 380    // Hex::Elements*  grid =
 381    doc->makeCartesianTop (size_x,size_y,size_z);
 382
 383    Hex::Law* law1 = doc->addLaw("loi1", 1);
 384    Hex::Law* law2 = doc->addLaw("loi2", 2);
 385    Hex::Law* law3 = doc->addLaw("loi3", 3);
 386
 387    law1->setKind (Hex::Uniform);
 388    law2->setKind (Hex::Arithmetic);
 389    law3->setKind (Hex::Geometric);
 390
 391    Hex::Propagation* prop1 = doc->getPropagation (0);
 392    Hex::Propagation* prop2 = doc->getPropagation (1);
 393    Hex::Propagation* prop3 = doc->getPropagation (2);
 394
 395    prop1->setLaw (law1);
 396    prop2->setLaw (law2);
 397    prop3->setLaw (law3);
 398
 399    prop1->setWay (true);
 400    prop2->setWay (false);
 401    prop3->setWay (true);
 402
 403    doc ->saveVtk ("mini.vtk");
 404    doc ->save ("Essai");
 405
 406    cpchar flux = doc ->getXml ();
 407    Hex::Document* docbis = mon_ex.addDocument ();
 408    docbis->setXml  (flux);
 409    docbis->saveVtk ("clone.vtk");
 410
 411    return HOK;
 412 }
 413 // ======================================================== test_relecture
 414 int test_relecture (int nbargs, cpchar tabargs[])
 415 {
 416    Hex::Hex mon_ex;
 417    cpchar nomdoc = "Essai";
 418    if (nbargs>1)
 419       nomdoc = tabargs[1];
 420    Hex::Document* doc = mon_ex.loadDocument (nomdoc);
 421
 422 /*********************
 423    Hex::Vertex* v4 = doc->findVertex (80.0, 0.0,  0.0);
 424    Hex::Vertex* v5 = doc->findVertex (80.0, 0.0, 40.0);
 425    Hex::Edge*   e4 = doc->findEdge   (v4, v5);
 426
 427    HexDump (v4);
 428    HexDump (v5);
 429    HexDump (e4);
 430
 431    e4->setScalar (5);
 432 ***********************/
 433
 434    doc ->dump ();
 435    doc ->saveVtk ("restore.vtk");
 436    doc ->save ("restore");
 437
 438    // doc ->reorderFaces ();
 439    // doc ->dump ();
 440
 441    // Hex::Elements* grid2 = doc->cut (e4, 2);
 442    return HOK;
 443 }
 444 // ======================================================== test_clone
 445 int test_clone ()
 446 {
 447    const int size_x = 2;
 448    const int size_y = 2;
 449    const int size_z = 2;
 450
 451    Hex::Hex mon_ex;
 452    Hex::Document* doc = mon_ex.addDocument ();
 453
 454    Hex::Elements* grid  = doc->makeCartesianTop (size_x,size_y,size_z);
 455    Hex::Vector*   bond  = doc->addVector (0, 0, 7);
 456    Hex::Elements* grid2 = doc->makeTranslation (grid, bond);
 457
 458    doc ->saveVtk ("clonage.vtk");
 459    doc ->dump();
 460
 461    HexDump (grid2->getHexa      (1));
 462    HexDump (grid2->getHexaIJK   (1,1,1));
 463    HexDump (grid2->getVertexIJK (1,1,1));
 464
 465    return HOK;
 466 }
 467 // ======================================================== test_separ
 468 int test_separ ()
 469 {
 470    const int size_x = 2;
 471    const int size_y = 2;
 472    const int size_z = 2;
 473
 474    Hex::Hex mon_ex;
 475    Hex::Document* doc = mon_ex.addDocument ();
 476
 477    doc->makeCartesianTop (size_x,size_y,size_z);
 478
 479    doc ->saveVtk ("separ.vtk");
 480    doc ->dump();
 481
 482    return HOK;
 483 }
 484 // ======================================================== test_shperical
 485 int test_spherical (int nbargs, const char* tabargs[])
 486 {
 487    Hex::Hex mon_ex;
 488    Hex::Document* doc = mon_ex.addDocument ();
 489
 490    int          nbr    = 3;
 491
 492    Hex::Elements* grid = doc->makeSphericalTop (nbr, 1);
 493
 494    int nbhexas = grid->countHexa ();
 495    HexDisplay (nbhexas);
 496    for (int nro=3 ; nro<nbhexas ; nro +=3)
 497        grid->getHexa(nro)->remove();
 498    HexDisplay (doc->countHexa ());
 499    doc->saveVtk ("shperical.vtk");
 500    // doc->dump ();
 501
 502    return HOK;
 503 }
 504 // ================================================== test_grille_cyl
 505 int test_grille_cyl (int nbargs, cpchar tabargs[])
 506 {
 507    Hex::Hex mon_ex;
 508    Hex::Document* doc = mon_ex.addDocument ();
 509
 510    Hex::Vertex* orig1 = doc->addVertex ( 0, 0,0);
 511    Hex::Vertex* orig2 = doc->addVertex (10, 0,0);
 512
 513    Hex::Vertex* orig3 = doc->addVertex ( 0,10,0);
 514    Hex::Vertex* orig4 = doc->addVertex (10,10,0);
 515
 516    Hex::Vertex* orig5 = doc->addVertex ( 0,20,0);
 517    Hex::Vertex* orig6 = doc->addVertex (10,20,0);
 518
 519    Hex::Vector* vz = doc->addVector (0,0,1);
 520    Hex::Vector* vx = doc->addVector (1,0,0);
 521
 522    double rint = 5;
 523    double rext = 1;
 524    double dl = 1;
 525    int    nr = 2;
 526    int    nl = 3;
 527
 528    doc->makeCylinderUni (orig1, vx,vz, rint,rext, 360,dl,nr, 4, nl);
 529    Hex::Elements* c2 =
 530    doc->makeCylinderUni (orig2, vx,vz, rint,rext, 360,dl,nr, 8, nl);
 531    doc->makeCylinderUni (orig3, vx,vz, rint,rext, 270,dl,nr, 8, nl);
 532    doc->makeCylinderUni (orig4, vx,vz, rint,rext, 270,dl,nr, 7, nl);
 533    doc->makeCylinderUni (orig5, vx,vz, rint,rext, 360,dl,nr, 5, nl);
 534    doc->makeCylinderUni (orig6, vx,vz, rint,rext, 360,dl,nr, 6, nl);
 535
 536    int base2 = nr*nl*8;
 537    c2->getHexa(base2 + 0)->setScalar (5);
 538    c2->getHexa(base2 + 1)->setScalar (5);
 539    c2->getHexa(base2 + 2)->setScalar (5);
 540
 541    doc->saveVtk ("cylindres.vtk");
 542    // doc->dump ();
 543
 544    return HOK;
 545 }
 546 // ===================================================== test_cylindrical
 547 int test_cylindrical (int nbargs, cpchar tabargs[])
 548 {
 549    cpchar fic_vtk = "cylindre";
 550
 551    Hex::Hex mon_ex;
 552    Hex::Document* doc = mon_ex.addDocument ();
 553
 554    Hex::Vertex* orig = doc->addVertex (0, 0, 0);
 555    Hex::Vector* vx   = doc->addVector (1, 0, 0);
 556    Hex::Vector* vz   = doc->addVector (0, 0, 1);
 557
 558    double da = 360;
 559    double dl = 1;
 560
 561    int nr = 1;
 562    int na = 8;
 563    int nl = 2;
 564
 565    if (nbargs>1)
 566       {
 567       na = atoi (tabargs[1]);
 568       HexDisplay (na);
 569       if (na <= 2)
 570           na = 8;
 571       }
 572
 573    if (nbargs>2)
 574       {
 575       da = atof (tabargs[2]);
 576       HexDisplay (da);
 577       }
 578
 579    doc->makeCylinderUni (orig,vx, vz, 5.0,1.0,da,dl, nr,na, nl);
 580    doc->saveVtk (fic_vtk, na);
 581    return HOK;
 582 }
 583 // ======================================================== test_joint2
 584 int test_joint2 (int nbargs, cpchar tabargs[])
 585 {
 586    Hex::Hex mon_ex;
 587    docu = mon_ex.addDocument ();
 588    case_name = "pb_joint";
 589
 590    Hex::Vector* vx    = docu->addVector (1, 0, 0);
 591    Hex::Vector* vz    = docu->addVector (0, 0, 1);
 592    Hex::Vertex* hori  = docu->addVertex (0, 0, 0);
 593
 594    double da = 360;
 595    double dl = 1;
 596    int    nr = 1;
 597    int    na = 8;
 598    int    nl = 1;
 599
 600    Hex::Elements *bgrid=NULL, *hgrid=NULL;
 601
 602    hgrid = docu->makeCylinderUni (hori, vx,vz, 5.0,1.0,da,dl, nr,na,nl);
 603    docu->dump ();
 604    save_vtk ();
 605
 606    Hex::Vertex* bori  = docu->addVertex (0, 0, -5);
 607    bgrid = docu->makeCylinderUni (bori, vx,vz, 5.0,1.0,da,dl, nr,na,nl);
 608    save_vtk ();
 609
 610    Hex::Quads qsource, qdest;
 611    printf (" Source = ");
 612    for (int ny=0 ; ny< na ; ny++)
 613        {
 614        Hex::Quad* quad = hgrid->getQuadIJ (0, ny, 0);
 615        PrintName (quad);
 616        qsource.push_back (quad);
 617        }
 618    printf ("\n          ");
 619    for (int ny=0 ; ny<4 ; ny++)
 620        {
 621        Hex::Quad* quad = hgrid->getKerHQuad (ny);
 622        PrintName (quad);
 623        qsource.push_back (quad);
 624        }
 625    printf ("\n");
 626
 627
 628    printf (" Cible  = ");
 629    for (int ny=0 ; ny< na ; ny++)
 630        {
 631        Hex::Quad* quad = bgrid->getQuadIJ (0, ny, 1);
 632        PrintName (quad);
 633        qdest.push_back (quad);
 634        }
 635    printf ("\n          ");
 636    for (int ny=0 ; ny<4 ; ny++)
 637        {
 638        Hex::Quad* quad = bgrid->getKerHQuad (ny+4);
 639        PrintName (quad);
 640        qdest.push_back (quad);
 641        }
 642    printf ("\n");
 643    docu ->setLevel (1);
 644    Hex::Quad*   cible = bgrid->getQuadIJ    (0, 1, 1);
 645    Hex::Vertex* vc1   = bgrid->getVertexIJK (0, 1, 1);
 646    Hex::Vertex* vc2   = bgrid->getVertexIJK (1, 1, 1);
 647
 648    Hex::Vertex* vs1 = hgrid->getVertexIJK (0, 0, 0);
 649    Hex::Vertex* vs2 = hgrid->getVertexIJK (1, 0, 0);
 650
 651    docu->joinQuadsUni (qsource, cible, vs1, vc1, vs2, vc2, 1);
 652    save_vtk ();
 653
 654    return HOK;
 655 }
 656 // ======================================================== test_disconnect2
 657 // === Disconnect Edge seul
 658 int test_disconnect2 (int nbargs, cpchar tabargs[])
 659 {
 660    const int size_x = 2;
 661    const int size_y = 2;
 662    const int size_z = 1;
 663
 664    Hex::Hex mon_ex;
 665    Hex::Document* doc   = mon_ex.addDocument ();
 666    Hex::Elements* grid2 = doc->makeCartesianTop (size_x,size_y,size_z);
 667
 668    int nvtk = 0;
 669    doc->setLevel (1);
 670
 671    Hex::Hexa* hexa2 = grid2->getHexaIJK (1,1,0);
 672    Hex::Edge* edge  = grid2->getEdgeK   (1,2,0);
 673
 674    hexa2->setColor  (2);
 675    edge->setColor   (5);
 676
 677    doc->saveVtk ("test_disco", nvtk);
 678
 679    doc->setLevel (4);
 680
 681    Hex::Elements* disco_edges =  doc->disconnectEdge (hexa2, edge);
 682    HexDisplay (disco_edges->countVertex());
 683    HexDisplay (disco_edges->countEdge());
 684    HexDisplay (disco_edges->countQuad());
 685    HexDisplay (disco_edges->countHexa());
 686
 687    doc->saveVtk ("test_disco", nvtk);
 688    doc->save ("test_disco");
 689    // doc->dump ();
 690    // hexa2->dumpFull ();
 691
 692    // doc->setLevel (4);
 693    return HOK;
 694 }
 695 // ======================================================== test_disconnect4
 696 // === Disconnect Edges
 697 int test_disconnect4 (int nbargs, cpchar tabargs[])
 698 {
 699    const int size_x = 1;
 700    const int size_y = 2;
 701    const int size_z = 2;
 702
 703    Hex::Hex mon_ex;
 704    Hex::Document* doc = mon_ex.addDocument ();
 705
 706    Hex::Elements* grid2 = doc->makeCartesianTop  (size_x,size_y,size_z);
 707
 708    int nvtk = 0;
 709    doc->setLevel (1);
 710
 711    Hex::Hexas t_hexas;
 712    Hex::Edges t_edges;
 713    for (int nk=0 ; nk< size_z; nk++)
 714        {
 715        Hex::Hexa* hexa2 = grid2->getHexaIJK (0,0,nk);
 716        Hex::Edge* edge  = grid2->getEdgeK   (0,1,nk);
 717
 718        hexa2->setScalar  (2);
 719        edge->setScalar   (5);
 720        t_hexas.push_back (hexa2);
 721        t_edges.push_back (edge);
 722
 723        doc->setLevel (4);
 724        }
 725
 726    doc->dump ();
 727    doc->saveVtk ("test_disco", nvtk);
 728    /* Hex::Elements* disco_edges = */  doc->disconnectEdges (t_hexas, t_edges);
 729    doc->saveVtk ("test_disco", nvtk);
 730    doc->dump ();
 731    // hexa2->dumpFull ();
 732
 733    doc->setLevel (4);
 734    return HOK;
 735 }
 736 // ======================================================== test_disconnect1
 737 // ==== Disconnect Quad
 738 int test_disconnect1 (int nbargs, cpchar tabargs[])
 739 {
 740    const int size_x = 2;
 741    const int size_y = 2;
 742    const int size_z = 1;
 743
 744    Hex::Hex mon_ex;
 745    Hex::Document* doc = mon_ex.addDocument ();
 746
 747    Hex::Elements* grid1 = doc->makeCartesianTop (size_x,size_y,size_z);
 748
 749    int nvtk = 0;
 750    doc->setLevel (1);
 751    Hex::Matrix  matrice;
 752    Hex::Vector* ecart  = doc->addVector (0.5,0.5,0);
 753    matrice.defTranslation (ecart);
 754
 755    Hex::Hexa* hexa1 = grid1->getHexaIJK (1,1,0);
 756    Hex::Quad* quad  = grid1->getQuadJK  (1,1,0);
 757
 758    quad->setScalar   (5);
 759
 760    doc->saveVtk ("test_disco", nvtk);
 761    doc->disconnectQuad (hexa1, quad);
 762    // hexa1 ->transform (&matrice);
 763    doc->saveVtk ("test_disco", nvtk);
 764
 765    // doc->dumpPropagation ();
 766    // doc->dump  ();
 767
 768    doc->save  ("disco_all");
 769    return HOK;
 770 }
 771 // ======================================================== test_disconnect3
 772 // ==== disconnectVertex
 773 int test_disconnect3 (int nbargs, cpchar tabargs[])
 774 {
 775    const int size_x = 2;
 776    const int size_y = 2;
 777    const int size_z = 1;
 778
 779    Hex::Hex mon_ex;
 780    Hex::Document* doc = mon_ex.addDocument ();
 781
 782    Hex::Elements* grid1 = doc->makeCartesianTop (size_x,size_y,size_z);
 783
 784    int nvtk = 0;
 785    doc->setLevel (1);
 786    Hex::Matrix  matrice;
 787    Hex::Vector* ecart  = doc->addVector (0.5,0.5,0);
 788    matrice.defTranslation (ecart);
 789
 790    Hex::Hexa* hexa1 = grid1->getHexaIJK (1,1,0);
 791    Hex::Vertex* vertex = grid1->getVertexIJK (1,1,1);
 792
 793    vertex->setScalar (5);
 794
 795    doc->saveVtk ("test_disco", nvtk);
 796
 797    doc->disconnectVertex (hexa1, vertex);
 798    // hexa1->transform (&matrice);
 799    doc->saveVtk ("test_disco", nvtk);
 800
 801    // doc->dumpPropagation ();
 802    // doc->dump  ();
 803
 804    doc->save  ("disco_all");
 805    return HOK;
 806 }
 807 // ======================================================== contraction
 808 void contraction (Hex::Hexa* hexa, Hex::Elements* grid)
 809 {
 810    return;
 811    Hex::Real3 cg = { 0, 0, 0 };
 812
 813    for (int nro=0; nro<Hex::HV_MAXI ; nro++)
 814        {
 815        cg [0] += hexa->getVertex(nro)->getX()/Hex::HV_MAXI;
 816        cg [1] += hexa->getVertex(nro)->getY()/Hex::HV_MAXI;
 817        cg [2] += hexa->getVertex(nro)->getZ()/Hex::HV_MAXI;
 818        }
 819
 820    int nbvertex = grid->countVertex();
 821    const double coeff = 0.5;
 822    for (int nro=0; nro<nbvertex ; nro++)
 823        {
 824        Hex::Vertex* pv = grid->getVertex(nro);
 825        Hex::Real3 pold = { pv->getX(), pv->getY(), pv->getZ() };
 826        Hex::Real3 pnew;
 827        for (int dd=0; dd<3 ; dd++)
 828            pnew [dd] = cg[dd] + coeff * (pold[dd]-cg[dd]);
 829
 830        pv->setX (pnew[0]);
 831        pv->setY (pnew[1]);
 832        pv->setZ (pnew[2]);
 833        }
 834 }
 835 // ======================================================== test_disconnect
 836 // ==== Les 3 disconnect
 837 int test_disconnect (int nbargs, cpchar tabargs[])
 838 {
 839    const int size_x = 2;
 840    const int size_y = 2;
 841    const int size_z = 1;
 842
 843    Hex::Hex mon_ex;
 844    Hex::Document* doc = mon_ex.addDocument ();
 845
 846    Hex::Vertex*   orig1 = doc->addVertex (0,0,0);
 847    Hex::Vertex*   orig2 = doc->addVertex (4,0,0);
 848    Hex::Vertex*   orig3 = doc->addVertex (8,0,0);
 849
 850    Hex::Vector*   vx   = doc->addVector (1,0,0);
 851    Hex::Vector*   vy   = doc->addVector (0,1,0);
 852    Hex::Vector*   vz   = doc->addVector (0,0,1);
 853    double dx=1, dy=1, dz=1;
 854
 855    Hex::Elements* grid1 = doc->makeCartesianUni (orig1, vx,vy,vz, dx,dy,dz,
 856                                                  size_x,size_y,size_z);
 857    Hex::Elements* grid2 = doc->makeCartesianUni (orig2, vx,vy,vz, dx,dy,dz,
 858                                                  size_x,size_y,size_z);
 859    Hex::Elements* grid3 = doc->makeCartesianUni (orig3, vx,vy,vz, dx,dy,dz,
 860                                                  size_x,size_y,size_z);
 861
 862    int nvtk = 0;
 863    doc->setLevel (1);
 864    Hex::Matrix  matrice;
 865    Hex::Vector* ecart  = doc->addVector (0.5,0.5,0);
 866    matrice.defTranslation (ecart);
 867
 868    Hex::Hexa* hexa1 = grid1->getHexaIJK (1,1,0);
 869    Hex::Hexa* hexa2 = grid2->getHexaIJK (1,1,0);
 870    Hex::Hexa* hexa3 = grid3->getHexaIJK (1,1,0);
 871
 872    Hex::Quad* quad  = grid1->getQuadJK  (1,1,0);
 873    Hex::Edge* edge  = grid2->getEdgeK   (1,2,0);
 874    Hex::Vertex* vertex = grid3->getVertexIJK (1,1,1);
 875
 876    quad->setScalar   (5);
 877    edge->setScalar   (5);
 878    vertex->setScalar (5);
 879
 880    doc->saveVtk ("test_disco", nvtk);
 881    doc->disconnectQuad (hexa1, quad);
 882    doc->saveVtk ("test_disco", nvtk);
 883
 884    doc->disconnectEdge (hexa2, edge);
 885    doc->saveVtk ("test_disco", nvtk);
 886
 887    doc->disconnectVertex (hexa3, vertex);
 888    doc->saveVtk ("test_disco", nvtk);
 889
 890    // doc->dumpPropagation ();
 891    // doc->dump  ();
 892
 893    doc->save  ("disco_all");
 894    return HOK;
 895 }
 896 // ======================================================== test_propagation
 897 int test_propagation ()
 898 {
 899    const int size_x = 2;
 900    const int size_y = 1;
 901    const int size_z = 1;
 902
 903    Hex::Hex mon_ex;
 904    Hex::Document* doc = mon_ex.addDocument ();
 905
 906    doc->makeCartesianTop (size_x,size_y,size_z);
 907
 908    int nb = doc->countPropagation ();
 909    for (int nro=0 ; nro<nb ; nro++)
 910        {
 911        Hex::Propagation*  prop  = doc ->getPropagation (nro);
 912        const Hex::Edges&  table = prop->getEdges ();
 913        printf (" ____________________________________ Prop nro %d\n", nro);
 914        for (int ned=0 ; ned<(int)table.size() ; ned++)
 915            {
 916            bool way = table [ned]->getWay ();
 917
 918            if (way)
 919               {
 920               printf ("     (");
 921               table [ned]->getVertex (0)->printName (", ");
 922               table [ned]->getVertex (1)->printName (")\n");
 923               }
 924           else
 925               {
 926               printf ("     (");
 927               table [ned]->getVertex (1)->printName (", ");
 928               table [ned]->getVertex (0)->printName (")\n");
 929               }
 930            }
 931        }
 932
 933    doc->dump  ();
 934    doc->saveVtk ("test_propagation.vtk");
 935    doc->save ("test_propagation");
 936
 937    return HOK;
 938 }
 939 // ======================================================== test_move
 940 int test_move ()
 941 {
 942    const int size_x = 1;
 943    const int size_y = 1;
 944    const int size_z = 2;
 945
 946    Hex::Hex mon_ex;
 947    Hex::Document* doc = mon_ex.addDocument ();
 948
 949    Hex::Elements*  grid = doc->makeCartesianTop (size_x,size_y,size_z);
 950
 951    Hex::Vector*   enhaut  = doc->addVector (0, 0, 5);
 952    Hex::Vector*   devant  = doc->addVector (5, 0, 0);
 953    // Hex::Vector*   agauche = doc->addVector (0, 5, 0);
 954
 955    Hex::Matrix matrice;
 956    matrice.defTranslation (enhaut);
 957
 958    Hex::Hexa* cube    = grid->getHexa (1);
 959    Hex::Quad* dessous = cube->getQuad (Hex::Q_A);
 960    dessous->dump();
 961
 962    Hex::Elements* grid2 = doc->makeTranslation (grid, devant);
 963    /* Hex::Elements* grid3 = doc->makeTranslation (grid, agauche); */
 964    Hex::Hexa* cube2     = grid2->getHexa (1);
 965
 966    doc ->saveVtk ("move0.vtk");
 967
 968    cube ->disconnectQuad (dessous);
 969    cube ->transform (&matrice);
 970    cube2->transform (&matrice);
 971
 972    doc ->saveVtk ("move1.vtk");
 973    doc ->dump();
 974
 975    return HOK;
 976 }
 977 // ======================================================== test_transfo2
 978 int test_transfo2 (int nbargs, cpchar tabargs[])
 979 {
 980    const int size_x = 1;
 981    const int size_y = 1;
 982    const int size_z = 2;
 983
 984    int    nvtk    = 0;
 985    cpchar fic_vtk = "transfo";
 986
 987    Hex::Hex mon_ex;
 988    Hex::Document* doc = mon_ex.addDocument ();
 989    doc ->setLevel (1);
 990
 991    Hex::Vertex* orig = doc->addVertex (0,0,0);
 992    Hex::Vector* dir  = doc->addVector (1,1,1);
 993    Hex::Elements* grid = doc->makeCartesianTop (size_x, size_y, size_z);
 994    if (grid==NULL)
 995       return HERR;
 996
 997    orig->setScalar(2);
 998
 999    doc ->saveVtk (fic_vtk, nvtk);
1000
1001    Hex::Vector*   devant  = doc->addVector (5, 0, 0);
1002
1003    Hex::Elements* grid2 = doc->makeTranslation (grid, devant);
1004    if (grid2==NULL)
1005       return HERR;
1006    doc ->saveVtk (fic_vtk, nvtk);
1007
1008    Hex::Elements* grid3  = doc->makeScale (grid2, orig, 2);
1009    if (grid3==NULL)
1010       return HERR;
1011    doc ->saveVtk (fic_vtk, nvtk);
1012
1013    Hex::Elements* grid4 = doc->makeRotation (grid2, orig, dir, 45);
1014    if (grid4==NULL)
1015       return HERR;
1016    doc ->saveVtk (fic_vtk, nvtk);
1017
1018    Hex::Elements* grid5 = doc->makeSymmetryPoint (grid4, orig);
1019    if (grid5==NULL)
1020       return HERR;
1021
1022    doc ->saveVtk (fic_vtk, nvtk);
1023
1024    Hex::Vector* dir1  = doc->addVector (1,0,0);
1025    Hex::Elements* grid6 = doc->makeSymmetryLine (grid4, orig, dir1);
1026    if (grid6==NULL)
1027       return HERR;
1028
1029    grid4->getHexa(0)->getVertex(0)->setScalar(3);
1030    grid6->getHexa(0)->getVertex(0)->setScalar(3);
1031    doc ->saveVtk (fic_vtk, nvtk);
1032
1033    grid4->getHexa(0)->getVertex(0)->setScalar(0);
1034    grid6->getHexa(0)->getVertex(0)->setScalar(0);
1035
1036    Hex::Elements* grid7 = doc->makeSymmetryLine (grid2, orig, dir1);
1037    if (grid7==NULL)
1038       return HERR;
1039
1040    grid2->getHexa(0)->getVertex(0)->setScalar(4);
1041    grid7->getHexa(0)->getVertex(0)->setScalar(4);
1042    doc ->saveVtk (fic_vtk, nvtk);
1043
1044    grid2->getHexa(0)->getVertex(0)->setScalar(0);
1045    grid7->getHexa(0)->getVertex(0)->setScalar(0);
1046
1047    Hex::Elements* grid8 = doc->makeSymmetryPlane (grid2, orig, dir1);
1048    if (grid8==NULL)
1049       return HERR;
1050
1051    grid2->getHexa(0)->getVertex(0)->setScalar(4);
1052    grid8->getHexa(0)->getVertex(0)->setScalar(4);
1053    doc ->saveVtk (fic_vtk, nvtk);
1054    grid2->getHexa(0)->getVertex(0)->setScalar(0);
1055    grid8->getHexa(0)->getVertex(0)->setScalar(0);
1056
1057    Hex::Elements* grid9 = doc->makeSymmetryPlane (grid3, orig, dir);
1058    if (grid9==NULL)
1059       return HERR;
1060
1061    grid3->getHexa(0)->getVertex(0)->setScalar(4);
1062    grid9->getHexa(0)->getVertex(0)->setScalar(4);
1063    doc ->saveVtk (fic_vtk, nvtk);
1064
1065    grid9->getHexa(0)->removeConnected ();
1066    doc ->saveVtk (fic_vtk, nvtk);
1067
1068    return HOK;
1069 }
1070 // ======================================================== test_copy_document
1071 int test_copy_document (int nbargs, cpchar tabargs[])
1072 {
1073    Hex::Hex mon_ex;
1074    Hex::Document* doc = mon_ex.loadDocument ("Essai");
1075    doc ->saveVtk ("restore1.vtk");
1076
1077    Hex::Document* clone = doc->copyDocument();
1078    clone->saveVtk ("restore2.vtk");
1079
1080    return HOK;
1081 }
1082 // ======================================================== test_remove
1083 int test_remove ()
1084 {
1085    const int size_x = 2;
1086    const int size_y = 2;
1087    const int size_z = 2;
1088
1089    Hex::Hex mon_ex;
1090    Hex::Document* doc = mon_ex.addDocument ();
1091
1092    Hex::Elements* grid  = doc->makeCartesianTop (size_x,size_y,size_z);
1093    doc->saveVtk ("removeConn1.vtk");
1094
1095    Echo ("--------- Avant destruction");
1096    HexDisplay (doc->countVertex ());
1097    HexDisplay (doc->countEdge ());
1098    HexDisplay (doc->countQuad ());
1099    HexDisplay (doc->countHexa ());
1100    HexDisplay (doc->countUsedVertex ());
1101    HexDisplay (doc->countUsedEdge ());
1102    HexDisplay (doc->countUsedQuad ());
1103    HexDisplay (doc->countUsedHexa ());
1104
1105
1106    doc->removeConnectedHexa (grid->getHexaIJK (0,0,0));
1107
1108    Echo ("--------- Apres destruction");
1109    HexDisplay (doc->countVertex ());
1110    HexDisplay (doc->countEdge ());
1111    HexDisplay (doc->countQuad ());
1112    HexDisplay (doc->countHexa ());
1113
1114    HexDisplay (doc->countUsedVertex ());
1115    HexDisplay (doc->countUsedEdge ());
1116    HexDisplay (doc->countUsedQuad ());
1117    HexDisplay (doc->countUsedHexa ());
1118    doc->saveVtk ("removeConn2.vtk");
1119
1120    return HOK;
1121 }
1122 // ================================================== init_vec
1123 void init_vec (Hex::RealVector& tab, double n0=0, double n1=0, double n2=0,
1124                double n3=0, double n4=0, double n5=0, double n6=0,
1125                double n7=0, double n8=0, double n9=0, double n10=0,
1126                double n11=0, double n12=0, double n13=0, double n14=0,
1127                double n15=0, double n16=0)
1128 {
1129    if (n0>0.0) tab.push_back (n0);
1130    if (n1>0.0) tab.push_back (n1);
1131    if (n2>0.0) tab.push_back (n2);
1132    if (n3>0.0) tab.push_back (n3);
1133    if (n4>0.0) tab.push_back (n4);
1134    if (n5>0.0) tab.push_back (n5);
1135    if (n6>0.0) tab.push_back (n6);
1136    if (n7>0.0) tab.push_back (n7);
1137    if (n8>0.0) tab.push_back (n8);
1138    if (n9>0.0) tab.push_back (n9);
1139
1140    if (n10>0.0) tab.push_back (n10);
1141    if (n11>0.0) tab.push_back (n11);
1142    if (n12>0.0) tab.push_back (n12);
1143    if (n13>0.0) tab.push_back (n13);
1144    if (n14>0.0) tab.push_back (n14);
1145    if (n15>0.0) tab.push_back (n15);
1146    if (n16>0.0) tab.push_back (n16);
1147 }
1148 // ======================================================== test_edge
1149 int test_edge (int nbargs, cpchar tabargs[])
1150 {
1151    Hex::Hex mon_ex;
1152    Hex::Document* doc = mon_ex.addDocument ();
1153
1154    Hex::Vertex* orig = doc->addVertex (0, 0, 0);
1155    Hex::Vector* vx   = doc->addVector (1 ,0, 0);
1156    doc->addEdgeVector   (orig, vx);
1157
1158    HexDisplay (doc->countVertex());
1159    HexDisplay (doc->countEdge());
1160    doc->dump ();
1161    return HOK;
1162 }
1163 // ======================================================== test_hexa
1164 int test_hexa (int nbargs, cpchar tabargs[])
1165 {
1166    int ier = 0;
1167    Hex::Hex mon_ex;
1168    Hex::Document* doc1 = mon_ex.loadDocument ("bielle_triang");
1169    Hex::Document* doc2 = mon_ex.loadDocument ("bielle_triang");
1170    PutData (doc1->getName ());
1171    PutData (doc2->getName ());
1172    return ier;
1173
1174
1175    goto_workspace ();
1176    free_workspace ();
1177
1178    return ier;
1179 }