doc/pyplots/bride.py

   1 # -*- coding: latin-1 -*-
   2 # Copyright (C) 2009-2023  CEA, EDF
   3 #
   4 # This library is free software; you can redistribute it and/or
   5 # modify it under the terms of the GNU Lesser General Public
   6 # License as published by the Free Software Foundation; either
   7 # version 2.1 of the License, or (at your option) any later version.
   8 #
   9 # This library is distributed in the hope that it will be useful,
  10 # but WITHOUT ANY WARRANTY; without even the implied warranty of
  11 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  12 # Lesser General Public License for more details.
  13 #
  14 # You should have received a copy of the GNU Lesser General Public
  15 # License along with this library; if not, write to the Free Software
  16 # Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  17 #
  18 # See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  19 #
  20
  21 import os
  22 import GEOM
  23 import geompy
  24 import smesh
  25 import hexablock
  26 import math
  27 import SALOMEDS
  28
  29 k1 = 1
  30
  31 OPT_QUAD_IK = 1
  32 OPT_FIRST = 2
  33
  34 count = 1
  35 def save_schema(doc):
  36     """
  37     sauvegarde vtk du modele de bloc
  38     """
  39     global count
  40     file_name = os.path.join(os.environ['TMP'], 'bride' + str(count) + '.vtk')
  41     doc.saveVtk(file_name)
  42     count += 1
  43     pass
  44
  45 def merge_quads(doc, quart, demi, ni1, nj1, ni2, nj2, option=0) :
  46     """
  47     fusion des quadrangles entre les 2 grilles cylindriques :
  48     """
  49
  50     prems = False
  51     if option == OPT_FIRST:
  52         prems = True
  53
  54     quad_ik = False
  55     if option == OPT_QUAD_IK:
  56         quad_ik = True
  57
  58     orig = None
  59     if quad_ik:
  60         orig = grille_cyl_quart.getQuadIK(ni1, nj1, k1)
  61     else:
  62         orig = grille_cyl_quart.getQuadJK(ni1, nj1, k1)
  63
  64     dest = grille_cyl_demi.getQuadJK(ni2, nj2, k1)
  65
  66
  67 # JPL le 10/05/2011 :
  68 # closeQuads() n'est pas accessible en python
  69 # a priori fonctionne avec mergeQuads()
  70 ##     if prems:
  71     if True:
  72         iq1 = 0
  73         if quad_ik:
  74             iq1 = 1
  75         iq3 = 1 - iq1
  76         v1 = dest.getVertex(iq1)
  77         v3 = dest.getVertex(iq3)
  78
  79         v2 = orig.getVertex(0)
  80         v4 = orig.getVertex(1)
  81
  82         doc.mergeQuads(dest, orig, v1, v2, v3, v4)
  83         pass
  84     else:
  85 ##         doc.closeQuads(dest, orig)
  86         print("closeQuads() : not yet implemented")
  87         pass
  88
  89     return None
  90
  91 BREP_PATH = os.path.expandvars("$HEXABLOCK_ROOT_DIR/bin/salome/bride.brep")
  92
  93 #=============================
  94 # CREATION DOCUMENT
  95 #=============================
  96
  97 doc = hexablock.addDocument()
  98
  99 #=============================
 100 # PARAMETRES
 101 #=============================
 102
 103 height = 1.0
 104
 105 # cylinder grid 1 :
 106 nr1 = 8
 107 na1 = 4
 108 nl1 = 5
 109 dr1 = 1.0
 110 da1 = 45.0  # angle
 111 dl1 = height
 112
 113 # cylinder grid 2 :
 114 nr2 = 3
 115 na2 = 8
 116 nl2 = nl1
 117 dr2 = 0.5
 118 da2 = 180.0  # angle
 119 dl2 = dl1
 120
 121 #=============================
 122 # Creation du modele de blocs
 123 #=============================
 124
 125 # JPL (le 09/05/2011)
 126 # repris de test_bride_abu.cxx (version la plus a jour dans ~/IHMHEXA/Alain/models):
 127 # et de BRIDE.py (Karima) pour les "vraies" coordonn�es :
 128
 129 #=================================================
 130 # Creation des centres des grilles cylindriques
 131 #=================================================
 132
 133 center1 = doc.addVertex(0, 0, height)
 134 center2 = doc.addVertex(6, 0, height)
 135
 136 dx = doc.addVector(height, 0, 0)
 137 dz = doc.addVector(0, 0, height)
 138
 139 # Creation des grilles cylindriques initiales
 140 #============================================
 141
 142 # 1. 1 ere grille (quart) :
 143 #==========================
 144
 145 # JPL (le 10/05/2011) : vu la geometrie, on ne remplit pas le centre :
 146 ## grille_cyl_quart = doc.makeCylindrical(orig1, dx, dz, dr_q, 45.0, dl,
 147 ##                                        nr_q, na_q, dim_z, True)
 148 grille_cyl_quart = doc.makeCylindrical(center1, dx, dz, dr1, da1, dl1,
 149                                        nr1, na1, nl1, False)
 150
 151
 152 # temporaire : sauvegarde du modele de blocs :
 153 save_schema(doc)
 154 # fin temporaire
 155
 156 # Elagage  :
 157 for nk in range(2, nl1):
 158     for nj in range(na1):
 159         ideb = 2
 160         if nk == nl1 - 1:
 161             ideb = 1
 162         for ni in range(ideb, nr1):
 163             doc.removeHexa(grille_cyl_quart.getHexaIJK(ni, nj, nk))
 164             pass
 165         pass
 166     pass
 167
 168 # temporaire : sauvegarde du modele de blocs :
 169 save_schema(doc)
 170 # fin temporaire
 171
 172 # Semelle :
 173 k0 = 0
 174 for nj in range(na1):
 175     for ni in range(2, nr1):
 176         doc.removeHexa(grille_cyl_quart.getHexaIJK(ni, nj, k0))
 177         pass
 178     pass
 179
 180 # temporaire : sauvegarde du modele de blocs :
 181 save_schema(doc)
 182 # fin temporaire
 183
 184 # @todo JPL : peut-on fusionner les edges du haut des hexaedres du haut du deuxieme
 185 # rang ? (cf. GEOM). Si oui revoir aussi l'association
 186
 187 # 2. 2�me grille (demi) :
 188 #========================
 189 grille_cyl_demi = doc.makeCylindrical(center2, dx, dz, dr2, da2, dl2,
 190                                       nr2, na2, nl2, True)
 191
 192 # temporaire : sauvegarde du modele de blocs :
 193 save_schema(doc)
 194 # fin temporaire
 195
 196 ni0 = [0, nr2, 2, 1, 0]  # en fonction de z (ie : nk)
 197 for nk in range(0, nl2):
 198     for nj in range(na2):  # elagage suivant toute la demi-circonference
 199         for ni in range(ni0[nk], nr2):
 200             doc.removeHexa(grille_cyl_demi.getHexaIJK(ni, nj, nk))
 201             pass
 202         pass
 203     pass
 204
 205
 206 # temporaire : sauvegarde du modele de blocs :
 207 save_schema(doc)
 208 # fin temporaire
 209
 210 # 3. creusement des fondations de demi dans quart :
 211 #==================================================
 212 for nj in range(2):
 213     for ni in range(3, nr1 - 1):
 214         doc.removeHexa(grille_cyl_quart.getHexaIJK(ni, nj, k1))
 215         pass
 216     pass
 217
 218 # temporaire : sauvegarde du modele de blocs :
 219 save_schema(doc)
 220 # fin temporaire
 221
 222 # 4. Fusion des bords :
 223 #======================
 224 merge_quads(doc, grille_cyl_quart, grille_cyl_demi, 7, 0, nr2, 0, OPT_FIRST)
 225 merge_quads(doc, grille_cyl_quart, grille_cyl_demi, 7, 1, nr2, 1)
 226 for ni1 in range(2, 6):
 227     merge_quads(doc, grille_cyl_quart, grille_cyl_demi, 8 - ni1, 2, nr2, ni1, OPT_QUAD_IK)
 228     pass
 229 merge_quads(doc, grille_cyl_quart, grille_cyl_demi, 3, 1, nr2, 6)
 230 merge_quads(doc, grille_cyl_quart, grille_cyl_demi, 3, 0, nr2, 7)
 231
 232 # temporaire : sauvegarde du modele de blocs :
 233 save_schema(doc)
 234 # fin temporaire
 235
 236 ###########
 237 # Geometry
 238 ###########
 239
 240 bride_geom = geompy.ImportFile(BREP_PATH, "BREP")
 241
 242 geompy.addToStudy(bride_geom, "bride_geom")
 243
 244
 245 # parametres de la geometrie :
 246
 247 r1 = 12.0
 248 r1_t = 7.88
 249 r2 = 20.0
 250 r2_t = 2.0
 251
 252 ##############
 253 # Association
 254 ##############
 255
 256 # association vertex/points de la grille 1
 257 # (tous les vertex qui ne sont pas fusionnes avec ceux de la grille #
 258 # 2)
 259
 260 dz_geom = geompy.MakeVectorDXDYDZ(0., 0., 1.)
 261
 262 # les vertex du cylindre 1 sont crees de bas en haut, en tournant dans
 263 # le sens trigonometrique :
 264
 265 # 6 vertex sont necessaires / axe z (pour associer aux 6 vertices du
 266 # modele) :
 267
 268 z_val = [-1, 1.5, 21.5, 34, 46.5]  # nl1 + 1 valeurs
 269 for ni in range(nr1 + 1):
 270     # suivant ni, valeurs des x pour les (nl1 + 1) points (selon l'axe z)
 271     x = []
 272     z = []
 273     nb_z = 0  # nombre de points suivant l'axe z
 274     if ni == 0:
 275         z = z_val
 276         x = [r1_t] * len(z)
 277         pass
 278     elif ni == 1:
 279         z = z_val
 280         x = [r1_t + 2.77] * len(z)
 281         pass
 282     elif ni == 2:
 283         z = z_val[0:-1]  # tout sauf le dernier
 284         x_last = r1_t + 2.77
 285         x = [24.0, 24.0, 19.0, (19.0 - x_last)/2 + x_last, x_last]  # pour le 4 eme point, moyenne
 286         # entre le 3 eme et le 5 eme
 287         pass
 288     elif ni == 3:
 289         z = z_val[1:3]
 290         x = [24.0, 19.0]
 291         pass
 292     elif ni == 4:
 293         z = z_val[1:3]
 294         x = [26.5, 21.0]  # a revoir pour le premier point ??
 295         pass
 296     elif ni == 8:
 297         z = z_val[1:3]
 298         x = [47.5] * 2
 299         pass
 300     else:  # ni = 5, 6, 7
 301         z = z_val[1:3]
 302         x = [26.5 + (47.5 - 26.5)/4*(ni - 4)] * 2
 303         pass
 304     pass
 305
 306     nb_z = len(z)
 307
 308     # creation des points pour y = 0 :
 309     vert_grid1_xi = [geompy.MakeVertex(xi, 0, zi) for (xi, zi) in \
 310                      zip(x, z)]
 311
 312     # les points suivants sont crees par rotation de PI/16 suivant
 313     # l'axe z / aux precedents :
 314     angle = math.pi/4.0/na1  # PI/4 (45 degres), divise par 4
 315     for j in range(na1):
 316         li = [geompy.MakeRotation(v, dz_geom, angle) for v in vert_grid1_xi[-nb_z:]]
 317         vert_grid1_xi.extend(li)
 318         pass
 319
 320     # ajout des points a l'etude et association :
 321     # les vertex fusionnes ou correspondant a des hexaedres effaces ne
 322     # sont pas pris en compte.
 323     for nj in range(na1 + 1):
 324         for nk in range(nb_z):
 325             if (ni <= 2) or (3 <= ni <= 7 and nj >= na1 - 1) or \
 326                    (ni == 8 and (nj == 0  or nj >= na1 - 1)):
 327                 v_mod = grille_cyl_quart.getVertexIJK(ni, nj, nk)
 328                 v_geo = vert_grid1_xi[nk + nj*nb_z]
 329                 geompy.addToStudy(v_geo, "vert_grid1_x" + str(ni) + "_y" + \
 330                                   str(nj) + "_z" + str(nk))
 331                 v_mod.setAssociation(v_geo)
 332                 pass
 333             pass
 334         pass
 335
 336     pass
 337
 338 # association vertex/points de la grille 2
 339 # (tous les vertex qui ne sont pas fusionnes avec ceux de la grille #
 340 # 1)
 341 ## dz_geom2 = geompy.MakeVectorDXDYDZ(33.5, 0., 1.)
 342 pt_a = geompy.MakeVertex(33.5, 0, 0)
 343 pt_b = geompy.MakeVertex(33.5, 0, 1.)
 344 dz_geom2 = geompy.MakeVector(pt_a, pt_b)
 345
 346 # les vertex du cylindre 2 sont crees de bas en haut, en tournant dans
 347 # le sens trigonometrique :
 348
 349 # REM : pour l'instant on met de cote la partie centrale du cylindre
 350 # (6 vertex selon z) => a faire a la fin. Ils sont d'ailleurs ranges
 351 # apres les autres dans le modele
 352
 353
 354 # 6 vertex sont necessaires / axe z (pour associer aux 6 vertices du
 355 # modele) :
 356
 357 z_val = [-1, 1.5, 21.5, 24, 36, 41.5]  # nl2 + 1 valeurs
 358 for ni in range(nr2 + 1):
 359     # suivant ni, valeurs des x pour les (nl1 + 1) points (selon l'axe z)
 360     x = []
 361     z = []
 362     nb_z = 0  # nombre de points suivant l'axe z
 363     if ni == 0:
 364         z = z_val
 365         x = [39.5] * len(z)
 366         pass
 367     elif ni == 1:
 368         z = z_val[1:-1]  # tout sauf le dernier et le premier
 369         x = [42.5] * len(z)
 370         pass
 371     elif ni == 2:
 372         z = z_val[1:-2]  # tout sauf les 2 derniers et le premier
 373         x = [46.] * len(z)
 374         pass
 375     elif ni == 3:
 376         z = z_val[1:3]
 377         x = [46.7] * len(z)  # valeur a revoir ??
 378         pass
 379     pass
 380
 381     nb_z = len(z)
 382
 383     # creation des points pour y = 0 :
 384     vert_grid2_xi = [geompy.MakeVertex(xi, 0, zi) for (xi, zi) in \
 385                      zip(x, z)]
 386
 387     # les points suivants sont crees par rotation de PI/16 suivant
 388     # l'axe z / aux precedents :
 389     angle = math.pi/na2  # PI (180 degres), divise par 8
 390     for j in range(na2):
 391         li = [geompy.MakeRotation(v, dz_geom2, angle) for v in vert_grid2_xi[-nb_z:]]
 392         vert_grid2_xi.extend(li)
 393         pass
 394
 395     # ajout des points a l'etude et association :
 396     for nj in range(na2 + 1):
 397         for nk in range(nb_z):
 398             v_mod = grille_cyl_demi.getVertexIJK(ni, nj, nk)
 399             v_geo = vert_grid2_xi[nk + nj*nb_z]
 400             geompy.addToStudy(v_geo, "vert_grid2_x" + str(ni) + "_y" + \
 401                               str(nj) + "_z" + str(nk))
 402             v_mod.setAssociation(v_geo)
 403             pass
 404         pass
 405
 406     pass
 407
 408 # association des vertex communs grille1/grille2
 409 # REM : cette etape n'est pas necessaire ? En effet, les vertex ayant
 410 # ete fusionnes, l'association a ete faite avec la grille 2