Salome HOME
ménage
[modules/smesh.git] / src / Tools / blocFissure / gmu / insereFissureGenerale.py
1 # -*- coding: utf-8 -*-
2 # Copyright (C) 2014-2021  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 """procédure complète fissure générale"""
21
22 import os
23 import math
24
25 import logging
26 from . import initLog
27
28 import salome
29 from salome.smesh import smeshBuilder
30 import GEOM
31 import SMESH
32
33 from .geomsmesh import geompy
34 from .geomsmesh import geomPublish
35 from .geomsmesh import geomPublishInFather
36 from .geomsmesh import smesh
37
38 from .extractionOrientee import extractionOrientee
39 from .extractionOrienteeMulti import extractionOrienteeMulti
40 from .sortFaces import sortFaces
41 from .sortEdges import sortEdges
42 from .substractSubShapes import substractSubShapes
43 from .produitMixte import produitMixte
44 from .findWireEndVertices import findWireEndVertices
45 from .findWireIntermediateVertices import findWireIntermediateVertices
46 from .orderEdgesFromWire import orderEdgesFromWire
47 from .putName import putName
48 from .enleveDefaut import enleveDefaut
49 from .shapeSurFissure import shapeSurFissure
50 from .regroupeSainEtDefaut import RegroupeSainEtDefaut
51 from .triedreBase import triedreBase
52 from .checkDecoupePartition import checkDecoupePartition
53 from .whichSide import whichSide
54 from .whichSideVertex import whichSideVertex
55 from .projettePointSurCourbe import projettePointSurCourbe
56 from .prolongeWire import prolongeWire
57
58 def insereFissureGenerale(maillagesSains,
59                           shapesFissure, shapeFissureParams,
60                           maillageFissureParams, elementsDefaut, \
61                           step=-1, mailleur="MeshGems"):
62   """ TODO: a completer"""
63   logging.info('start')
64
65   shapeDefaut       = shapesFissure[0] # faces de fissure, débordant
66   fondFiss          = shapesFissure[4] # groupe d'edges de fond de fissure
67
68   rayonPipe = shapeFissureParams['rayonPipe']
69   if 'lenSegPipe' in shapeFissureParams:
70     lenSegPipe = shapeFissureParams['lenSegPipe']
71   else:
72     lenSegPipe = rayonPipe
73
74   nomRep            = maillageFissureParams['nomRep']
75   nomFicSain        = maillageFissureParams['nomFicSain']
76   nomFicFissure     = maillageFissureParams['nomFicFissure']
77
78   nbsegRad          = maillageFissureParams['nbsegRad']      # nombre de couches selon un rayon du pipe
79   nbsegCercle       = maillageFissureParams['nbsegCercle']   # nombre de secteur dans un cercle du pipe
80   areteFaceFissure  = maillageFissureParams['areteFaceFissure']
81
82   pointIn_x = 0.0
83   pointIn_y = 0.0
84   pointIn_z = 0.0
85   isPointInterne = False
86   if 'pointIn_x' in shapeFissureParams:
87     pointIn_x = shapeFissureParams['pointIn_x']
88     isPointInterne = True
89   if 'pointIn_y' in shapeFissureParams:
90     pointIn_y = shapeFissureParams['pointIn_y']
91     isPointInterne = True
92   if 'pointIn_z' in shapeFissureParams:
93     pointIn_z = shapeFissureParams['pointIn_z']
94     isPointInterne = True
95   if isPointInterne:
96     pointInterne = geompy.MakeVertex(pointIn_x, pointIn_y, pointIn_z)
97
98   #fichierMaillageSain = os.path.join(nomRep, '{}.med'.format(nomFicSain))
99   fichierMaillageFissure = os.path.join(nomRep, '{}.med'.format(nomFicFissure))
100
101   # fillings des faces en peau
102   facesDefaut = elementsDefaut[0]
103   #centresDefaut            = elementsDefaut[1]
104   #normalsDefaut            = elementsDefaut[2]
105   #extrusionsDefaut         = elementsDefaut[3]
106   dmoyen                   = elementsDefaut[4]
107   bordsPartages = elementsDefaut[5]
108   #fillconts                = elementsDefaut[6]
109   #idFilToCont              = elementsDefaut[7]
110   maillageSain             = elementsDefaut[8]
111   internalBoundary         = elementsDefaut[9]
112   zoneDefaut               = elementsDefaut[10]
113   zoneDefaut_skin          = elementsDefaut[11]
114   zoneDefaut_internalFaces = elementsDefaut[12]
115   zoneDefaut_internalEdges = elementsDefaut[13]
116   #edgeFondExt              = elementsDefaut[14]
117   centreFondFiss           = elementsDefaut[15]
118   tgtCentre                = elementsDefaut[16]
119
120   # --- restriction de la face de fissure au domaine solide :
121   #     partition face fissure étendue par fillings, on garde la plus grande face
122
123   partShapeDefaut = geompy.MakePartition([shapeDefaut], facesDefaut, list(), list(), geompy.ShapeType["FACE"], 0, [], 0)
124   geomPublish(initLog.debug, partShapeDefaut, 'partShapeDefaut')
125   facesPartShapeDefaut = geompy.ExtractShapes(partShapeDefaut, geompy.ShapeType["FACE"], False)
126   if isPointInterne:
127     distfaces = [(geompy.MinDistance(face,pointInterne), i, face) for i, face in enumerate(facesPartShapeDefaut)]
128     distfaces.sort()
129     logging.debug("selection de la face la plus proche du point interne, distance=%s",distfaces[0][0])
130     facesPortFissure = distfaces[0][2]
131   else:
132     facesPartShapeDefautSorted, minSurf, maxSurf = sortFaces(facesPartShapeDefaut) # la face de fissure dans le volume doit être la plus grande
133     logging.debug("surfaces faces fissure étendue, min %s, max %s", minSurf, maxSurf)
134     facesPortFissure = facesPartShapeDefautSorted[-1] #= global
135
136   geomPublish(initLog.debug, facesPortFissure, "facesPortFissure")
137
138   O, _, _, _ = triedreBase()
139
140   # -----------------------------------------------------------------------------
141   # --- pipe de fond de fissure, prolongé, partition face fissure par pipe
142   #     identification des edges communes pipe et face fissure
143
144   if geompy.NumberOfFaces(shapeDefaut) == 1:
145     plan = geompy.MakePlane(centreFondFiss, tgtCentre, 10000)
146     shapeDefaut = geompy.MakePartition([shapeDefaut], [plan], [], [], geompy.ShapeType["FACE"], 0, [], 0) #= local
147     #fondFissCoupe = geompy.GetInPlaceByHistory(shapeDefaut, fondFiss) #= inutile
148     geomPublish(initLog.debug, shapeDefaut, 'shapeDefaut_coupe')
149     #geomPublishInFather(initLog.debug,shapeDefaut, fondFissCoupe, 'fondFiss_coupe')
150
151   extrem, norms = findWireEndVertices(fondFiss, True)
152   logging.debug("extrem: %s, norm: %s",extrem, norms)
153   cercle = geompy.MakeCircle(extrem[0], norms[0], rayonPipe)
154   cercle = geompy.MakeRotation(cercle, norms[0], math.pi/3.0 ) # éviter d'avoir l'arête de couture du pipe presque confondue avec la face fissure
155   geomPublish(initLog.debug, cercle, 'cercle')
156   fondFissProlonge = prolongeWire(fondFiss, extrem, norms, 2*rayonPipe)
157   pipeFiss = geompy.MakePipe(cercle, fondFissProlonge)
158   geomPublish(initLog.debug, pipeFiss, 'pipeFiss')
159   partFissPipe = geompy.MakePartition([shapeDefaut, pipeFiss], [], [], [], geompy.ShapeType["FACE"], 0, [], 1)
160   geomPublish(initLog.debug, partFissPipe, 'partFissPipe')
161   fissPipe = geompy.GetInPlaceByHistory(partFissPipe, shapeDefaut) #= global
162   geomPublish(initLog.debug, fissPipe, 'fissPipe')
163   partPipe = geompy.GetInPlaceByHistory(partFissPipe, pipeFiss) #= local
164   geomPublish(initLog.debug, partPipe, 'partPipe')
165
166   edgesPipeFiss = geompy.GetSharedShapesMulti([fissPipe, partPipe], geompy.ShapeType["EDGE"]) #= global
167   for i, edge in enumerate(edgesPipeFiss):
168     name = "edgePipe%d"%i
169     geomPublishInFather(initLog.debug,fissPipe, edge, name)
170   try:
171     wirePipeFiss = geompy.MakeWire(edgesPipeFiss) #= global
172   except:
173     wirePipeFiss = geompy.MakeCompound(edgesPipeFiss)
174     logging.debug("wirePipeFiss construit sous forme de compound")
175   geomPublish(initLog.debug, wirePipeFiss, "wirePipeFiss")
176
177   wireFondFiss = geompy.GetInPlace(partFissPipe,fondFiss)
178   edgesFondFiss = geompy.GetSharedShapesMulti([fissPipe, wireFondFiss], geompy.ShapeType["EDGE"])
179   for i, edge in enumerate(edgesFondFiss):
180     name = "edgeFondFiss%d"%i
181     geomPublishInFather(initLog.debug,fissPipe, edge, name)
182   wireFondFiss = geompy.MakeWire(edgesFondFiss) #= global
183   geomPublish(initLog.debug, wireFondFiss,"wireFondFiss")
184
185   # -----------------------------------------------------------------------------
186   # --- peau et face de fissure
187   #
188   # --- partition peau défaut - face de fissure prolongée - wire de fond de fissure prolongée
189   #     il peut y avoir plusieurs faces externes, dont certaines sont découpées par la fissure
190   #     liste de faces externes : facesDefaut
191   #     liste de partitions face externe - fissure : partitionPeauFissFond (None quand pas d'intersection)
192
193   partitionsPeauFissFond = list() #= global
194   ipart = 0
195   for filling in facesDefaut:
196     part = geompy.MakePartition([fissPipe, filling], list(), list(), list(), geompy.ShapeType["FACE"], 0, list(), 0)
197     isPart = checkDecoupePartition([fissPipe, filling], part)
198     if isPart: # on recrée la partition avec toutes les faces filling en outil pour avoir une face de fissure correcte
199       otherFD = [fd for fd in facesDefaut if fd != filling]
200       if len(otherFD) > 0:
201         fissPipePart = geompy.MakePartition([fissPipe], otherFD, list(), list(), geompy.ShapeType["FACE"], 0, list(), 0)
202       else:
203         fissPipePart = fissPipe
204       part = geompy.MakePartition([fissPipePart, filling], list(), list(), list(), geompy.ShapeType["FACE"], 0, list(), 0)
205       partitionsPeauFissFond.append(part)
206       geomPublish(initLog.debug,  part, 'partitionPeauFissFond%d'%ipart )
207     else:
208       partitionsPeauFissFond.append(None)
209     ipart = ipart +1
210
211
212   # --- arêtes vives détectées (dans quadranglesToShapeNoCorner
213   #                             et quadranglesToShapeWithCorner)
214
215   aretesVives = list()
216   aretesVivesCoupees = list()  #= global
217   ia = 0
218   for a in bordsPartages:
219     if not isinstance(a, list):
220       aretesVives.append(a)
221       name = "areteVive%d"%ia
222       geomPublish(initLog.debug, a, name)
223       ia += 1
224     else:
225       if a[0] is not None:
226         aretesVives.append(a[0])
227         name = "areteVive%d"%ia
228         geomPublish(initLog.debug, a[0], name)
229         ia += 1
230
231   aretesVivesC = None #= global
232   if len(aretesVives) > 0:
233     aretesVivesC =geompy.MakeCompound(aretesVives)
234
235   # -------------------------------------------------------
236   # --- inventaire des faces de peau coupées par la fissure
237   #     pour chaque face de peau : 0, 1 ou 2 faces débouchante du fond de fissure
238   #                                0, 1 ou plus edges de la face de fissure externe au pipe
239
240   nbFacesFilling = len(partitionsPeauFissFond)
241   ptEdgeFond = [ list()  for i in range(nbFacesFilling)] # pour chaque face [points edge fond de fissure aux débouchés du pipe]
242   fsPipePeau = [ list()  for i in range(nbFacesFilling)] # pour chaque face [faces du pipe débouchantes]
243   edRadFPiPo = [ list()  for i in range(nbFacesFilling)] # pour chaque face [edge radiale des faces du pipe débouchantes ]
244   fsFissuExt = [ list()  for i in range(nbFacesFilling)] # pour chaque face [faces de fissure externes au pipe]
245   edFisExtPe = [ list()  for i in range(nbFacesFilling)] # pour chaque face [edge en peau des faces de fissure externes (pas subshape facePeau)]
246   edFisExtPi = [ list()  for i in range(nbFacesFilling)] # pour chaque face [edge commun au pipe des faces de fissure externes]
247   facesPeaux = [None for i in range(nbFacesFilling)] # pour chaque face : la face de peau finale a mailler (percée des faces débouchantes)
248   edCircPeau = [ list()  for i in range(nbFacesFilling)] # pour chaque face de peau : [subshape edge circulaire aux débouchés du pipe]
249   ptCircPeau = [ list()  for i in range(nbFacesFilling)] # pour chaque face de peau : [subshape point sur edge circulaire aux débouchés du pipe]
250   gpedgeBord = [None for i in range(nbFacesFilling)] # pour chaque face de peau : groupe subshape des edges aux bords liés à la partie saine
251   gpedgeVifs = [None for i in range(nbFacesFilling)] # pour chaque face de peau : groupes subshape des edges aux arêtes vives entre fillings
252   edFissPeau = [ list()  for i in range(nbFacesFilling)] # pour chaque face de peau : [subshape edge en peau des faces de fissure externes]
253   ptFisExtPi = [ list()  for i in range(nbFacesFilling)] # pour chaque face de peau : [point commun edFissPeau edCircPeau]
254
255   for ifil, partitionPeauFissFond in enumerate(partitionsPeauFissFond):
256     if partitionPeauFissFond is not None:
257       fillingFaceExterne = facesDefaut[ifil]
258       #fillingSansDecoupe = fillconts[idFilToCont[ifil]]
259       logging.debug("traitement partitionPeauFissFond %s", ifil)
260       # -----------------------------------------------------------------------
261       # --- identification edges fond de fissure, edges pipe sur la face de fissure,
262       #     edges prolongées
263
264       edgesPipeC = geompy.GetInPlace(partitionPeauFissFond, geompy.MakeCompound(edgesPipeFiss)) #= local
265       geomPublishInFather(initLog.debug,partitionPeauFissFond, edgesPipeC, "edgesPipeFiss")
266       edgesFondC = geompy.GetInPlace(partitionPeauFissFond, geompy.MakeCompound(edgesFondFiss)) #= local
267       geomPublishInFather(initLog.debug,partitionPeauFissFond, edgesFondC, "edgesFondFiss")
268
269       if aretesVivesC is None: #= global facesInside facesOnside
270         [edgesInside, _, _] = extractionOrientee(fillingFaceExterne, partitionPeauFissFond, centreFondFiss, "EDGE", 1.e-3)
271         [facesInside, _, facesOnside] = extractionOrientee(fillingFaceExterne, partitionPeauFissFond, centreFondFiss, "FACE", 1.e-3)
272       else:
273         [edgesInside, _, _] = extractionOrienteeMulti(facesDefaut, ifil, partitionPeauFissFond, centreFondFiss, "EDGE", 1.e-3)
274         [facesInside, _, facesOnside] = extractionOrienteeMulti(facesDefaut, ifil, partitionPeauFissFond, centreFondFiss, "FACE", 1.e-3)
275
276       edgesPipeIn = geompy.GetSharedShapesMulti([edgesPipeC, geompy.MakeCompound(edgesInside)], geompy.ShapeType["EDGE"]) #= global
277       verticesPipePeau = list() #= global
278
279       for i, edge in enumerate(edgesPipeIn):
280         try:
281           vertices = geompy.GetSharedShapesMulti([edge, geompy.MakeCompound(facesOnside)], geompy.ShapeType["VERTEX"])
282           verticesPipePeau.append(vertices[0])
283           name = "edgePipeIn%d"%i
284           geomPublishInFather(initLog.debug,partitionPeauFissFond, edge, name)
285           name = "verticePipePeau%d"%i
286           geomPublishInFather(initLog.debug,partitionPeauFissFond, vertices[0], name)
287           logging.debug("edgePipeIn%s coupe les faces OnSide", i)
288         except:
289           logging.debug("edgePipeIn%s ne coupe pas les faces OnSide", i)
290       #edgesFondOut = list()     #= inutile
291       edgesFondIn =list() #= global
292       if len(verticesPipePeau) > 0: # au moins une extrémité du pipe sur cette face de peau
293         #tmp = geompy.GetSharedShapesMulti([edgesFondC, geompy.MakeCompound(edgesOutside)], geompy.ShapeType["EDGE"])
294         #edgesFondOut = [ ed for ed in tmp if geompy.MinDistance(ed, geompy.MakeCompound(facesOnside)) < 1.e-3]
295         tmp = geompy.GetSharedShapesMulti([edgesFondC, geompy.MakeCompound(edgesInside)], geompy.ShapeType["EDGE"])
296         edgesFondIn = [ ed for ed in tmp if geompy.MinDistance(ed, geompy.MakeCompound(facesOnside)) < 1.e-3]
297
298       verticesEdgesFondIn = list() # les points du fond de fissure au débouché du pipe sur la peau (indice de edgesFondIn)
299       pipexts = list()             # les segments de pipe associés au points de fond de fissure débouchants (même indice)
300       cercles = list()             # les cercles de generation des pipes débouchant (même indice)
301       facesFissExt = list()        # les faces de la fissure externe associés au points de fond de fissure débouchants (même indice)
302       edgesFissExtPeau = list()    # edges des faces de fissure externe sur la peau (même indice)
303       edgesFissExtPipe = list()    # edges des faces de fissure externe sur le pipe (même indice)
304       #logging.debug("edgesFondIn %s", edgesFondIn)
305
306       edgesFondFiss, edgesIdByOrientation = orderEdgesFromWire(wireFondFiss)
307       for i,edge in enumerate(edgesFondFiss):
308         geomPublishInFather(initLog.debug,wireFondFiss, edge, "edgeFondFiss%d"%i)
309
310       for iedf, edge in enumerate(edgesFondIn):
311         name = "edgeFondIn%d"%iedf
312         geomPublishInFather(initLog.debug,partitionPeauFissFond, edge, name)
313         dist = [ geompy.MinDistance(pt, edge) for pt in verticesPipePeau]
314         ptPeau = verticesPipePeau[dist.index(min(dist))] # le point de verticesPipePeau a distance minimale de l'edge
315         [parametre, PointOnEdge, EdgeInWireIndex]  = geompy.MakeProjectionOnWire(ptPeau, wireFondFiss)
316         logging.debug("parametre:%s, EdgeInWireIndex: %s, len(edgesFondFiss): %s", parametre, EdgeInWireIndex, len(edgesFondFiss))
317         localEdgeInFondFiss = edgesFondFiss[EdgeInWireIndex]
318         centre = PointOnEdge
319         centre2 = geompy.MakeVertexOnCurve(localEdgeInFondFiss, parametre)
320         geomPublishInFather(initLog.debug,partitionPeauFissFond, centre2, "centre2_%d"%iedf)
321         verticesEdgesFondIn.append(centre)
322         name = "verticeEdgesFondIn%d"%iedf
323         geomPublishInFather(initLog.debug,partitionPeauFissFond, centre, name)
324         norm = geompy.MakeTangentOnCurve(localEdgeInFondFiss, parametre)
325         geomPublishInFather(initLog.debug,partitionPeauFissFond, centre, "norm%d"%iedf)
326         cercle = geompy.MakeCircle(centre, norm, rayonPipe)
327         geomPublishInFather(initLog.debug,partitionPeauFissFond, cercle, "cerclorig%d"%iedf)
328         [vertex] = geompy.ExtractShapes(cercle, geompy.ShapeType["VERTEX"], False)
329         vec1 = geompy.MakeVector(centre, vertex)
330         vec2 = geompy.MakeVector(centre, ptPeau)
331         angle = geompy.GetAngleRadians(vec1, vec2)
332         # cas général : on reconstitue une portion de pipe, avec l'arête de couture qui coincide
333         #   avec la face de fissure, au niveau du débouché sur la face externe
334         # cas dégénéré : le pipe débouche perpendiculairement à une surface plane à l'origine.
335         #   La partition filling / pipe reconstruit échoue.
336         #   - Si on partitionne le filling avec un simple pipe obtenu par extrusion droite du cercle,
337         #     cela donne un point en trop sur le cercle.
338         #   - Si on prend une vraie surface plane (pas un filling), on peut faire la partition avec
339         #     les pipes reconstruits
340         logging.debug("angle=%s", angle)
341         #if abs(angle) > 1.e-7:
342         sommetAxe = geompy.MakeTranslationVector(centre, norm)
343         pm = produitMixte(centre, vertex, ptPeau, sommetAxe)
344         if pm > 0:  # ajout de pi a (-)angle pour éviter des points confondus (partition échoue) dans les cas dégénérés
345           cercle = geompy.MakeRotation(cercle, norm, angle + math.pi)
346         else:
347           cercle = geompy.MakeRotation(cercle, norm, -angle + math.pi)
348         name = "cercle%d"%iedf
349         geomPublishInFather(initLog.debug,partitionPeauFissFond, cercle, name)
350         cercles.append(cercle)
351
352         # --- estimation de la longueur du pipe necessaire de part et d'autre du point de sortie
353         if aretesVivesC is None:
354           faceTestPeau = fillingFaceExterne
355         else:
356           faceTestPeau = facesDefaut[ifil]
357         sideCentre = whichSide(faceTestPeau, centre)
358         locPt0 = geompy.MakeVertexOnCurve(localEdgeInFondFiss, 0.0)
359         locPt1 = geompy.MakeVertexOnCurve(localEdgeInFondFiss, 1.0)
360         sidePt0 = whichSide(faceTestPeau, locPt0)
361         sidePt1 = whichSide(faceTestPeau, locPt1)
362         logging.debug("position centre cercle: %s, extremité edge u0: %s, u1: %s", sideCentre, sidePt0, sidePt1)
363         normFace = geompy.GetNormal(faceTestPeau, ptPeau)
364         inclPipe = abs(geompy.GetAngleRadians(norm, normFace))
365         lgp = max(rayonPipe/2., abs(3*rayonPipe*math.tan(inclPipe)))
366         logging.debug("angle inclinaison Pipe en sortie: %s degres, lgp: %s", inclPipe*180/math.pi, lgp)
367
368         # --- position des points extremite du pipe sur l'edge debouchante
369         #     il faut la distance curviligne ofp du point central par rapport à une extrémité de l'edge débouchante
370         locEdgePart = geompy.MakePartition([localEdgeInFondFiss],[centre], list(), list(), geompy.ShapeType["EDGE"], 0, list(), 0)
371         edgesLoc = geompy.ExtractShapes(locEdgePart, geompy.ShapeType["EDGE"], False)
372         edgesLocSorted =[(geompy.MinDistance(edge, locPt0), kk, edge) for kk, edge in enumerate(edgesLoc)]
373         edgesLocSorted.sort()
374         ofp = geompy.BasicProperties(edgesLocSorted[0][2])[0] # distance curviligne centre locPt0
375         logging.debug("distance curviligne centre extremite0: %s", ofp)
376         p1 = geompy.MakeVertexOnCurveByLength(localEdgeInFondFiss, ofp +lgp, locPt0)
377         p2 = geompy.MakeVertexOnCurveByLength(localEdgeInFondFiss, ofp -lgp, locPt0)
378         geomPublishInFather(initLog.debug,wireFondFiss, p1, "p1_%d"%iedf)
379         geomPublishInFather(initLog.debug,wireFondFiss, p2, "p2_%d"%iedf)
380
381         edgePart = geompy.MakePartition([localEdgeInFondFiss], [p1,p2], list(), list(), geompy.ShapeType["EDGE"], 0, list(), 0)
382         edps = geompy.ExtractShapes(edgePart, geompy.ShapeType["EDGE"], True)
383         for edp in edps:
384           if geompy.MinDistance(centre, edp) < 1.e-3:
385             pipext = geompy.MakePipe(cercle, edp)
386             name = "pipeExt%d"%iedf
387             geomPublishInFather(initLog.debug,partitionPeauFissFond, pipext, name)
388             pipexts.append(pipext)
389
390         for face in facesInside:
391           logging.debug("recherche edges communes entre une face inside et (faces onside, edges pipe et fond débouchante)")
392           edgesPeauFis = list()
393           edgesPipeFis = list()
394           edgesPipeFnd = list()
395           try:
396             edgesPeauFis = geompy.GetSharedShapesMulti([geompy.MakeCompound(facesOnside), face], geompy.ShapeType["EDGE"])
397             logging.debug("    faces onside %s",edgesPeauFis)
398             edgesPipeFis = geompy.GetSharedShapesMulti([geompy.MakeCompound(edgesPipeIn), face], geompy.ShapeType["EDGE"])
399             logging.debug("    edgesPipeIn %s", edgesPipeFis)
400             edgesPipeFnd = geompy.GetSharedShapesMulti([geompy.MakeCompound(edgesFondIn), face], geompy.ShapeType["EDGE"])
401             logging.debug("    edgesFondIn %s ", edgesPipeFnd)
402           except:
403             logging.debug("  pb edges communes %s %s %s",edgesPeauFis, edgesPipeFis, edgesPipeFnd)
404           if (len(edgesPeauFis) > 0) and (len(edgesPipeFis) > 0) and (len(edgesPipeFnd) == 0):
405             dist = geompy.MinDistance(geompy.MakeCompound(edgesPeauFis), ptPeau)
406             logging.debug("    test distance extrémité reference %s", dist)
407             if dist < 1.e-3: # c'est la face de fissure externe associée
408               logging.debug("    face %s inside ajoutée", i)
409               facesFissExt.append(face)
410               name="faceFissExt%d"%iedf
411               geomPublishInFather(initLog.debug,partitionPeauFissFond, face, name)
412               dist = 1.
413               for _, edpe in enumerate(edgesPeauFis):
414                 for _, edpi in enumerate(edgesPipeFis):
415                   dist = geompy.MinDistance(edpe, edpi)
416                   if dist < 1.e-3:
417                     edgesFissExtPeau.append(edpe)
418                     name="edgesFissExtPeau%d"%iedf
419                     geomPublishInFather(initLog.debug,partitionPeauFissFond, edpe, name)
420                     edgesFissExtPipe.append(edpi)
421                     name="edgesFissExtPipe%d"%iedf
422                     geomPublishInFather(initLog.debug,partitionPeauFissFond, edpi, name)
423                     break
424                 if dist < 1.e-3:
425                   break
426
427       if len(verticesPipePeau) == 0: # aucune extrémité du pipe sur cette face de peau
428                                      # il faut recenser les edges de fissure sur la face de peau
429         j = 0
430         for face in facesInside:
431           edgesPeauFis = list()
432           edgesPipeFis = list()
433           edgesPipeFnd = list()
434           try:
435             edgesPeauFis = geompy.GetSharedShapesMulti([geompy.MakeCompound(facesOnside), face], geompy.ShapeType["EDGE"])
436             edgesPipeFis = geompy.GetSharedShapesMulti([geompy.MakeCompound(edgesPipeIn), face], geompy.ShapeType["EDGE"])
437             edgesPipeFnd = geompy.GetSharedShapesMulti([geompy.MakeCompound(edgesFondIn), face], geompy.ShapeType["EDGE"])
438           except:
439             pass
440           if (len(edgesPeauFis) > 0) and (len(edgesPipeFis) > 0) and (len(edgesPipeFnd) == 0):
441             edgesFissExtPeau.append(edgesPeauFis[0])
442             name="edgesFissExtPeau%d"%j
443             geomPublishInFather(initLog.debug,partitionPeauFissFond, edgesPeauFis[0], name)
444             j += 1
445
446       # -----------------------------------------------------------------------
447       # --- identification faces de peau : face de peau percée du pipe, extrémités du pipe
448       #     La partition avec le pipe peut créer un vertex (et un edge) de trop sur le cercle projeté,
449       #     quand le cercle est très proche de la face.
450       #     dans ce cas, la projection du cercle sur la face suivie d'une partition permet
451       #     d'éviter le point en trop
452
453       facesAndFond = facesOnside
454       facesAndFond.append(wireFondFiss)
455       try:
456         partitionPeauByPipe = geompy.MakePartition(facesAndFond, pipexts, list(), list(), geompy.ShapeType["FACE"], 0, list(), 1)
457       except:
458         logging.debug("probleme partition face pipe, contournement avec MakeSection")
459         sections = list()
460         for pipext in pipexts:
461           sections.append(geompy.MakeSection(facesOnside[0], pipext))
462         partitionPeauByPipe = geompy.MakePartition(facesAndFond, sections, list(), list(), geompy.ShapeType["FACE"], 0, list(), 1)
463
464       # contrôle edge en trop sur edges circulaires
465       if len(verticesPipePeau) > 0: # --- au moins une extrémité du pipe sur cette face de peau
466         edgeEnTrop = list()
467         outilPart = pipexts
468         facesPeau = geompy.ExtractShapes(partitionPeauByPipe, geompy.ShapeType["FACE"], False)
469         facesPeauSorted, _, _ = sortFaces(facesPeau)
470         for i, face in enumerate(facesPeauSorted[:-1]): # on ne teste que la ou les petites faces "circulaires"
471           nbv = geompy.NumberOfEdges(face)
472           logging.debug("nombre d'edges sur face circulaire: %s", nbv)
473           if nbv > 3:
474             edgeEnTrop.append(True) # TODO : distinguer les cas avec deux faces circulaires dont l'une est correcte
475           else:
476             edgeEnTrop.append(False)
477         refaire = sum(edgeEnTrop)
478         if refaire > 0:
479           dc = [(geompy.MinDistance(verticesEdgesFondIn[0], fac), i)  for i, fac in enumerate(facesPeauSorted[:-1])]
480           dc.sort()
481           logging.debug("dc sorted: %s", dc)
482           i0 = dc[0][1] # indice de facesPeauSorted qui correspond à verticesEdgesFondIn[0], donc 0 pour cercles
483           direct = (i0 == 0)
484           for i, bad in enumerate(edgeEnTrop):
485             if direct:
486               j = i
487             else:
488               j = 1-i
489             if bad:
490               outilPart[j] = geompy.MakeProjection(cercles[j],facesOnside[0])
491           partitionPeauByPipe = geompy.MakePartition(facesAndFond, outilPart, list(), list(), geompy.ShapeType["FACE"], 0, list(), 1)
492
493       name="partitionPeauByPipe%d"%ifil
494       geomPublish(initLog.debug, partitionPeauByPipe, name)
495       [edgesPeauFondIn, _, _] = extractionOrientee(fillingFaceExterne, partitionPeauByPipe, centreFondFiss, "EDGE", 1.e-3)
496       [_, _, facesPeauFondOn] = extractionOrientee(fillingFaceExterne, partitionPeauByPipe, centreFondFiss, "FACE", 1.e-3)
497
498       if len(verticesPipePeau) > 0: # --- au moins une extrémité du pipe sur cette face de peau
499         facesPeauSorted, _, _ = sortFaces(facesPeauFondOn)
500         facePeau = facesPeauSorted[-1] # la plus grande face
501       else:
502         facePeau =geompy.MakePartition(facesPeauFondOn, list(), list(), list(), geompy.ShapeType["FACE"], 0, list(), 1)
503       name="facePeau%d"%ifil
504       geomPublish(initLog.debug, facePeau, name)
505
506       facesPipePeau = [None for i in range(len(edgesFissExtPipe))]
507       endsEdgeFond = [None for i in range(len(edgesFissExtPipe))]
508       edgeRadFacePipePeau = [None for i in range(len(edgesFissExtPipe))]
509
510       edgesListees = list()
511       edgesCircPeau = list()
512       verticesCircPeau = list()
513       if len(verticesPipePeau) > 0: # --- au moins une extrémité du pipe sur cette face de peau
514
515         for face in facesPeauSorted[:-1]: # la ou les faces débouchantes, pas la grande face de peau
516           logging.debug("examen face debouchante circulaire")
517           for i,efep in enumerate(edgesFissExtPipe):
518             dist = geompy.MinDistance(face, efep)
519             logging.debug("  distance face circulaire edge %s", dist)
520             if dist < 1e-3:
521               for ik, edpfi in enumerate(edgesPeauFondIn):
522                 if geompy.MinDistance(face, edpfi) < 1e-3:
523                   ikok = ik
524                   break
525               sharedVertices = geompy.GetSharedShapesMulti([face, edgesPeauFondIn[ikok]], geompy.ShapeType["VERTEX"])
526               nameFace = "facePipePeau%d"%i
527               nameVert = "endEdgeFond%d"%i
528               nameEdge = "edgeRadFacePipePeau%d"%i
529               facesPipePeau[i] = face
530               endsEdgeFond[i] = sharedVertices[0]
531               geomPublish(initLog.debug, face, nameFace)
532               geomPublish(initLog.debug, sharedVertices[0], nameVert)
533               edgesFace = geompy.ExtractShapes(face, geompy.ShapeType["EDGE"], True)
534               for edge in edgesFace:
535                 if geompy.MinDistance(edge, sharedVertices[0]) < 1e-3:
536                   edgeRadFacePipePeau[i] = edge
537                   geomPublish(initLog.debug, edge, nameEdge)
538                   break
539
540         # --- edges circulaires de la face de peau et points de jonction de la face externe de fissure
541         logging.debug("facesPipePeau: %s", facesPipePeau)
542         edgesCircPeau = [None for i in range(len(facesPipePeau))]
543         verticesCircPeau = [None for i in range(len(facesPipePeau))]
544         for i,fcirc in enumerate(facesPipePeau):
545           edges = geompy.GetSharedShapesMulti([facePeau, fcirc], geompy.ShapeType["EDGE"])
546           grpEdgesCirc = geompy.CreateGroup(facePeau, geompy.ShapeType["EDGE"])
547           geompy.UnionList(grpEdgesCirc, edges)
548           edgesCircPeau[i] = grpEdgesCirc
549           name = "edgeCirc%d"%i
550           geomPublishInFather(initLog.debug,facePeau, grpEdgesCirc, name)
551           edgesListees = edgesListees + edges
552           vertices = geompy.GetSharedShapesMulti([facePeau, fcirc], geompy.ShapeType["VERTEX"])
553           grpVertCircPeau = geompy.CreateGroup(facePeau, geompy.ShapeType["VERTEX"])
554           geompy.UnionList(grpVertCircPeau, vertices)
555           verticesCircPeau[i] = grpVertCircPeau
556           name = "pointEdgeCirc%d"%i
557           geomPublishInFather(initLog.debug,facePeau, grpVertCircPeau, name)
558         # --- au moins une extrémité du pipe sur cette face de peau
559
560       # --- edges de bord de la face de peau
561
562       edgesFilling = geompy.ExtractShapes(fillingFaceExterne, geompy.ShapeType["EDGE"], False)
563       edgesBords = list()
564       for i, edge in enumerate(edgesFilling):
565         edgepeau = geompy.GetInPlace(facePeau, edge)
566         name = "edgepeau%d"%i
567         geomPublishInFather(initLog.debug,facePeau,edgepeau, name)
568         logging.debug("edgepeau %s", geompy.ShapeInfo(edgepeau))
569         if geompy.ShapeInfo(edgepeau)['EDGE'] > 1:
570           logging.debug("  EDGES multiples")
571           edgs = geompy.ExtractShapes(edgepeau, geompy.ShapeType["EDGE"], False)
572           edgesBords += edgs
573           edgesListees += edgs
574         else:
575           logging.debug("  EDGE")
576           edgesBords.append(edgepeau)
577           edgesListees.append(edgepeau)
578       groupEdgesBordPeau = geompy.CreateGroup(facePeau, geompy.ShapeType["EDGE"])
579       geompy.UnionList(groupEdgesBordPeau, edgesBords)
580       bordsVifs = None
581       if aretesVivesC is not None:
582         bordsVifs = geompy.GetInPlace(facePeau, aretesVivesC)
583       if bordsVifs is not None:
584         geomPublishInFather(initLog.debug,facePeau, bordsVifs, "bordsVifs")
585         groupEdgesBordPeau = geompy.CutGroups(groupEdgesBordPeau, bordsVifs)
586         grptmp = None
587         if len(aretesVivesCoupees) > 0:
588           grpC = geompy.MakeCompound(aretesVivesCoupees)
589           grptmp = geompy.GetInPlace(facePeau, grpC)
590         if grptmp is not None:
591           grpnew = geompy.CutGroups(bordsVifs, grptmp) # ce qui est nouveau dans bordsVifs
592         else:
593           grpnew = bordsVifs
594         if grpnew is not None:
595           edv = geompy.ExtractShapes(grpnew, geompy.ShapeType["EDGE"], False)
596           aretesVivesCoupees += edv
597       logging.debug("aretesVivesCoupees %s",aretesVivesCoupees)
598       geomPublishInFather(initLog.debug,facePeau, groupEdgesBordPeau , "EdgesBords")
599
600       # ---  edges de la face de peau partagées avec la face de fissure
601
602       edgesPeau = geompy.ExtractShapes(facePeau, geompy.ShapeType["EDGE"], False)
603       edges = substractSubShapes(facePeau, edgesPeau, edgesListees)
604       edgesFissurePeau = list()
605       if len(verticesPipePeau) > 0: # --- au moins une extrémité du pipe sur cette face de peau
606         edgesFissurePeau = [None for i in range(len(verticesCircPeau))] # edges associés aux extrémités du pipe, en premier
607         for edge in edges:
608           for i, grpVert in enumerate(verticesCircPeau):
609             if (geompy.MinDistance(grpVert, edge) < 1.e-3) and (edge not in edgesFissurePeau):
610               edgesFissurePeau[i] = edge
611               name = "edgeFissurePeau%d"%i
612               geomPublishInFather(initLog.debug,facePeau,  edge, name)
613         for edge in edges: # on ajoute après les edges manquantes
614           if edge not in edgesFissurePeau:
615             edgesFissurePeau.append(edge)
616       else:
617         for i, edge in enumerate(edges):
618           edgesFissurePeau.append(edge)
619           name = "edgeFissurePeau%d"%i
620           geomPublishInFather(initLog.debug,facePeau,  edge, name)
621
622       ptEdgeFond[ifil] = endsEdgeFond        # pour chaque face [points edge fond de fissure aux débouchés du pipe]
623       fsPipePeau[ifil] = facesPipePeau       # pour chaque face [faces du pipe débouchantes]
624       edRadFPiPo[ifil] = edgeRadFacePipePeau # pour chaque face [edge radiale des faces du pipe débouchantes ]
625       fsFissuExt[ifil] = facesFissExt        # pour chaque face [faces de fissure externes au pipe]
626       edFisExtPe[ifil] = edgesFissExtPeau    # pour chaque face [edge en peau des faces de fissure externes (pas subshape facePeau)]
627       edFisExtPi[ifil] = edgesFissExtPipe    # pour chaque face [edge commun au pipe des faces de fissure externes]
628       facesPeaux[ifil] = facePeau            # pour chaque face : la face de peau finale a mailler (percee des faces débouchantes)
629       edCircPeau[ifil] = edgesCircPeau       # pour chaque face de peau : [groupe subshapes edges circulaires aux débouchés du pipe]
630       ptCircPeau[ifil] = verticesCircPeau    # pour chaque face de peau : [groupe subshapes points sur edges circulaires aux débouchés du pipe]
631       gpedgeBord[ifil] = groupEdgesBordPeau  # pour chaque face de peau : groupe subshape des edges aux bords liés à la partie saine
632       gpedgeVifs[ifil] = bordsVifs           # pour chaque face de peau : groupe subshape des edges aux bords correspondant à des arêtes vives
633       edFissPeau[ifil] = edgesFissurePeau    # pour chaque face de peau : [subshape edge en peau des faces de fissure externes]
634       ptFisExtPi[ifil] = verticesPipePeau    # pour chaque face de peau : [point commun edFissPeau edCircPeau]
635
636   # -----------------------------------------------------------------------
637   # fin de la boucle sur les faces de filling
638   # -----------------------------------------------------------------------
639
640   for i, avc in enumerate(aretesVivesCoupees):
641     name = "areteViveCoupee%d"%i
642     geomPublish(initLog.debug, avc, name)
643
644   # --- identification des faces et edges de fissure externe pour maillage
645
646   facesFissExt = list()
647   edgesFissExtPeau = list()
648   edgesFissExtPipe = list()
649   for ifil in range(nbFacesFilling): # TODO: éliminer les doublons (comparer tous les vertices triés, avec mesure de distance ?)
650     facesFissExt += fsFissuExt[ifil]
651     edgesFissExtPeau += edFisExtPe[ifil]
652     edgesFissExtPipe += edFisExtPi[ifil]
653   logging.debug("---------------------------- identification faces de fissure externes au pipe :%s ", len(facesFissExt))
654   # regroupement des faces de fissure externes au pipe.
655
656   if len(facesFissExt) > 1:
657     faceFissureExterne = geompy.MakePartition(facesFissExt, list(), list(), list(), geompy.ShapeType["FACE"], 0, list(), 0)
658     edgesPipeFissureExterneC = geompy.GetInPlace(faceFissureExterne, geompy.MakeCompound(edgesPipeFiss))    # edgesFissExtPipe peut ne pas couvrir toute la longueur
659     # edgesPeauFissureExterneC = geompy.GetInPlace(faceFissureExterne, geompy.MakeCompound(edgesFissExtPeau))
660     # il peut manquer des edges de faceFissureExterne en contact avec la peau dans edgesFissExtPeau
661     (_, closedFreeBoundaries, _) = geompy.GetFreeBoundary(faceFissureExterne)
662     edgesBordFFE = list()
663     for bound in closedFreeBoundaries:
664       edgesBordFFE += geompy.ExtractShapes(bound, geompy.ShapeType["EDGE"], False)
665     edgesBordFFEid = [ (ed,geompy.GetSubShapeID(faceFissureExterne, ed)) for ed in edgesBordFFE]
666     logging.debug("edgesBordFFEid %s", edgesBordFFEid)
667     edgesPPE = geompy.ExtractShapes(edgesPipeFissureExterneC, geompy.ShapeType["EDGE"], False)
668     edgesPPEid = [ geompy.GetSubShapeID(faceFissureExterne, ed) for ed in edgesPPE]
669     logging.debug("edgesPPEid %s", edgesPPEid)
670     edgesPFE = [ edid[0] for edid in edgesBordFFEid if edid[1] not in edgesPPEid] # on garde toutes les edges de bord non en contact avec le pipe
671     logging.debug("edgesPFE %s", edgesPFE)
672     edgesPeauFissureExterneC = geompy.MakeCompound(edgesPFE)
673   else:
674     faceFissureExterne = facesFissExt[0]
675     edgesPeauFissureExterneC = geompy.MakeCompound(edgesFissExtPeau)
676     edgesPipeFissureExterneC = geompy.MakeCompound(edgesFissExtPipe)
677   wirePipeFissureExterne = geompy.MakeWire(geompy.ExtractShapes(edgesPipeFissureExterneC, geompy.ShapeType["EDGE"], False))
678   geomPublish(initLog.debug, faceFissureExterne, "faceFissureExterne")
679   geomPublishInFather(initLog.debug,faceFissureExterne, edgesPeauFissureExterneC, "edgesPeauFissureExterne")
680   geomPublishInFather(initLog.debug,faceFissureExterne, edgesPipeFissureExterneC, "edgesPipeFissureExterne")
681
682   logging.debug("---------------------------- Preparation Maillage du Pipe --------------")
683   # -----------------------------------------------------------------------
684   # --- preparation maillage du pipe :
685   #     - détections des points a respecter : jonction des edges/faces constituant
686   #       la face de fissure externe au pipe
687   #     - points sur les edges de fond de fissure et edges pipe/face fissure,
688   #     - vecteurs tangents au fond de fissure (normal au disque maillé)
689
690   # --- option de maillage selon le rayon de courbure du fond de fissure
691   lenEdgeFondExt = 0
692   for edff in edgesFondFiss:
693     lenEdgeFondExt += geompy.BasicProperties(edff)[0]
694
695   disfond = list()
696   for filling in facesDefaut:
697     disfond.append(geompy.MinDistance(centreFondFiss, filling))
698   disfond.sort()
699   rcourb = disfond[0]
700   nbSegQuart = 5 # on veut 5 segments min sur un quart de cercle
701   alpha = math.pi/(4*nbSegQuart)
702   deflexion = rcourb*(1.0 -math.cos(alpha))
703   lgmin = lenSegPipe*0.25
704   lgmax = lenSegPipe*1.5
705   logging.debug("rcourb: %s, lenFond:%s, deflexion: %s, lgmin: %s, lgmax: %s", rcourb, lenEdgeFondExt, deflexion, lgmin, lgmax)
706
707   meshFondExt = smesh.Mesh(wireFondFiss)
708   algo1d = meshFondExt.Segment()
709   hypo1d = algo1d.Adaptive(lgmin, lgmax, deflexion) # a ajuster selon la profondeur de la fissure
710
711   is_done = meshFondExt.Compute()
712   text = "meshFondExt.Compute"
713   if is_done:
714     logging.info(text+" OK")
715   else:
716     text = "Erreur au calcul du maillage.\n" + text
717     logging.info(text)
718     raise Exception(text)
719
720   ptGSdic = {} # dictionnaire [paramètre sur la courbe] --> point géométrique
721   allNodeIds = meshFondExt.GetNodesId()
722   for nodeId in allNodeIds:
723     xyz = meshFondExt.GetNodeXYZ(nodeId)
724     #logging.debug("nodeId %s, coords %s", nodeId, str(xyz))
725     pt = geompy.MakeVertex(xyz[0], xyz[1], xyz[2])
726     parametre, PointOnEdge, EdgeInWireIndex = geompy.MakeProjectionOnWire(pt, wireFondFiss) # parametre compris entre 0 et 1
727     edgeOrder = edgesIdByOrientation[EdgeInWireIndex]
728     ptGSdic[(edgeOrder, EdgeInWireIndex, parametre)] = pt
729     #logging.debug("nodeId %s, parametre %s", nodeId, str(parametre))
730   usort = sorted(ptGSdic)
731   logging.debug("nombre de points obtenus par deflexion %s",len(usort))
732
733   centres = list()
734   origins = list()
735   normals = list()
736   for edu in usort:
737     vertcx = ptGSdic[edu]
738     norm = geompy.MakeTangentOnCurve(edgesFondFiss[edu[1]], edu[2])
739     plan = geompy.MakePlane(vertcx, norm, 3*rayonPipe)
740     part = geompy.MakePartition([plan], [wirePipeFiss], list(), list(), geompy.ShapeType["VERTEX"], 0, list(), 0)
741     liste = geompy.ExtractShapes(part, geompy.ShapeType["VERTEX"], True)
742     if len(liste) == 5: # 4 coins du plan plus intersection recherchée
743       for point in liste:
744         if geompy.MinDistance(point, vertcx) < 1.1*rayonPipe: # les quatre coins sont plus loin
745           vertpx = point
746           break
747       centres.append(vertcx)
748       origins.append(vertpx)
749       normals.append(norm)
750 #      name = "vertcx%d"%i
751 #      geomPublishInFather(initLog.debug,wireFondFiss, vertcx, name)
752 #      name = "vertpx%d"%i
753 #      geomPublishInFather(initLog.debug,wireFondFiss, vertpx, name)
754 #      name = "plan%d"%i
755 #      geomPublishInFather(initLog.debug,wireFondFiss, plan, name)
756
757   # --- maillage du pipe étendu, sans tenir compte de l'intersection avec la face de peau
758
759   logging.debug("nbsegCercle %s", nbsegCercle)
760
761   # -----------------------------------------------------------------------
762   # --- points géométriques
763
764   gptsdisks = list() # vertices géométrie de tous les disques
765   raydisks = [list() for i in range(nbsegCercle)]
766   for i, centres_i in enumerate(centres): # boucle sur les disques
767     gptdsk = list() # vertices géométrie d'un disque
768     vertcx = centres_i
769     vertpx = origins[i]
770     normal = normals[i]
771     vec1 = geompy.MakeVector(vertcx, vertpx)
772
773     points = [vertcx] # les points du rayon de référence
774     for j in range(nbsegRad):
775       pt = geompy.MakeTranslationVectorDistance(vertcx, vec1, (j+1)*float(rayonPipe)/nbsegRad)
776       points.append(pt)
777     gptdsk.append(points)
778     pt = geompy.MakeTranslationVectorDistance(vertcx, vec1, 1.5*rayonPipe)
779     rayon = geompy.MakeLineTwoPnt(vertcx, pt)
780     raydisks[0].append(rayon)
781
782     for k in range(nbsegCercle-1):
783       angle = (k+1)*2*math.pi/nbsegCercle
784       pts = [vertcx] # les points d'un rayon obtenu par rotation
785       for j in range(nbsegRad):
786         pt = geompy.MakeRotation(points[j+1], normal, angle)
787         pts.append(pt)
788       gptdsk.append(pts)
789       ray = geompy.MakeRotation(rayon, normal, angle)
790       raydisks[k+1].append(ray)
791
792     gptsdisks.append(gptdsk)
793
794   # -----------------------------------------------------------------------
795   # --- recherche des points en trop (externes au volume à remailler)
796   #     - on associe chaque extrémité du pipe à une face filling
797   #     - on part des disques aux extrémités du pipe
798   #     - pour chaque disque, on prend les vertices de géométrie,
799   #       on marque leur position relative à la face.
800   #     - on s'arrete quand tous les noeuds sont dedans
801
802   logging.debug("---------------------------- recherche des points du pipe a éliminer --------------")
803
804   pt0 = centres[0]
805   pt1 = centres[-1]
806   idFillingFromBout = [None, None]                 # contiendra l'index du filling pour les extrémités 0 et 1
807   for ifil in range(nbFacesFilling):
808     for _, pt in enumerate(ptEdgeFond[ifil]): # il y a un ou deux points débouchant sur cette face
809       if geompy.MinDistance(pt,pt0) < geompy.MinDistance(pt,pt1): # TODO: trouver plus fiable pour les cas tordus...
810         idFillingFromBout[0] = ifil
811       else:
812         idFillingFromBout[1] = ifil
813   logging.debug("association bouts du pipe - faces de filling: %s", idFillingFromBout)
814
815   facesPipePeau = list()
816   edgeRadFacePipePeau = list()
817   for ifil in range(nbFacesFilling):
818     facesPipePeau += fsPipePeau[ifil]
819     edgeRadFacePipePeau += edRadFPiPo[ifil]
820
821   logging.debug("recherche des disques de noeuds complètement internes")
822   idisklim = list() # indices des premier et dernier disques internes
823   idiskout = list() # indices des premier et dernier disques externes
824   for bout in range(2):
825     if bout == 0:
826       idisk = -1
827       inc = 1
828       numout = -1
829     else:
830       idisk = len(gptsdisks)
831       inc = -1
832       numout = len(gptsdisks)
833     inside = False
834     outside = True
835     while not inside:
836       idisk = idisk + inc
837       logging.debug("examen disque %s", idisk)
838       gptdsk = gptsdisks[idisk]
839       inside = True
840       for k in range(nbsegCercle):
841         points = gptdsk[k]
842         for j, pt in enumerate(points):
843           side = whichSideVertex(facesDefaut[idFillingFromBout[bout]], pt)
844           if side < 0:
845             if outside: # premier point detecté dedans
846               outside = False
847               numout = idisk -inc # le disque précédent était dehors
848           else:
849             inside = False # ce point est dehors
850         if not inside and not outside:
851           break
852     idisklim.append(idisk)  # premier et dernier disques internes
853     idiskout.append(numout) # premier et dernier disques externes
854
855   # --- listes de nappes radiales en filling à chaque extrémité débouchante
856   facesDebouchantes = [False, False]
857   idFacesDebouchantes = [-1, -1] # contiendra les indices des faces disque débouchantes (facesPipePeau)
858   listNappes =list()
859   for i, idisk in enumerate(idisklim):
860     numout = idiskout[i]
861     logging.debug("extremité %s, indices disques interne %s, externe %s",i, idisk, numout)
862     nappes = list()
863     if  (idisk != 0) and (idisk != len(gptsdisks)-1): # si extrémité débouchante
864       for k in range(nbsegCercle):
865         if i == 0:
866           iddeb = max(0, numout)
867           idfin = max(iddeb+3,idisk+1) # il faut 3 rayons pour faire un filling qui suive le fond de fissure
868           #logging.debug("extremité %s, indices retenus interne %s, externe %s",i, idfin, iddeb)
869           comp = geompy.MakeCompound(raydisks[k][iddeb:idfin])
870           name='compoundRay%d'%k
871           geomPublish(initLog.debug, comp, name)
872         else:
873           idfin = min(len(gptsdisks), numout+1)
874           iddeb = min(idfin-3, idisk) # il faut 3 rayons pour faire un filling qui suive le fond de fissure
875           #logging.debug("extremité %s, indices retenus interne %s, externe %s",i, idfin, iddeb)
876           comp = geompy.MakeCompound(raydisks[k][iddeb:idfin])
877           name='compoundRay%d'%k
878           geomPublish(initLog.debug, comp, name)
879         nappe = geompy.MakeFilling(comp, 2, 5, 0.0001, 0.0001, 0, GEOM.FOM_Default)
880         nappes.append(nappe)
881         name='nappe%d'%k
882         geomPublish(initLog.debug, nappe, name)
883         facesDebouchantes[i] = True
884     listNappes.append(nappes)
885
886   # --- mise en correspondance avec les indices des faces disque débouchantes (facesPipePeau)
887   for i, nappes in enumerate(listNappes):
888     if facesDebouchantes[i]:
889       for k, face in enumerate(facesPipePeau):
890         edge = geompy.MakeSection(face, nappes[0])
891         if geompy.NbShapes(edge, geompy.ShapeType["EDGE"]) > 0:
892           idFacesDebouchantes[i] = k
893           break
894   logging.debug("idFacesDebouchantes: %s", idFacesDebouchantes)
895
896   # --- construction des listes d'edges radiales sur chaque extrémité débouchante
897   listEdges = list()
898   for i, nappes in enumerate(listNappes):
899     indice = idFacesDebouchantes[i] # indice de face débouchante (facesPipePeau)
900     if indice < 0:
901       listEdges.append(list())
902     else:
903       face = facesPipePeau[indice]
904       edges = [edgeRadFacePipePeau[indice]]
905       for k, nappe in enumerate(nappes):
906         if k > 0:
907           obj = geompy.MakeSection(face, nappes[k]) # normalement une edge, parfois un compound d'edges dont un tout petit
908           edge = obj
909           vs = geompy.ExtractShapes(obj, geompy.ShapeType["VERTEX"], False)
910           if len(vs) > 2:
911             eds = geompy.ExtractShapes(obj, geompy.ShapeType["EDGE"], False)
912             [edsorted, _,maxl] = sortEdges(eds)
913             edge = edsorted[-1]
914           else:
915             maxl = geompy.BasicProperties(edge)[0]
916           if maxl < 0.01: # problème MakeSection
917             logging.debug("problème MakeSection recherche edge radiale %s, longueur trop faible: %s, utilisation partition", k, maxl)
918             partNappeFace = geompy.MakePartition([face, nappes[k]], list() , list(), list(), geompy.ShapeType["FACE"], 0, list(), 0)
919             edps= geompy.ExtractShapes(partNappeFace, geompy.ShapeType["EDGE"], False)
920             ednouv = list()
921             for ii, ed in enumerate(edps):
922               vxs = geompy.ExtractShapes(ed, geompy.ShapeType["VERTEX"], False)
923               distx = [geompy.MinDistance(vx, face) for vx in vxs]
924               distx += [geompy.MinDistance(vx, nappes[k]) for vx in vxs]
925               dmax = max(distx)
926               logging.debug("  dmax %s",dmax)
927               if dmax < 0.01:
928                 ednouv.append(ed)
929             logging.debug("  edges issues de la partition: %s", ednouv)
930             for ii, ed in enumerate(ednouv):
931               geomPublish(initLog.debug, ed, "ednouv%d"%ii)
932             [edsorted, _,maxl] = sortEdges(ednouv)
933             logging.debug("  longueur edge trouvée: %s", maxl)
934             edge = edsorted[-1]
935           edges.append(edge)
936           name = 'edgeEndPipe%d'%k
937           geomPublish(initLog.debug, edge, name)
938       listEdges.append(edges)
939
940   # --- création des points du maillage du pipe sur la face de peau
941   for i, edges in enumerate(listEdges):
942     indice = idFacesDebouchantes[i] # indice de face débouchante (facesPipePeau)
943     if indice >= 0:
944       gptdsk = list()
945       if indice > 0: # indice vaut 0 ou 1
946         indice = -1  # si indice vaut 1, on prend le dernier élément de la liste (1 ou 2 extrémités débouchent sur la face)
947       centre = ptEdgeFond[idFillingFromBout[i]][indice]
948       name = "centre%d"%indice
949       geomPublish(initLog.debug, centre, name)
950       vertPipePeau = ptFisExtPi[idFillingFromBout[i]][indice]
951       geomPublishInFather(initLog.debug,centre, vertPipePeau, "vertPipePeau")
952       grpsEdgesCirc = edCircPeau[idFillingFromBout[i]] # liste de groupes
953       edgesCirc = list()
954       for grpEdgesCirc in grpsEdgesCirc:
955         edgesCirc += geompy.ExtractShapes(grpEdgesCirc, geompy.ShapeType["EDGE"], False)
956       for k, edge in enumerate(edges):
957         extrems = geompy.ExtractShapes(edge, geompy.ShapeType["VERTEX"], True)
958         if geompy.MinDistance(centre, extrems[0]) < geompy.MinDistance(centre, extrems[1]):
959           bout = extrems[1]
960         else:
961           bout = extrems[0]
962         # ajustement du point extrémité (bout) sur l'edge circulaire en face de peau
963         logging.debug("edgesCirc: %s", edgesCirc)
964         distEdgeCirc = [(geompy.MinDistance(bout, edgeCirc), k2, edgeCirc) for k2, edgeCirc in enumerate(edgesCirc)]
965         distEdgeCirc.sort()
966         logging.debug("distEdgeCirc: %s", distEdgeCirc)
967         u = projettePointSurCourbe(bout, distEdgeCirc[0][2])
968         if (abs(u) < 0.02) or (abs(1-u) < 0.02): # les points très proches d'une extrémité doivent y être mis précisément.
969           extrCircs = geompy.ExtractShapes(distEdgeCirc[0][2], geompy.ShapeType["VERTEX"], True)
970           if geompy.MinDistance(bout, extrCircs[0]) < geompy.MinDistance(bout, extrCircs[1]):
971             bout = extrCircs[0]
972           else:
973             bout = extrCircs[1]
974         else:
975           bout = geompy.MakeVertexOnCurve(distEdgeCirc[0][2], u)
976         name ="bout%d"%k
977         geomPublishInFather(initLog.debug,centre, bout, name)
978         # enregistrement des points dans la structure
979         points = list()
980         for j in range(nbsegRad +1):
981           u = j/float(nbsegRad)
982           points.append(geompy.MakeVertexOnCurve(edge, u))
983         if geompy.MinDistance(bout, points[0]) < geompy.MinDistance(centre, points[0]):
984           points.reverse()
985         points[0] = centre
986         points[-1] = bout
987         gptdsk.append(points)
988       if i == 0:
989         gptsdisks[idisklim[0] -1] = gptdsk
990         idisklim[0] = idisklim[0] -1
991       else:
992         gptsdisks[idisklim[1] +1] = gptdsk
993         idisklim[1] = idisklim[1] +1
994
995   # --- ajustement precis des points sur edgesPipeFissureExterneC
996
997   edgesPFE = geompy.ExtractShapes(edgesPipeFissureExterneC, geompy.ShapeType["EDGE"], False)
998   verticesPFE, _ = findWireIntermediateVertices(wirePipeFissureExterne)  # vertices intermédiaires (des points en trop dans ptsInWireFissExtPipe)
999   idiskmin = idisklim[0] + 1 # on ne prend pas le disque sur la peau, déjà ajusté
1000   idiskmax = idisklim[1]     # on ne prend pas le disque sur la peau, déjà ajusté
1001   idiskint = list()
1002   for vtx in verticesPFE:
1003     distPtVt = list()
1004     for idisk in range(idiskmin, idiskmax):
1005       gptdsk = gptsdisks[idisk]
1006       pt = gptdsk[0][-1]       # le point sur l'edge de la fissure externe au pipe
1007       distPtVt.append((geompy.MinDistance(pt, vtx), idisk))
1008     distPtVt.sort()
1009     idiskint.append(distPtVt[0][1])
1010     gptsdisks[idiskint[-1]][0][-1] = vtx
1011     logging.debug("ajustement point sur edgePipeFissureExterne, vertex: %s %s", idiskint[-1], distPtVt[0][0])
1012   for idisk in range(idiskmin, idiskmax):
1013     if idisk in idiskint:
1014       break
1015     logging.debug("ajustement point sur edgePipeFissureExterne: %s", idisk)
1016     gptdsk = gptsdisks[idisk]
1017     pt = gptdsk[0][-1]       # le point sur l'edge de la fissure externe au pipe
1018     distPtEd = [(geompy.MinDistance(pt, edgePFE), k, edgePFE) for k, edgePFE in enumerate(edgesPFE)]
1019     distPtEd.sort()
1020     edgePFE = distPtEd[0][2]
1021     u = projettePointSurCourbe(pt, edgePFE)
1022     ptproj = geompy.MakeVertexOnCurve(edgePFE, u)
1023     gptsdisks[idisk][0][-1] = ptproj
1024
1025   # -----------------------------------------------------------------------
1026   # --- maillage effectif du pipe
1027
1028   logging.debug("---------------------------- maillage effectif du pipe --------------")
1029   meshPipe = smesh.Mesh(None, "meshPipe")
1030   fondFissGroup = meshPipe.CreateEmptyGroup(SMESH.EDGE, "FONDFISS")
1031   nodesFondFissGroup = meshPipe.CreateEmptyGroup(SMESH.NODE, "nfondfis")
1032   faceFissGroup = meshPipe.CreateEmptyGroup(SMESH.FACE, "fisInPi")
1033   edgeFaceFissGroup = meshPipe.CreateEmptyGroup(SMESH.EDGE, "edgeFaceFiss")
1034   edgeCircPipe0Group = meshPipe.CreateEmptyGroup(SMESH.EDGE, "edgeCircPipe0")
1035   edgeCircPipe1Group = meshPipe.CreateEmptyGroup(SMESH.EDGE, "edgeCircPipe1")
1036   faceCircPipe0Group = meshPipe.CreateEmptyGroup(SMESH.FACE, "faceCircPipe0")
1037   faceCircPipe1Group = meshPipe.CreateEmptyGroup(SMESH.FACE, "faceCircPipe1")
1038   mptsdisks = list()  # vertices maillage de tous les disques
1039   mEdges = list()     # identifiants edges maillage fond de fissure
1040   mEdgeFaces = list() # identifiants edges maillage edge face de fissure externe
1041   mFaces = list()     # identifiants faces maillage fissure
1042   mVols  = list()     # identifiants volumes maillage pipe
1043
1044   mptdsk = list()
1045   for idisk in range(idisklim[0], idisklim[1]+1): # boucle sur les disques internes
1046
1047     # -----------------------------------------------------------------------
1048     # --- points
1049
1050     gptdsk = gptsdisks[idisk]
1051     if idisk > idisklim[0]:
1052       oldmpts = mptdsk
1053     mptdsk = list() # vertices maillage d'un disque
1054     for k in range(nbsegCercle):
1055       points = gptdsk[k]
1056       mptids = list()
1057       for j, pt in enumerate(points):
1058         if j == 0 and k > 0:
1059           indice = mptdsk[0][0]
1060         else:
1061           coords = geompy.PointCoordinates(pt)
1062           indice = meshPipe.AddNode(coords[0], coords[1], coords[2])
1063         mptids.append(indice)
1064       mptdsk.append(mptids)
1065     mptsdisks.append(mptdsk)
1066
1067     # -----------------------------------------------------------------------
1068     # --- groupes edges cercles debouchants
1069
1070     if idisk == idisklim[0]:
1071       pts = list()
1072       for k in range(nbsegCercle):
1073         pts.append(mptdsk[k][-1])
1074       edges = list()
1075       for k, pts_k in enumerate(pts):
1076         k1 = (k+1)%len(pts)
1077         idEdge = meshPipe.AddEdge([pts_k, pts[k1]])
1078         edges.append(idEdge)
1079       edgeCircPipe0Group.Add(edges)
1080
1081     if idisk == idisklim[1]:
1082       pts = list()
1083       for k in range(nbsegCercle):
1084         pts.append(mptdsk[k][-1])
1085       edges = list()
1086       for k, pts_k in enumerate(pts):
1087         k1 = (k+1)%len(pts)
1088         idEdge = meshPipe.AddEdge([pts_k, pts[k1]])
1089         edges.append(idEdge)
1090       edgeCircPipe1Group.Add(edges)
1091
1092     # -----------------------------------------------------------------------
1093     # --- groupes faces  debouchantes
1094
1095     if idisk == idisklim[0]:
1096       faces = list()
1097       for j in range(nbsegRad):
1098         for k in range(nbsegCercle):
1099           k1 = k+1
1100           if k ==  nbsegCercle-1:
1101             k1 = 0
1102           if j == 0:
1103             idf = meshPipe.AddFace([mptdsk[k][0], mptdsk[k][1], mptdsk[k1][1]]) # triangle
1104           else:
1105             idf = meshPipe.AddFace([mptdsk[k][j], mptdsk[k][j+1], mptdsk[k1][j+1], mptdsk[k1][j]]) # quadrangle
1106           faces.append(idf)
1107       faceCircPipe0Group.Add(faces)
1108
1109     if idisk == idisklim[1]:
1110       faces = list()
1111       for j in range(nbsegRad):
1112         for k in range(nbsegCercle):
1113           k1 = k+1
1114           if k ==  nbsegCercle-1:
1115             k1 = 0
1116           if j == 0:
1117             idf = meshPipe.AddFace([mptdsk[k][0], mptdsk[k][1], mptdsk[k1][1]]) # triangle
1118           else:
1119             idf = meshPipe.AddFace([mptdsk[k][j], mptdsk[k][j+1], mptdsk[k1][j+1], mptdsk[k1][j]]) # quadrangle
1120           faces.append(idf)
1121       faceCircPipe1Group.Add(faces)
1122
1123     # -----------------------------------------------------------------------
1124     # --- mailles volumiques, groupes noeuds et edges de fond de fissure, groupe de face de fissure
1125
1126     if idisk == idisklim[0]:
1127       mEdges.append(0)
1128       mEdgeFaces.append(0)
1129       mFaces.append([0])
1130       mVols.append([[0]])
1131       nodesFondFissGroup.Add([mptdsk[0][0]])
1132     else:
1133       ide = meshPipe.AddEdge([oldmpts[0][0], mptdsk[0][0]])
1134       mEdges.append(ide)
1135       fondFissGroup.Add([ide])
1136       nodesFondFissGroup.Add([mptdsk[0][0]])
1137       ide2 = meshPipe.AddEdge([oldmpts[0][-1], mptdsk[0][-1]])
1138       mEdgeFaces.append(ide2)
1139       edgeFaceFissGroup.Add([ide2])
1140       idFaces = list()
1141       idVols = list()
1142
1143       for j in range(nbsegRad):
1144         idf = meshPipe.AddFace([oldmpts[0][j], mptdsk[0][j], mptdsk[0][j+1], oldmpts[0][j+1]])
1145         faceFissGroup.Add([idf])
1146         idFaces.append(idf)
1147
1148         idVolCercle = list()
1149         for k in range(nbsegCercle):
1150           k1 = k+1
1151           if k ==  nbsegCercle-1:
1152             k1 = 0
1153           if j == 0:
1154             idv = meshPipe.AddVolume([mptdsk[k][j], mptdsk[k][j+1], mptdsk[k1][j+1],
1155                                       oldmpts[k][j], oldmpts[k][j+1], oldmpts[k1][j+1]])
1156           else:
1157             idv = meshPipe.AddVolume([mptdsk[k][j], mptdsk[k][j+1], mptdsk[k1][j+1], mptdsk[k1][j],
1158                                       oldmpts[k][j], oldmpts[k][j+1], oldmpts[k1][j+1], oldmpts[k1][j]])
1159           idVolCercle.append(idv)
1160         idVols.append(idVolCercle)
1161
1162       mFaces.append(idFaces)
1163       mVols.append(idVols)
1164
1165   pipeFissGroup = meshPipe.CreateEmptyGroup( SMESH.VOLUME, 'PIPEFISS' )
1166   _ = pipeFissGroup.AddFrom( meshPipe.GetMesh() )
1167
1168   _, _, _ = meshPipe.MakeBoundaryElements(SMESH.BND_2DFROM3D, "pipeBoundaries")
1169   edgesCircPipeGroup = [edgeCircPipe0Group, edgeCircPipe1Group]
1170
1171   # --- fin du maillage du pipe
1172   # -----------------------------------------------------------------------
1173   # --- edges de bord, faces défaut à respecter
1174
1175   _ = smesh.CreateFilterManager()
1176   _, internalBoundary, _NoneGroup = internalBoundary.MakeBoundaryElements( SMESH.BND_1DFROM2D, '', '', 0, [  ])
1177   criteres = list()
1178   unCritere = smesh.GetCriterion(SMESH.EDGE,SMESH.FT_FreeBorders,SMESH.FT_Undefined,0)
1179   criteres.append(unCritere)
1180   filtre = smesh.GetFilterFromCriteria(criteres)
1181   bordsLibres = internalBoundary.MakeGroupByFilter( 'bords', filtre )
1182   smesh.SetName(bordsLibres, 'bordsLibres')
1183
1184   # --- pour aider l'algo hexa-tetra à ne pas mettre de pyramides à l'exterieur des volumes repliés sur eux-mêmes
1185   #     on désigne les faces de peau en quadrangles par le groupe "skinFaces"
1186
1187   skinFaces = internalBoundary.CreateEmptyGroup( SMESH.FACE, 'skinFaces' )
1188   _ = skinFaces.AddFrom( internalBoundary.GetMesh() )
1189
1190   # --- maillage des éventuelles arêtes vives entre faces reconstruites
1191
1192   if aretesVivesCoupees:
1193
1194     aretesVivesC = geompy.MakeCompound(aretesVivesCoupees)
1195     meshAretesVives = smesh.Mesh(aretesVivesC)
1196     algo1d = meshAretesVives.Segment()
1197     hypo1d = algo1d.LocalLength(dmoyen,list(),1e-07)
1198     putName(algo1d.GetSubMesh(), "aretesVives")
1199     putName(algo1d, "algo1d_aretesVives")
1200     putName(hypo1d, "hypo1d_aretesVives")
1201
1202     is_done = meshAretesVives.Compute()
1203     text = "meshAretesVives.Compute"
1204     if is_done:
1205       logging.info(text+" OK")
1206     else:
1207       text = "Erreur au calcul du maillage.\n" + text
1208       logging.info(text)
1209       raise Exception(text)
1210     grpAretesVives = meshAretesVives.CreateEmptyGroup( SMESH.EDGE, 'grpAretesVives' )
1211     _ = grpAretesVives.AddFrom( meshAretesVives.GetMesh() )
1212
1213   # -----------------------------------------------------------------------
1214   # --- maillage faces de fissure
1215
1216   logging.debug("---------------------------- maillage faces de fissure externes au pipe :%s --------------", len(facesFissExt))
1217
1218   meshFaceFiss = smesh.Mesh(faceFissureExterne)
1219   logging.info("Maillage avec %s", mailleur)
1220   if ( mailleur == "MeshGems"):
1221     algo2d = meshFaceFiss.Triangle(algo=smeshBuilder.MG_CADSurf)
1222     hypo2d = algo2d.Parameters()
1223     hypo2d.SetPhySize( areteFaceFissure )
1224     hypo2d.SetMinSize( rayonPipe/float(nbsegRad) )
1225     hypo2d.SetMaxSize( areteFaceFissure*3. )
1226     hypo2d.SetChordalError( areteFaceFissure*0.25 )
1227     hypo2d.SetVerbosity( 0 )
1228   else:
1229     algo2d = meshFaceFiss.Triangle(algo=smeshBuilder.NETGEN_1D2D)
1230     hypo2d = algo2d.Parameters()
1231     hypo2d.SetMaxSize( areteFaceFissure )
1232     hypo2d.SetSecondOrder( 0 )
1233     hypo2d.SetOptimize( 1 )
1234     hypo2d.SetFineness( 2 )
1235     hypo2d.SetMinSize( rayonPipe/float(nbsegRad) )
1236     hypo2d.SetQuadAllowed( 0 )
1237   putName(algo2d.GetSubMesh(), "faceFiss")
1238   putName(algo2d, "algo2d_faceFiss")
1239   putName(hypo2d, "hypo2d_faceFiss")
1240
1241   algo1d = meshFaceFiss.UseExisting1DElements(geom=edgesPipeFissureExterneC)
1242   hypo1d = algo1d.SourceEdges([ edgeFaceFissGroup ],0,0)
1243   putName(algo1d.GetSubMesh(), "edgeFissPeau")
1244   putName(algo1d, "algo1d_edgeFissPeau")
1245   putName(hypo1d, "hypo1d_edgeFissPeau")
1246
1247   _ = meshFaceFiss.GroupOnGeom(faceFissureExterne, "fisOutPi", SMESH.FACE)
1248   grpEdgesPeauFissureExterne = meshFaceFiss.GroupOnGeom(edgesPeauFissureExterneC,'edgesPeauFissureExterne',SMESH.EDGE)
1249   _ = meshFaceFiss.GroupOnGeom(edgesPipeFissureExterneC,'edgesPipeFissureExterne',SMESH.EDGE)
1250
1251   is_done = meshFaceFiss.Compute()
1252   text = "meshFaceFiss.Compute"
1253   if is_done:
1254     logging.info(text+" OK")
1255   else:
1256     text = "Erreur au calcul du maillage.\n" + text
1257     logging.info(text)
1258     raise Exception(text)
1259
1260   # --- maillage faces de peau
1261
1262   boutFromIfil = [None for i in range(nbFacesFilling)]
1263   if idFillingFromBout[0] != idFillingFromBout[1]: # repérage des extremites du pipe quand elles débouchent sur des faces différentes
1264     boutFromIfil[idFillingFromBout[0]] = 0
1265     boutFromIfil[idFillingFromBout[1]] = 1
1266
1267   logging.debug("---------------------------- maillage faces de peau --------------")
1268   meshesFacesPeau = list()
1269   for ifil in range(nbFacesFilling):
1270     meshFacePeau = None
1271     if partitionsPeauFissFond[ifil] is None: # face de peau maillage sain intacte
1272
1273       # --- edges de bord de la face de filling
1274       filling = facesDefaut[ifil]
1275       edgesFilling = geompy.ExtractShapes(filling, geompy.ShapeType["EDGE"], False)
1276       groupEdgesBordPeau = geompy.CreateGroup(filling, geompy.ShapeType["EDGE"])
1277       geompy.UnionList(groupEdgesBordPeau, edgesFilling)
1278       geomPublishInFather(initLog.debug,filling, groupEdgesBordPeau , "EdgesBords")
1279
1280       meshFacePeau = smesh.Mesh(facesDefaut[ifil])
1281
1282       algo1d = meshFacePeau.UseExisting1DElements(geom=groupEdgesBordPeau)
1283       hypo1d = algo1d.SourceEdges([ bordsLibres ],0,0)
1284       putName(algo1d.GetSubMesh(), "bordsLibres", ifil)
1285       putName(algo1d, "algo1d_bordsLibres", ifil)
1286       putName(hypo1d, "hypo1d_bordsLibres", ifil)
1287
1288     else:
1289
1290       facePeau           = facesPeaux[ifil] # pour chaque face : la face de peau finale a mailler (percée des faces débouchantes)
1291       edgesCircPeau      = edCircPeau[ifil] # pour chaque face de peau : [subshape edge circulaire aux débouchés du pipe]
1292       verticesCircPeau   = ptCircPeau[ifil] # pour chaque face de peau : [subshape point sur edge circulaire aux débouchés du pipe]
1293       groupEdgesBordPeau = gpedgeBord[ifil] # pour chaque face de peau : groupe subshape des edges aux bords liés à la partie saine
1294       bordsVifs          = gpedgeVifs[ifil] # pour chaque face de peau : groupe subshape des edges aux bords correspondant à des arêtes vives
1295       edgesFissurePeau   = edFissPeau[ifil] # pour chaque face de peau : [subshape edge en peau des faces de fissure externes]
1296
1297       meshFacePeau = smesh.Mesh(facePeau)
1298
1299       algo1d = meshFacePeau.UseExisting1DElements(geom=groupEdgesBordPeau)
1300       hypo1d = algo1d.SourceEdges([ bordsLibres ],0,0)
1301       putName(algo1d.GetSubMesh(), "bordsLibres", ifil)
1302       putName(algo1d, "algo1d_bordsLibres", ifil)
1303       putName(hypo1d, "hypo1d_bordsLibres", ifil)
1304
1305       algo1d = meshFacePeau.UseExisting1DElements(geom=geompy.MakeCompound(edgesFissurePeau))
1306       hypo1d = algo1d.SourceEdges([ grpEdgesPeauFissureExterne ],0,0)
1307       putName(algo1d.GetSubMesh(), "edgePeauFiss", ifil)
1308       putName(algo1d, "algo1d_edgePeauFiss", ifil)
1309       putName(hypo1d, "hypo1d_edgePeauFiss", ifil)
1310
1311       if bordsVifs is not None:
1312         algo1d = meshFacePeau.UseExisting1DElements(geom=bordsVifs)
1313         hypo1d = algo1d.SourceEdges([ grpAretesVives ],0,0)
1314         putName(algo1d.GetSubMesh(), "bordsVifs", ifil)
1315         putName(algo1d, "algo1d_bordsVifs", ifil)
1316         putName(hypo1d, "hypo1d_bordsVifs", ifil)
1317
1318       for i, edgeCirc in enumerate(edgesCircPeau):
1319         if edgeCirc is not None:
1320           algo1d = meshFacePeau.UseExisting1DElements(geom=edgeCirc)
1321           if boutFromIfil[ifil] is None:
1322             hypo1d = algo1d.SourceEdges([ edgesCircPipeGroup[i] ],0,0)
1323           else:
1324             hypo1d = algo1d.SourceEdges([ edgesCircPipeGroup[boutFromIfil[ifil]] ],0,0)
1325           name = "cercle%d"%i
1326           putName(algo1d.GetSubMesh(), name, ifil)
1327           putName(algo1d, "algo1d_" + name, ifil)
1328           putName(hypo1d, "hypo1d_" + name, ifil)
1329
1330     logging.info("Maillage avec %s", mailleur)
1331     if ( mailleur == "MeshGems"):
1332       algo2d = meshFacePeau.Triangle(algo=smeshBuilder.MG_CADSurf)
1333       hypo2d = algo2d.Parameters()
1334       hypo2d.SetPhySize( dmoyen )
1335       hypo2d.SetMinSize( rayonPipe/float(nbsegRad) )
1336       hypo2d.SetMaxSize( dmoyen*3. )
1337       hypo2d.SetChordalError( dmoyen*0.25 )
1338       hypo2d.SetVerbosity( 0 )
1339     else:
1340       algo2d = meshFacePeau.Triangle(algo=smeshBuilder.NETGEN_1D2D)
1341       hypo2d = algo2d.Parameters()
1342       hypo2d.SetMaxSize( dmoyen*0.75 )
1343       hypo2d.SetOptimize( 1 )
1344       hypo2d.SetFineness( 2 )
1345       hypo2d.SetMinSize( rayonPipe/float(nbsegRad) )
1346       hypo2d.SetQuadAllowed( 0 )
1347     putName(algo2d.GetSubMesh(), "facePeau", ifil)
1348     putName(algo2d, "algo2d_facePeau", ifil)
1349     putName(hypo2d, "hypo2d_facePeau", ifil)
1350
1351     is_done = meshFacePeau.Compute()
1352     text = "meshFacePeau {} Compute".format(ifil)
1353     if is_done:
1354       logging.info(text+" OK")
1355     else:
1356       text = "Erreur au calcul du maillage.\n" + text
1357       logging.info(text)
1358       raise Exception(text)
1359
1360     GroupFaces = meshFacePeau.CreateEmptyGroup( SMESH.FACE, "facePeau%d"%ifil )
1361     _ = GroupFaces.AddFrom( meshFacePeau.GetMesh() )
1362     meshesFacesPeau.append(meshFacePeau)
1363
1364   # --- regroupement des maillages du défaut
1365
1366   listMeshes = [internalBoundary.GetMesh(),
1367                 meshPipe.GetMesh(),
1368                 meshFaceFiss.GetMesh()]
1369   for mp in meshesFacesPeau:
1370     listMeshes.append(mp.GetMesh())
1371
1372   meshBoiteDefaut = smesh.Concatenate(listMeshes, 1, 1, 1e-05,False)
1373   # pour aider l'algo hexa-tetra à ne pas mettre de pyramides à l'exterieur des volumes repliés sur eux-mêmes
1374   # on désigne les faces de peau en quadrangles par le groupe "skinFaces"
1375   group_faceFissOutPipe = None
1376   group_faceFissInPipe = None
1377   groups = meshBoiteDefaut.GetGroups()
1378   for grp in groups:
1379     if grp.GetType() == SMESH.FACE:
1380       #if "internalBoundary" in grp.GetName():
1381       #  grp.SetName("skinFaces")
1382       if grp.GetName() == "fisOutPi":
1383         group_faceFissOutPipe = grp
1384       elif grp.GetName() == "fisInPi":
1385         group_faceFissInPipe = grp
1386
1387   # le maillage NETGEN ne passe pas toujours ==> on force l'usage de MG_Tetra
1388   mailleur = "MeshGems"
1389   logging.info("Maillage avec %s", mailleur)
1390   if ( mailleur == "MeshGems"):
1391     algo3d = meshBoiteDefaut.Tetrahedron(algo=smeshBuilder.MG_Tetra)
1392   else:
1393     algo3d = meshBoiteDefaut.Tetrahedron(algo=smeshBuilder.NETGEN)
1394     hypo3d = algo3d.MaxElementVolume(1000.0)
1395     hypo3d.SetVerboseLevel( 0 )
1396     hypo3d.SetStandardOutputLog( 0 )
1397     hypo3d.SetRemoveLogOnSuccess( 1 )
1398   putName(algo3d.GetSubMesh(), "boiteDefaut")
1399   putName(algo3d, "algo3d_boiteDefaut")
1400   putName(meshBoiteDefaut, "boiteDefaut")
1401
1402   is_done = meshBoiteDefaut.Compute()
1403   text = "meshBoiteDefaut.Compute"
1404   if is_done:
1405     logging.info(text+" OK")
1406   else:
1407     text = "Erreur au calcul du maillage.\n" + text
1408     logging.info(text)
1409     raise Exception(text)
1410
1411   _ = meshBoiteDefaut.GetMesh().UnionListOfGroups( [ group_faceFissOutPipe, group_faceFissInPipe ], 'FACE1' )
1412   maillageSain = enleveDefaut(maillageSain, zoneDefaut, zoneDefaut_skin,
1413                               zoneDefaut_internalFaces, zoneDefaut_internalEdges)
1414   putName(maillageSain, nomFicSain+"_coupe")
1415   _, normfiss = shapeSurFissure(facesPortFissure)
1416   maillageComplet = RegroupeSainEtDefaut(maillageSain, meshBoiteDefaut,
1417                                          None, None, 'COMPLET', normfiss)
1418
1419   logging.info("conversion quadratique")
1420   maillageComplet.ConvertToQuadratic( 1 )
1421   logging.info("groupes")
1422   groups = maillageComplet.GetGroups()
1423   grps = [ grp for grp in groups if grp.GetName() == 'FONDFISS']
1424   _ = maillageComplet.GetMesh().CreateDimGroup( grps, SMESH.NODE, 'FONDFISS' )
1425
1426   logging.info("réorientation face de fissure FACE1")
1427   grps = [ grp for grp in groups if grp.GetName() == 'FACE1']
1428   _ = maillageComplet.Reorient2D( grps[0], normfiss, grps[0].GetID(1))
1429
1430   logging.info("réorientation face de fissure FACE2")
1431   plansim = geompy.MakePlane(O, normfiss, 10000)
1432   fissnorm = geompy.MakeMirrorByPlane(normfiss, plansim)
1433   grps = [ grp for grp in groups if grp.GetName() == 'FACE2']
1434   _ = maillageComplet.Reorient2D( grps[0], fissnorm, grps[0].GetID(1))
1435   _ = maillageComplet.GetMesh().CreateDimGroup( grps, SMESH.NODE, 'FACE2' )
1436
1437   logging.info("export maillage fini")
1438   maillageComplet.ExportMED(fichierMaillageFissure)
1439   putName(maillageComplet, nomFicFissure)
1440   logging.info("fichier maillage fissure %s", fichierMaillageFissure)
1441
1442   if salome.sg.hasDesktop():
1443     salome.sg.updateObjBrowser()
1444
1445   logging.info("maillage fissure fini")
1446
1447   return maillageComplet