src/Tools/blocFissure/gmu/fissureCoude.py

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