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