src/Tools/blocFissure/gmu/fissureCoude.py

   1 # -*- coding: utf-8 -*-
   2 # Copyright (C) 2014-2019  CEA/DEN, EDF R&D
   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 from .geomsmesh import geompy, smesh
  22 from .geomsmesh import geomPublish
  23 from .geomsmesh import geomPublishInFather
  24 from . import initLog
  25
  26 import math
  27 import GEOM
  28 import SALOMEDS
  29 import SMESH
  30 #import StdMeshers
  31 #import GHS3DPlugin
  32 #import NETGENPlugin
  33 import logging
  34
  35 from .fissureGenerique import fissureGenerique
  36
  37 from .triedreBase import triedreBase
  38 from .genereMeshCalculZoneDefaut import genereMeshCalculZoneDefaut
  39 from .creeZoneDefautDansObjetSain import creeZoneDefautDansObjetSain
  40 from .construitFissureGenerale import construitFissureGenerale
  41 from .sortEdges import sortEdges
  42
  43 O, OX, OY, OZ = triedreBase()
  44
  45 class fissureCoude(fissureGenerique):
  46   """
  47   problème de fissure du Coude : version de base
  48   maillage hexa
  49   """
  50
  51   nomProbleme = "tuyau_Coude"
  52
  53   # ---------------------------------------------------------------------------
  54   def setParamGeometrieSaine(self):
  55     """
  56     Paramètres géométriques du tuyau coudé sain:
  57     angleCoude
  58     r_cintr
  59     l_tube_p1
  60     l_tube_p2
  61     epais
  62     de
  63     """
  64     self.geomParams = dict(angleCoude = 60,
  65                            r_cintr    = 1200,
  66                            l_tube_p1  = 1600,
  67                            l_tube_p2  = 1200,
  68                            epais      = 40,
  69                            de         = 760)
  70
  71   # ---------------------------------------------------------------------------
  72   def genereGeometrieSaine(self, geomParams):
  73     logging.info("genereGeometrieSaine %s", self.nomCas)
  74
  75     angleCoude = geomParams['angleCoude']
  76     r_cintr    = geomParams['r_cintr']
  77     l_tube_p1  = geomParams['l_tube_p1']
  78     l_tube_p2  = geomParams['l_tube_p2']
  79     epais      = geomParams['epais']
  80     de         = geomParams['de']
  81
  82     centre = geompy.MakeVertex(0, 0, -l_tube_p1)
  83     diskext = geompy.MakeDiskPntVecR(centre, OZ, de/2.)
  84     diskint = geompy.MakeDiskPntVecR(centre, OZ, de/2. -epais)
  85     couronne = geompy.MakeCut(diskext, diskint)
  86     tube_1 = geompy.MakePrismVecH(couronne, OZ, l_tube_p1)
  87     axe = geompy.MakeTranslation(OY, -r_cintr, 0, -l_tube_p1)
  88     coude = geompy.MakeRevolution(couronne, axe, angleCoude*math.pi/180.0)
  89     Rotation_1 = geompy.MakeRotation(couronne, axe, angleCoude*math.pi/180.0)
  90     Rotation_2 = geompy.MakeRotation(OZ, OY, angleCoude*math.pi/180.0)
  91     tube_2 = geompy.MakePrismVecH(Rotation_1, Rotation_2, -l_tube_p2)
  92     plan_y = geompy.MakePlaneLCS(None, 100000, 3)
  93     geomPublish(initLog.debug,  plan_y, "plan_y" )
  94     geomPublish(initLog.debug,  tube_1, "tube_1" )
  95     geomPublish(initLog.debug,  coude, "coude" )
  96     geomPublish(initLog.debug,  tube_2, "tube_2" )
  97
  98     P1 = O
  99     geompy.addToStudy(P1, "P1" )
 100     op2 = geompy.MakeVertex(0, 0, -l_tube_p1)
 101     P2 = geompy.MakeRotation(op2, axe, angleCoude*math.pi/180.0)
 102     P2 = geompy.MakeTranslationVectorDistance(P2, Rotation_2, -l_tube_p2)
 103     geompy.addToStudy(P2, "P2" )
 104
 105     # --- tube coude sain
 106
 107     geometrieSaine = geompy.MakePartition([tube_1, coude, tube_2, P1, P2], [plan_y], [], [], geompy.ShapeType["SOLID"], 0, [], 1)
 108     geomPublish(initLog.debug,  geometrieSaine, self.nomCas )
 109     [P1, P2] = geompy.RestoreGivenSubShapes(geometrieSaine, [P1, P2], GEOM.FSM_GetInPlaceByHistory, False, True)
 110
 111     xmin = -de -r_cintr -l_tube_p2
 112     zmin = -l_tube_p1 -r_cintr -l_tube_p2 -de
 113     ymax = de +100.
 114     boxypos = geompy.MakeBox(xmin, 0, zmin, ymax, ymax, 100, "boxypos")
 115     boxyneg = geompy.MakeBox(xmin, 0, zmin, ymax, -ymax, 100, "boxyneg")
 116     edgesypos = geompy.GetShapesOnShape(boxypos, geometrieSaine, geompy.ShapeType["EDGE"], GEOM.ST_IN)
 117     edgesyneg = geompy.GetShapesOnShape(boxyneg, geometrieSaine, geompy.ShapeType["EDGE"], GEOM.ST_IN)
 118     circ_g = geompy.CreateGroup(geometrieSaine, geompy.ShapeType["EDGE"])
 119     geompy.UnionList(circ_g, edgesyneg)
 120     circ_d = geompy.CreateGroup(geometrieSaine, geompy.ShapeType["EDGE"])
 121     geompy.UnionList(circ_d, edgesypos)
 122     edgesy0pos = geompy.GetShapesOnShape(boxypos, geometrieSaine, geompy.ShapeType["EDGE"], GEOM.ST_ONIN)
 123     grpedpos = geompy.CreateGroup(geometrieSaine, geompy.ShapeType["EDGE"])
 124     geompy.UnionList(grpedpos, edgesy0pos)
 125     grpedy0 = geompy.CutGroups(grpedpos, circ_d, "edges_y0")
 126     boxtub1 = geompy.MakeBox(-de/2.0 -1, -1, -l_tube_p1, de, de, 0, "boxtub1")
 127     edgestub1 = geompy.GetShapesOnShape(boxtub1, geometrieSaine, geompy.ShapeType["EDGE"], GEOM.ST_IN)
 128     grped = geompy.CreateGroup(geometrieSaine, geompy.ShapeType["EDGE"])
 129     geompy.UnionList(grped, edgestub1)
 130     long_p1 = geompy.IntersectGroups(grped, grpedy0)
 131     boxtub  = geompy.MakeBox(-de/2.0 -1, -1, -l_tube_p1 -l_tube_p2, de, de, -l_tube_p1)
 132     boxtub2 = geompy.MakeRotation(boxtub, axe, angleCoude*math.pi/180.0, "boxttub2")
 133     edgestub2 = geompy.GetShapesOnShape(boxtub2, geometrieSaine, geompy.ShapeType["EDGE"], GEOM.ST_IN)
 134     grped = geompy.CreateGroup(geometrieSaine, geompy.ShapeType["EDGE"])
 135     geompy.UnionList(grped, edgestub2)
 136     long_p2 = geompy.IntersectGroups(grped, grpedy0)
 137     boxtub1t = geompy.MakeTranslationVectorDistance(boxtub1, OZ, -l_tube_p1)
 138     facer = geompy.GetShapesOnShape(boxtub1t, boxtub1, geompy.ShapeType["FACE"], GEOM.ST_ONIN, "facer")
 139     boxcoud = geompy.MakeRevolution(facer[0], axe, angleCoude*math.pi/180.0, "boxcoud")
 140     edgescoud = geompy.GetShapesOnShape(boxcoud, geometrieSaine, geompy.ShapeType["EDGE"], GEOM.ST_IN)
 141     grped = geompy.CreateGroup(geometrieSaine, geompy.ShapeType["EDGE"])
 142     geompy.UnionList(grped, edgescoud)
 143     long_coude = geompy.IntersectGroups(grped, grpedy0)
 144     grped = geompy.CutGroups(grpedy0, long_p1)
 145     grped = geompy.CutGroups(grped, long_p2)
 146     ep = geompy.CutGroups(grped, long_coude)
 147     geomPublishInFather(initLog.debug, geometrieSaine, long_p1, 'long_p1' )
 148     geomPublishInFather(initLog.debug, geometrieSaine, ep, 'ep' )
 149     geomPublishInFather(initLog.debug, geometrieSaine, long_coude, 'long_coude' )
 150     geomPublishInFather(initLog.debug, geometrieSaine, circ_g, 'circ_g' )
 151     geomPublishInFather(initLog.debug, geometrieSaine, circ_d, 'circ_d' )
 152     geomPublishInFather(initLog.debug, geometrieSaine, long_p2, 'long_p2' )
 153
 154     # --- face extremite tube (EXTUBE)
 155
 156     facesIds = geompy.GetShapesOnPlaneIDs(geometrieSaine, geompy.ShapeType["FACE"], OZ, GEOM.ST_ON)
 157     EXTUBE = geompy.CreateGroup(geometrieSaine, geompy.ShapeType["FACE"])
 158     geompy.UnionIDs(EXTUBE, facesIds)
 159     geomPublishInFather(initLog.debug, geometrieSaine, EXTUBE, 'EXTUBE' )
 160
 161     # --- edge bord extremite tube (BORDTU)
 162
 163     edge1Ids = geompy.GetShapesOnPlaneIDs(geometrieSaine, geompy.ShapeType["EDGE"], OZ, GEOM.ST_ON)
 164     edge2Ids = geompy.GetShapesOnCylinderIDs(geometrieSaine, geompy.ShapeType["EDGE"], OZ, de/2. -epais, GEOM.ST_ON)
 165     edgesIds = []
 166     for edge in edge1Ids:
 167       if edge in edge2Ids:
 168         edgesIds.append(edge)
 169     BORDTU = geompy.CreateGroup(geometrieSaine, geompy.ShapeType["EDGE"])
 170     geompy.UnionIDs(BORDTU, edgesIds)
 171     geomPublishInFather(initLog.debug, geometrieSaine, BORDTU, 'BORDTU' )
 172
 173     # --- face origine tube (CLGV)
 174
 175     pp2 = geompy.MakeTranslationVectorDistance(P2, Rotation_2, 10)
 176     vec2 = geompy.MakeVector(P2, pp2)
 177     #geomPublish(initLog.debug, vec2, 'vec2')
 178     facesIds = geompy.GetShapesOnPlaneIDs(geometrieSaine, geompy.ShapeType["FACE"], vec2, GEOM.ST_ON)
 179     CLGV = geompy.CreateGroup(geometrieSaine, geompy.ShapeType["FACE"])
 180     geompy.UnionIDs(CLGV, facesIds)
 181     geomPublishInFather(initLog.debug, geometrieSaine, CLGV, 'CLGV' )
 182
 183     # --- peau tube interieur (PEAUINT)
 184
 185     extru1 = geompy.MakePrismVecH(diskint, OZ, l_tube_p1)
 186     revol1 = geompy.MakeRevolution(diskint, axe, angleCoude*math.pi/180.0)
 187     rot1 = geompy.MakeRotation(diskint, axe, angleCoude*math.pi/180.0)
 188     extru2 = geompy.MakePrismVecH(rot1, Rotation_2, -l_tube_p2)
 189     interne = geompy.MakeFuse(extru1, revol1)
 190     interne = geompy.MakeFuse(extru2, interne)
 191     geomPublish(initLog.debug, interne, 'interne')
 192     facesIds = geompy.GetShapesOnShapeIDs(interne, geometrieSaine, geompy.ShapeType["FACE"], GEOM.ST_ONIN)
 193     PEAUINT = geompy.CreateGroup(geometrieSaine, geompy.ShapeType["FACE"])
 194     geompy.UnionIDs(PEAUINT, facesIds)
 195     geomPublishInFather(initLog.debug, geometrieSaine, PEAUINT, 'PEAUINT' )
 196
 197     # --- peau tube exterieur (PEAUEXT)
 198
 199     Disk_3 = geompy.MakeDiskPntVecR(centre, OZ, de/2. +epais)
 200     extru1 = geompy.MakePrismVecH(Disk_3, OZ, l_tube_p1)
 201     revol1 = geompy.MakeRevolution(Disk_3, axe, angleCoude*math.pi/180.0)
 202     rot1 = geompy.MakeRotation(Disk_3, axe, angleCoude*math.pi/180.0)
 203     extru2 = geompy.MakePrismVecH(rot1, Rotation_2, -l_tube_p2)
 204     externe = geompy.MakeFuse(extru1, revol1)
 205     externe = geompy.MakeFuse(extru2, externe)
 206     geomPublish(initLog.debug, externe, 'externe')
 207     facesIds = geompy.GetShapesOnShapeIDs(externe, geometrieSaine, geompy.ShapeType["FACE"], GEOM.ST_ON)
 208     PEAUEXT = geompy.CreateGroup(geometrieSaine, geompy.ShapeType["FACE"])
 209     geompy.UnionIDs(PEAUEXT, facesIds)
 210     geomPublishInFather(initLog.debug, geometrieSaine, PEAUEXT, 'PEAUEXT' )
 211
 212     # --- solide sain
 213
 214     volIds = geompy.SubShapeAllIDs(geometrieSaine, geompy.ShapeType["SOLID"])
 215     COUDE = geompy.CreateGroup(geometrieSaine, geompy.ShapeType["SOLID"])
 216     geompy.UnionIDs(COUDE, volIds)
 217     geomPublishInFather(initLog.debug, geometrieSaine, COUDE, 'COUDSAIN' )
 218
 219     geometriesSaines = [geometrieSaine, long_p1, ep, long_coude, circ_g, circ_d, long_p2, P1, P2, EXTUBE, BORDTU, CLGV, PEAUINT, PEAUEXT, COUDE]
 220
 221     return geometriesSaines
 222
 223   # ---------------------------------------------------------------------------
 224   def setParamMaillageSain(self):
 225     self.meshParams = dict(n_long_p1    = 16,
 226                            n_ep         = 3,
 227                            n_long_coude = 15,
 228                            n_circ_g     = 20,
 229                            n_circ_d     = 20,
 230                            n_long_p2    = 12)
 231
 232   # ---------------------------------------------------------------------------
 233   def genereMaillageSain(self, geometriesSaines, meshParams):
 234     logging.info("genereMaillageSain %s", self.nomCas)
 235
 236     geometrieSaine = geometriesSaines[0]
 237     long_p1        = geometriesSaines[1]
 238     ep             = geometriesSaines[2]
 239     long_coude     = geometriesSaines[3]
 240     circ_g         = geometriesSaines[4]
 241     circ_d         = geometriesSaines[5]
 242     long_p2        = geometriesSaines[6]
 243     P1             = geometriesSaines[7]
 244     P2             = geometriesSaines[8]
 245     EXTUBE         = geometriesSaines[9]
 246     BORDTU         = geometriesSaines[10]
 247     CLGV           = geometriesSaines[11]
 248     PEAUINT        = geometriesSaines[12]
 249     PEAUEXT        = geometriesSaines[13]
 250     COUDE          = geometriesSaines[14]
 251
 252     n_long_p1    = meshParams['n_long_p1']
 253     n_ep         = meshParams['n_ep']
 254     n_long_coude = meshParams['n_long_coude']
 255     n_circ_g     = meshParams['n_circ_g']
 256     n_circ_d     = meshParams['n_circ_d']
 257     n_long_p2    = meshParams['n_long_p2']
 258
 259     maillageSain = smesh.Mesh(geometrieSaine)
 260
 261     algo3d = maillageSain.Hexahedron()
 262     algo2d = maillageSain.Quadrangle()
 263     smesh.SetName(algo3d, "algo3d_maillageSain")
 264     smesh.SetName(algo2d, "algo2d_maillageSain")
 265
 266     algo1d_long_p1 = maillageSain.Segment(geom=long_p1)
 267     hypo1d_long_p1 = algo1d_long_p1.NumberOfSegments(n_long_p1)
 268     smesh.SetName(algo1d_long_p1, "algo1d_long_p1")
 269     smesh.SetName(hypo1d_long_p1, "hypo1d_long_p1")
 270
 271     algo1d_ep = maillageSain.Segment(geom=ep)
 272     hypo1d_ep = algo1d_ep.NumberOfSegments(n_ep)
 273     smesh.SetName(algo1d_ep, "algo1d_ep")
 274     smesh.SetName(hypo1d_ep, "hypo1d_ep")
 275
 276     algo1d_long_coude = maillageSain.Segment(geom=long_coude)
 277     hypo1d_long_coude = algo1d_long_coude.NumberOfSegments(n_long_coude)
 278     smesh.SetName(algo1d_long_coude, "algo1d_long_coude")
 279     smesh.SetName(hypo1d_long_coude, "hypo1d_long_coude")
 280
 281     algo1d_circ_g = maillageSain.Segment(geom=circ_g)
 282     hypo1d_circ_g = algo1d_circ_g.NumberOfSegments(n_circ_g)
 283     smesh.SetName(algo1d_circ_g, "algo1d_circ_g")
 284     smesh.SetName(hypo1d_circ_g, "hypo1d_circ_g")
 285
 286     algo1d_circ_d = maillageSain.Segment(geom=circ_d)
 287     hypo1d_circ_d = algo1d_circ_d.NumberOfSegments(n_circ_d)
 288     smesh.SetName(algo1d_circ_d, "algo1d_circ_d")
 289     smesh.SetName(hypo1d_circ_d, "hypo1d_circ_d")
 290
 291     algo1d_long_p2 = maillageSain.Segment(geom=long_p2)
 292     hypo1d_long_p2 = algo1d_long_p2.NumberOfSegments(n_long_p2)
 293     smesh.SetName(algo1d_long_p2, "algo1d_long_p2")
 294     smesh.SetName(hypo1d_long_p2, "hypo1d_long_p2")
 295
 296     isDone = maillageSain.Compute()
 297
 298     mp1 = maillageSain.GroupOnGeom(P1,'P1',SMESH.NODE)
 299     mp2 = maillageSain.GroupOnGeom(P2,'P2',SMESH.NODE)
 300     ext = maillageSain.GroupOnGeom(EXTUBE,'EXTUBE',SMESH.FACE)
 301     btu = maillageSain.GroupOnGeom(BORDTU,'BORDTU',SMESH.EDGE)
 302     clg = maillageSain.GroupOnGeom(CLGV,'CLGV',SMESH.FACE)
 303     pei = maillageSain.GroupOnGeom(PEAUINT,'PEAUINT',SMESH.FACE)
 304     pex = maillageSain.GroupOnGeom(PEAUEXT,'PEAUEXT',SMESH.FACE)
 305     cou = maillageSain.GroupOnGeom(COUDE,'COUDSAIN',SMESH.VOLUME)
 306
 307     return [maillageSain, True] # True : maillage hexa
 308
 309   # ---------------------------------------------------------------------------
 310   def setParamShapeFissure(self):
 311     """
 312     paramètres de la fissure pour le tuyau coude
 313     profondeur  : 0 < profondeur <= épaisseur
 314     rayonPipe   : rayon du pipe correspondant au maillage rayonnant
 315     lenSegPipe  : longueur des mailles rayonnantes le long du fond de fissure (= rayonPipe par défaut)
 316     azimut      : entre 0 et 360°
 317     alpha       : 0 < alpha < angleCoude
 318     longueur    : <=2*profondeur ==> force une fissure elliptique (longueur/profondeur = grand axe/petit axe).
 319     orientation : 0° : longitudinale, 90° : circonférentielle, autre : uniquement fissures elliptiques
 320     lgInfluence : distance autour de la shape de fissure a remailler (si 0, pris égal à profondeur. A ajuster selon le maillage)
 321     elliptique  : True : fissure elliptique (longueur/profondeur = grand axe/petit axe); False : fissure longue (fond de fissure de profondeur constante, demi-cercles aux extrémites)
 322     pointIn_x   : optionnel coordonnées x d'un point dans le solide, pas trop loin du centre du fond de fissure (idem y,z)
 323     externe     : True : fissure face externe, False : fissure face interne
 324     """
 325     logging.info("setParamShapeFissure %s", self.nomCas)
 326     self.shapeFissureParams = dict(profondeur  = 10,
 327                                    rayonPipe   = 2.5,
 328                                    lenSegPipe  = 2.5,
 329                                    azimut      = 160,
 330                                    alpha       = 20,
 331                                    longueur    = 400,
 332                                    orientation = 90,
 333                                    lgInfluence = 50,
 334                                    elliptique  = False,
 335                                    externe     = True)
 336
 337   # ---------------------------------------------------------------------------
 338   def genereShapeFissure( self, geometriesSaines, geomParams, shapeFissureParams):
 339     logging.info("genereShapeFissure %s", self.nomCas)
 340     logging.info("shapeFissureParams %s", shapeFissureParams)
 341
 342     angleCoude = geomParams['angleCoude']
 343     r_cintr    = geomParams['r_cintr']
 344     l_tube_p1  = geomParams['l_tube_p1']
 345     l_tube_p2  = geomParams['l_tube_p2']
 346     epais      = geomParams['epais']
 347     de         = geomParams['de']
 348
 349     profondeur  = shapeFissureParams['profondeur']
 350     azimut      = shapeFissureParams['azimut']
 351     alpha       = shapeFissureParams['alpha']
 352     longueur    = shapeFissureParams['longueur']
 353     orientation = shapeFissureParams['orientation']
 354     externe     = shapeFissureParams['externe']
 355     lgInfluence = shapeFissureParams['lgInfluence']
 356     self.elliptique  = False
 357     if 'elliptique' in shapeFissureParams:
 358       self.elliptique = shapeFissureParams['elliptique']
 359
 360
 361
 362     azimut = -azimut # axe inverse / ASCOUF
 363     axe = geompy.MakeTranslation(OY, -r_cintr, 0, -l_tube_p1)
 364     geomPublish(initLog.debug, axe,"axe")
 365
 366     if not lgInfluence:
 367       lgInfluence = profondeur
 368
 369     if longueur > 2*profondeur:
 370       self.fissureLongue=True
 371     else:
 372       self.fissureLongue=False
 373       self.elliptique = True
 374
 375     self.circonferentielle = False
 376     self.longitudinale = False
 377     if self.fissureLongue and not self.elliptique:
 378       if abs(orientation) < 45 :
 379         self.longitudinale = True
 380       else:
 381         self.circonferentielle = True
 382
 383     nbp1 = 10
 384     if self.circonferentielle:
 385       if externe:
 386         dp = -1.0
 387         raybor = de/2.
 388         rayint = raybor - profondeur
 389         rayext = raybor + profondeur/5.0
 390       else:
 391         dp = 1.0
 392         raybor = de/2. - epais
 393         rayint = raybor + profondeur
 394         rayext = raybor - profondeur/5.0
 395       lgfond = longueur -2*profondeur
 396       angle = lgfond/(2*raybor)
 397       pb = geompy.MakeVertex(raybor, 0, 0)
 398       pi = geompy.MakeVertex(rayint, 0, 0)
 399       pbl = geompy.MakeRotation(pb, OZ, angle)
 400       pbr = geompy.MakeRotation(pb, OZ, -angle)
 401       geomPublish(initLog.debug, pbl,"pbl")
 402       geomPublish(initLog.debug, pbr,"pbr")
 403       pal = geompy.MakeTranslationVector(pbl, OZ)
 404       par = geompy.MakeTranslationVector(pbr, OZ)
 405       axl = geompy.MakeVector(pbl,pal)
 406       axr = geompy.MakeVector(pbr,par)
 407       pil = geompy.MakeRotation(pi, OZ, angle)
 408       pir = geompy.MakeRotation(pi, OZ, -angle)
 409       points = []
 410       nbp = 3*nbp1
 411       for i in range(nbp):
 412         angi = dp*(nbp -i)*(2.0*math.pi/3.0)/nbp
 413         pt = geompy.MakeRotation(pil, axl, angi)
 414         points.append(pt)
 415       for i in range(nbp):
 416         angi = angle -2.0*i*angle/nbp
 417         pt = geompy.MakeRotation(pi, OZ, angi)
 418         points.append(pt)
 419       for i in range(nbp+1):
 420         angi = -dp*i*(2.0*math.pi/3.0)/nbp
 421         pt = geompy.MakeRotation(pir, axr, angi)
 422         points.append(pt)
 423       for i, pt in enumerate(points):
 424         pt = geompy.MakeRotation(pt, OZ, azimut*math.pi/180.)
 425         pt = geompy.MakeTranslation(pt, 0, 0, -l_tube_p1)
 426         pt = geompy.MakeRotation(pt, axe, alpha*math.pi/180.)
 427         points[i] = pt
 428       wire0 = geompy.MakeInterpol(points[0:nbp+1])
 429       wire1 = geompy.MakeInterpol(points[nbp:2*nbp+1])
 430       wire2 = geompy.MakeInterpol(points[2*nbp:3*nbp+1])
 431       #wiretube = geompy.MakeInterpol(points)
 432       wiretube=geompy.MakeWire([wire0,wire1,wire2])
 433       geomPublish(initLog.debug, wiretube,"wiretube")
 434
 435       pe = geompy.MakeVertex(rayext, 0, 0)
 436       pe = geompy.MakeRotation(pe, OZ, azimut*math.pi/180.)
 437       pe = geompy.MakeTranslation(pe, 0, 0, -l_tube_p1)
 438       pe = geompy.MakeRotation(pe, axe, alpha*math.pi/180.)
 439
 440       arce = geompy.MakeArc(points[0], pe, points[-1])
 441       geomPublish(initLog.debug, arce,"arce")
 442
 443       facefiss = geompy.MakeFaceWires([arce, wiretube], 1)
 444       geomPublish(initLog.debug,  facefiss, 'facefissPlace' )
 445
 446       pc = geompy.MakeVertex((raybor + rayint)/2.0, 0, 0)
 447       centre = geompy.MakeRotation(pc, OZ, azimut*math.pi/180.)
 448       centre = geompy.MakeTranslation(centre, 0, 0, -l_tube_p1)
 449       centre = geompy.MakeRotation(centre, axe, alpha*math.pi/180.)
 450       geomPublish(initLog.debug,  centre, 'centrefissPlace' )
 451
 452       wiretube = geompy.GetInPlace(facefiss, wiretube)
 453       geomPublish(initLog.debug, wiretube, 'wiretubePlace' )
 454       try:
 455         edgetube = geompy.MakeEdgeWire(wiretube)
 456         geomPublish(initLog.debug, edgetube,"edgetube")
 457       except:
 458         logging.debug("erreur MakeEdgeWire sur fond de fissure, on fait sans")
 459         edgetube = None
 460
 461     # ---------------------------------------------------------
 462
 463     elif self.longitudinale:
 464       if externe:
 465         raybor = de/2.
 466         dp = -1.0
 467       else:
 468         raybor = de/2. - epais
 469         dp = +1.0
 470       prof = dp * profondeur
 471       lgfond = longueur -2*profondeur
 472       cosaz = math.cos(azimut*math.pi/180.)
 473       sinaz = math.sin(azimut*math.pi/180.)
 474       alfrd = alpha*math.pi/180.
 475       rayxy = r_cintr + raybor*cosaz
 476       angle = lgfond/(2.*rayxy)
 477       logging.debug("longueur: %s, angle: %s, rayon: %s",lgfond, angle, rayxy)
 478       pb = geompy.MakeVertex(raybor*cosaz, raybor*sinaz, -l_tube_p1, "pb")
 479       pi = geompy.MakeTranslation(pb, prof*cosaz, prof*sinaz, 0., "pi")
 480       pbv = geompy.MakeTranslation(pb, -sinaz, cosaz, 0., "pbv")
 481       axb  = geompy.MakeVector(pb,pbv, "axb")
 482       pbl = geompy.MakeRotation(pb, axe, alfrd -angle, "pbl")
 483       pbr = geompy.MakeRotation(pb, axe, alfrd +angle, "pbr")
 484       axl = geompy.MakeRotation(axb, axe, alfrd -angle, "axl")
 485       axr = geompy.MakeRotation(axb, axe, alfrd +angle, "axr")
 486       pil = geompy.MakeRotation(pi, axe, alfrd -angle, "pil")
 487       pir = geompy.MakeRotation(pi, axe, alfrd +angle, "pir")
 488
 489       curves = []
 490
 491       points = []
 492       nbp = 3*nbp1
 493       xs = []
 494       totx = 0
 495       for i in range(nbp+2):
 496         x = math.sin(i*math.pi/(nbp+1)) # fonction de répartition des points : distance relative
 497         x2 = x*x
 498         totx += x2
 499         xs.append(totx)
 500         logging.debug("x2: %s, totx: %s", x2, totx)
 501       for i in range(nbp+1):
 502         #posi = nbp -i             # répartition équidistante des points sur la courbe
 503         posi = nbp*(1 -xs[i]/totx) # points plus resserrés aux extrémités de la courbe
 504         angi = -dp*posi*(5.0*math.pi/8.0)/nbp
 505         pt = geompy.MakeRotation(pil, axl, angi)
 506         points.append(pt)
 507       curves.append(geompy.MakeInterpol(points))
 508       point0 = points[0]
 509       geomPublish(initLog.debug, curves[-1],"curve0")
 510 #      for i, pt in enumerate(points):
 511 #        name = "point%d"%i
 512 #        geomPublishInFather(initLog.debug,curves[-1], pt, name)
 513
 514       points = []
 515       nbp = 3*nbp1
 516       xs =[]
 517       totx = 0
 518       for i in range(nbp+1):
 519         x = math.sin(i*math.pi/nbp)
 520         #x = 1.0 # répartition équidistante des points sur la courbe
 521         x2 = x*x # points plus resserrés aux extrémités de la courbe
 522         totx += x2
 523         xs.append(totx)
 524         logging.debug("x2: %s, totx: %s", x2, totx)
 525       for i in range(nbp):
 526         angi = alfrd -angle +2.0*angle*xs[i]/totx
 527         pt = geompy.MakeRotation(pi, axe, angi)
 528         points.append(pt)
 529       curves.append(geompy.MakeInterpol(points))
 530       geomPublish(initLog.debug, curves[-1],"curve1")
 531 #      for i, pt in enumerate(points):
 532 #        name = "point%d"%i
 533 #        geomPublishInFather(initLog.debug,curves[-1], pt, name)
 534
 535       points = []
 536       nbp = 3*nbp1
 537       xs = []
 538       totx = 0
 539       for i in range(nbp+2):
 540         x = math.sin(i*math.pi/(nbp+1))
 541         x2 = x*x
 542         totx += x2
 543         xs.append(totx)
 544         logging.debug("x2: %s, totx: %s", x2, totx)
 545       for i in range(nbp+1):
 546         #posi = nbp -i        # répartition équidistante des points sur la courbe
 547         posi = nbp*xs[i]/totx # points plus resserrés aux extrémités de la courbe
 548         angi = dp*posi*(5.0*math.pi/8.0)/nbp
 549         pt = geompy.MakeRotation(pir, axr, angi)
 550         points.append(pt)
 551       curves.append(geompy.MakeInterpol(points))
 552       point1 = points[-1]
 553       geomPublish(initLog.debug, curves[-1],"curve2")
 554 #      for i, pt in enumerate(points):
 555 #        name = "point%d"%i
 556 #        geomPublishInFather(initLog.debug,curves[-1], pt, name)
 557
 558       wiretube = geompy.MakeWire(curves)
 559       geomPublish(initLog.debug, wiretube,"wiretube")
 560       try:
 561         edgetube = geompy.MakeEdgeWire(wiretube)
 562         geomPublish(initLog.debug, edgetube,"edgetube")
 563       except:
 564         logging.debug("erreur MakeEdgeWire sur fond de fissure, on fait sans")
 565         edgetube = None
 566
 567       pts = []
 568       pts.append(point0)
 569       dpr = prof*math.cos(5.0*math.pi/8.0)
 570       pe = geompy.MakeTranslation(pb, dpr*cosaz, dpr*sinaz, 0., "pe")
 571       for i in range(nbp):
 572         angi = alfrd -angle +2.0*i*angle/nbp
 573         pt = geompy.MakeRotation(pe, axe, angi)
 574         pts.append(pt)
 575       pts.append(point1)
 576       arce = geompy.MakeInterpol(pts)
 577       geomPublish(initLog.debug, arce,"arce")
 578
 579       facefiss = geompy.MakeFaceWires([arce, wiretube], 0)
 580       geomPublish(initLog.debug,  facefiss, 'facefissPlace' )
 581
 582       pc = geompy.MakeTranslation(pb, 0.5*prof*cosaz, 0.5*prof*sinaz, 0.)
 583       centre = geompy.MakeRotation(pc, axe, alfrd)
 584       geomPublish(initLog.debug,  centre, 'centrefissPlace' )
 585
 586       edges = geompy.ExtractShapes(facefiss, geompy.ShapeType["EDGE"], True)
 587       edgesTriees, minl, maxl = sortEdges(edges)
 588       edges = edgesTriees[:-1] # la plus grande correspond à arce, on l'elimine
 589       wiretube = geompy.MakeWire(edges)
 590       #wiretube = edgesTriees[-1]
 591       geomPublish(initLog.debug, wiretube, 'wiretubePlace' )
 592
 593     # ---------------------------------------------------------
 594
 595     else: # fissure elliptique, longue ou courte
 596       if externe:
 597         raybor = de/2.
 598         dp = -1.0
 599       else:
 600         raybor = de/2. - epais
 601         dp = +1.0
 602       prof = dp * profondeur
 603       cosaz = math.cos(azimut*math.pi/180.)
 604       sinaz = math.sin(azimut*math.pi/180.)
 605       alfrd = alpha*math.pi/180.
 606       pb = geompy.MakeVertex(raybor*cosaz, raybor*sinaz, -l_tube_p1, "pb")
 607       pi = geompy.MakeTranslation(pb, prof*cosaz, prof*sinaz, 0., "pi")
 608       pbv = geompy.MakeTranslation(pb, -profondeur*sinaz, profondeur*cosaz, 0., "pbv")
 609       ayb  = geompy.MakeVector(pb,pbv, "ayb")
 610       pb0 = geompy.MakeRotation(pb, axe, alfrd, "pb0")
 611       ay0 = geompy.MakeRotation(ayb, axe, alfrd, "ay0")
 612       pi0 = geompy.MakeRotation(pi, axe, alfrd, "pi0")
 613       az_ = geompy.MakeVector(pi0, pb0, "az_")
 614       az0 = geompy.MakeTranslationVector(az_, az_, "az0") #normale sortante
 615       ax0 = geompy.MakeRotation(ay0, az0, -math.pi/2.0, "ax0")
 616       ax1 = geompy.MakeRotation(ax0, az0, orientation*math.pi/180., "ax1")
 617       ay1 = geompy.MakeRotation(ay0, az0, orientation*math.pi/180., "ay1")
 618       originLCS = geompy.MakeMarker(0, 0, 0, 1, 0, 0, 0, 1, 0, "originLCS")
 619       coo = geompy.PointCoordinates(pb0)
 620       cox = geompy.VectorCoordinates(ax1)
 621       coy = geompy.VectorCoordinates(ay1)
 622       localLCS = geompy.MakeMarker(coo[0], coo[1], coo[2], cox[0], cox[1], cox[2], coy[0], coy[1], coy[2], "localLCS")
 623
 624       pco = geompy.MakeVertex(0, 0, -profondeur, "pco")
 625       pao = geompy.MakeRotation(pco, OY, 0.6*math.pi, "pao")
 626       pbo = geompy.MakeRotation(pco, OY, -0.6*math.pi, "pbo")
 627       pce = geompy.MakeVertex(0, 0, 0.1*profondeur,"pce")
 628       arcoo = geompy.MakeArc(pao, pco, pbo, "arcoo")
 629       linoo = geompy.MakeArc(pao, pce, pbo, "linoo")
 630       scalex = longueur/profondeur
 631       arco =geompy.MakeScaleAlongAxes(arcoo, O, scalex, 1., 1., "arco")
 632       lino =geompy.MakeScaleAlongAxes(linoo, O, scalex, 1., 1., "lino")
 633       arci = geompy.MakePosition(arco, originLCS, localLCS, "arci")
 634       arce = geompy.MakePosition(lino, originLCS, localLCS, "arce")
 635       facefiss = geompy.MakeFaceWires([arce, arci], 0)
 636       geomPublish(initLog.debug,  facefiss, 'facefissPlace' )
 637       edges = geompy.ExtractShapes(facefiss, geompy.ShapeType["EDGE"], True)
 638       edgesTriees, minl, maxl = sortEdges(edges)
 639       edgetube = edgesTriees[-1] # la plus grande correspond à arci
 640       wiretube = edgetube
 641
 642       pc = geompy.MakeTranslation(pb, 0.5*prof*cosaz, 0.5*prof*sinaz, 0.)
 643       centre = geompy.MakeRotation(pc, axe, alfrd)
 644       geomPublish(initLog.debug,  centre, 'centrefissPlace' )
 645
 646     coordsNoeudsFissure = genereMeshCalculZoneDefaut(facefiss, profondeur/2. ,profondeur)
 647
 648     return [facefiss, centre, lgInfluence, coordsNoeudsFissure, wiretube, edgetube]
 649
 650   # ---------------------------------------------------------------------------
 651   def setParamMaillageFissure(self):
 652     """
 653     Paramètres du maillage de la fissure pour le tuyau coudé
 654     Voir également setParamShapeFissure, paramètres rayonPipe et lenSegPipe.
 655     nbSegRad = nombre de couronnes
 656     nbSegCercle = nombre de secteurs
 657     areteFaceFissure = taille cible de l'arête des triangles en face de fissure.
 658     """
 659     self.maillageFissureParams = dict(nomRep        = '.',
 660                                       nomFicSain    = self.nomCas,
 661                                       nomFicFissure = 'fissure_' + self.nomCas,
 662                                       nbsegRad      = 5,
 663                                       nbsegCercle   = 6,
 664                                       areteFaceFissure = 5)
 665
 666   # ---------------------------------------------------------------------------
 667   def genereZoneDefaut(self, geometriesSaines, maillagesSains, shapesFissure, shapeFissureParams, maillageFissureParams):
 668     elementsDefaut = creeZoneDefautDansObjetSain(geometriesSaines, maillagesSains, shapesFissure, shapeFissureParams, maillageFissureParams)
 669     return elementsDefaut
 670
 671   # ---------------------------------------------------------------------------
 672   def genereMaillageFissure(self, geometriesSaines, maillagesSains,
 673                             shapesFissure, shapeFissureParams,
 674                             maillageFissureParams, elementsDefaut, step):
 675     maillageFissure = construitFissureGenerale(maillagesSains,
 676                                                shapesFissure, shapeFissureParams,
 677                                                maillageFissureParams, elementsDefaut, step)
 678     return maillageFissure
 679
 680   # ---------------------------------------------------------------------------
 681   def setReferencesMaillageFissure(self):
 682     self.referencesMaillageFissure = dict(Entity_Node            = 77917,
 683                                           Entity_Quad_Edge       = 975,
 684                                           Entity_Quad_Triangle   = 2182,
 685                                           Entity_Quad_Quadrangle = 6842,
 686                                           Entity_Quad_Tetra      = 20135,
 687                                           Entity_Quad_Hexa       = 8994,
 688                                           Entity_Quad_Penta      = 972,
 689                                           Entity_Quad_Pyramid    = 1038)
 690