src/Tools/blocFissure/gmu/meshBlocPart.py

   1 # -*- coding: utf-8 -*-
   2 # Copyright (C) 2014-2022  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 """Maillage du bloc partitionné"""
  21
  22 import logging
  23
  24 import SMESH
  25 from salome.smesh import smeshBuilder
  26 from salome.StdMeshers import StdMeshersBuilder
  27
  28 from .geomsmesh import geompy
  29 from .geomsmesh import smesh
  30
  31 from .putName import putName
  32
  33 # -----------------------------------------------------------------------------
  34
  35 def meshBlocPart(blocPartition, faceFissure, tore, centres, edges, diams, circles, faces, \
  36                 gencnt, facefissoutore, edgeext, facesExternes, facesExtBloc, facesExtElli, \
  37                 aretesInternes, internalBoundary, ellipsoidep, sharedFaces, sharedEdges, edgesBords, \
  38                 nbsegExt, nbsegGen, nbsegRad, scaleRad, reverses, reverext, nbsegCercle, nbsegFis, dmoyen, lensegEllipsoide, \
  39                 mailleur="MeshGems", nro_cas=None):
  40   """Maillage du bloc partitionné"""
  41   logging.info('start')
  42   logging.info("Maillage avec %s pour le cas n°%s", mailleur, nro_cas)
  43
  44   # --- edges de bord à respecter
  45
  46   _ = smesh.CreateFilterManager()
  47   _, internalBoundary, _NoneGroup = internalBoundary.MakeBoundaryElements( SMESH.BND_1DFROM2D, '', '', 0, [ ])
  48   criteres = list()
  49   unCritere = smesh.GetCriterion(SMESH.EDGE,SMESH.FT_FreeBorders,SMESH.FT_Undefined,0)
  50   criteres.append(unCritere)
  51   filtre = smesh.GetFilterFromCriteria(criteres)
  52   bordsLibres = internalBoundary.MakeGroupByFilter( 'bords', filtre )
  53   putName(bordsLibres, 'bordsLibres', i_pref=nro_cas)
  54
  55   # --- maillage bloc
  56
  57   bloc1 = smesh.Mesh(blocPartition)
  58
  59   for i_aux, sharedFaces_i in enumerate(sharedFaces):
  60     algo2d = bloc1.Triangle(algo=smeshBuilder.NETGEN, geom=sharedFaces_i)
  61     putName(algo2d.GetSubMesh(), "sharedFaces", i_aux, nro_cas)
  62     hypo2d = algo2d.Parameters(which=smesh.SIMPLE)
  63     hypo2d.SetLocalLength(lensegEllipsoide)
  64     hypo2d.LengthFromEdges()
  65     hypo2d.SetAllowQuadrangles(0)
  66     putName(hypo2d, "sharedFaces", i_aux, nro_cas)
  67
  68   for i_aux, sharedEdges_i in enumerate(sharedEdges):
  69     algo1d = bloc1.Segment(geom=sharedEdges_i)
  70     putName(algo1d.GetSubMesh(), "sharedEdges", i_aux, nro_cas)
  71     hypo1d = algo1d.LocalLength(lensegEllipsoide)
  72     putName(hypo1d, "sharedEdges={}".format(lensegEllipsoide), i_aux, nro_cas)
  73
  74   declareAlgoEllipsoideFirst = False
  75   if declareAlgoEllipsoideFirst:
  76     algo3d = bloc1.Tetrahedron(algo=smeshBuilder.NETGEN,geom=ellipsoidep)
  77     putName(algo3d.GetSubMesh(), "ellipsoide", i_pref=nro_cas)
  78     hypo3d = algo3d.MaxElementVolume(1000.0)
  79     putName(hypo3d, "ellipsoide", i_pref=nro_cas)
  80
  81   algo3d = bloc1.Prism(geom=tore)
  82   putName(algo3d.GetSubMesh(), "tore", i_pref=nro_cas)
  83   algo2d = bloc1.Quadrangle(geom=tore)
  84   algo1d = bloc1.Segment(geom=tore)
  85   hypo1d = algo1d.NumberOfSegments(nbsegGen)
  86   putName(hypo1d, "tore={}".format(nbsegGen), i_pref=nro_cas)
  87
  88   for i_aux, faces_i in enumerate(faces):
  89     algo2d = bloc1.Quadrangle(geom=faces_i)
  90     putName(algo2d.GetSubMesh(), "faces", i_aux, nro_cas)
  91     hypo2d = smesh.CreateHypothesis('QuadrangleParams')
  92     hypo2d.SetTriaVertex( geompy.GetSubShapeID(blocPartition,centres[i_aux]) )
  93     hypo2d.SetQuadType( StdMeshersBuilder.QUAD_STANDARD )
  94     _ = bloc1.AddHypothesis(hypo2d,faces_i)
  95     putName(hypo2d, "faces", i_aux, nro_cas)
  96
  97   for i_aux, edges_i in enumerate(edges):
  98     algo1d = bloc1.Segment(geom=edges_i)
  99     putName(algo1d.GetSubMesh(), "edges", i_aux, nro_cas)
 100     if reverses[i_aux] > 0:
 101       hypo1d = algo1d.NumberOfSegments(nbsegRad, scaleRad,[ geompy.GetSubShapeID(blocPartition,edges_i) ])
 102     else:
 103       hypo1d = algo1d.NumberOfSegments(nbsegRad, scaleRad,[ ])
 104     putName(hypo1d, "edges", i_aux, nro_cas)
 105
 106   for i_aux, circles_i in enumerate(circles):
 107     algo1d = bloc1.Segment(geom=circles_i)
 108     putName(algo1d.GetSubMesh(), "circles", i_aux, nro_cas)
 109     hypo1d = algo1d.NumberOfSegments(nbsegCercle)
 110     putName(hypo1d, "circles={}".format(nbsegCercle), i_aux, nro_cas)
 111
 112   if len(edgeext) == 1:
 113     densite = int(round(nbsegFis/2))
 114     algo1d = bloc1.Segment(geom=edgeext[0])
 115     putName(algo1d.GetSubMesh(), "edgeext", i_pref=nro_cas)
 116     hypo1d = algo1d.NumberOfSegments(nbsegFis)
 117     hypo1d.SetDistrType( 2 )
 118     hypo1d.SetConversionMode( 1 )
 119     hypo1d.SetTableFunction( [ 0, densite, 0.4, 1, 0.6, 1, 1, densite ] )
 120     putName(hypo1d, "edgeext", i_pref=nro_cas)
 121   else:
 122     longTotal = 0
 123     longEdgeExts = list()
 124     for edgeext_i in edgeext:
 125       props = geompy.BasicProperties(edgeext_i)
 126       longEdgeExts.append(props[0])
 127       longTotal += props[0]
 128     for i_aux, edgeext_i in enumerate(edgeext):
 129       nbLocal = int(round(nbsegFis*longEdgeExts[i_aux]/longTotal))
 130       densite = int(round(nbLocal/2))
 131       algo1d = bloc1.Segment(geom=edgeext_i)
 132       putName(algo1d.GetSubMesh(), "edgeext", i_aux, nro_cas)
 133       hypo1d = algo1d.NumberOfSegments(nbLocal)
 134       hypo1d.SetDistrType( 2 )
 135       hypo1d.SetConversionMode( 1 )
 136       hypo1d.SetTableFunction( [ 0, densite, 0.8, 1, 1, 1 ] )
 137       if reverext[i_aux]:
 138         hypo1d.SetReversedEdges([ geompy.GetSubShapeID(blocPartition, edgeext_i) ])
 139       putName(hypo1d, "edgeext", i_aux, nro_cas)
 140
 141   algo2d = bloc1.Triangle(algo=smeshBuilder.NETGEN_2D, geom=facefissoutore)
 142   putName(algo2d.GetSubMesh(), "facefissoutore", i_pref=nro_cas)
 143   hypo2d = algo2d.LengthFromEdges()
 144   putName(hypo2d, "facefissoutore", i_pref=nro_cas)
 145
 146
 147   maxElemArea = 0.5*dmoyen*dmoyen
 148   logging.debug("dmoyen %s, maxElemArea %s", dmoyen, maxElemArea)
 149
 150   for i_aux, facesExternes_i in enumerate(facesExternes):
 151     algo2d = bloc1.Triangle(algo=smeshBuilder.NETGEN_2D, geom=facesExternes_i)
 152     putName(algo2d.GetSubMesh(), "facesExternes", i_aux, nro_cas)
 153     hypo2d = algo2d.MaxElementArea(maxElemArea)
 154     putName(hypo2d, "facesExternes={}".format(maxElemArea), i_aux, nro_cas)
 155     if edgesBords is None:
 156       algo1d = bloc1.Segment(geom=facesExternes_i)
 157       hypo1d = algo1d.NumberOfSegments(1)
 158       putName(hypo1d, "facesExternes", i_aux, nro_cas)
 159
 160   for i_aux, aretesInternes_i in enumerate(aretesInternes):
 161     algo1d = bloc1.Segment(geom=aretesInternes_i)
 162     putName(algo1d.GetSubMesh(), "aretesInternes", i_aux, nro_cas)
 163     hypo1d = algo1d.NumberOfSegments(nbsegExt)
 164     putName(hypo1d, "aretesInternes={}".format(nbsegExt), i_aux, nro_cas)
 165
 166   if edgesBords is not None:
 167     algo1d = bloc1.UseExisting1DElements(geom=edgesBords)
 168     putName(algo1d.GetSubMesh(), "bordsLibres", i_pref=nro_cas)
 169     hypo1d = algo1d.SourceEdges([ bordsLibres ],0,0)
 170     putName(hypo1d, "bordsLibres", i_pref=nro_cas)
 171
 172   if not declareAlgoEllipsoideFirst:
 173     algo3d = bloc1.Tetrahedron(algo=smeshBuilder.NETGEN,geom=ellipsoidep)
 174     putName(algo3d.GetSubMesh(), "ellipsoide", i_pref=nro_cas)
 175     hypo3d = algo3d.MaxElementVolume(1000.0)
 176     putName(hypo3d, "ellipsoide", i_pref=nro_cas)
 177
 178   _ = bloc1.GroupOnGeom(faceFissure,'FACE1',SMESH.FACE)
 179   _ = bloc1.GroupOnGeom(gencnt,'nfondfis',SMESH.NODE)
 180
 181   groups_faceCommuneEllipsoideBloc = list()
 182   for i_aux, sharedFaces_i in enumerate(sharedFaces):
 183     name = "faceCommuneEllipsoideBloc_{}".format(i_aux)
 184     groups_faceCommuneEllipsoideBloc.append(bloc1.GroupOnGeom(sharedFaces_i, name, SMESH.FACE))
 185   groups_faceExterneBloc = list()
 186   for i_aux, facesExtBloc_i in enumerate(facesExtBloc):
 187     name = "faceExterneBloc_{}".format(i_aux)
 188     groups_faceExterneBloc.append(bloc1.GroupOnGeom(facesExtBloc_i, name, SMESH.FACE))
 189
 190   is_done = bloc1.Compute()
 191   text = "bloc1.Compute"
 192   if is_done:
 193     logging.info(text+" OK")
 194   else:
 195     text = "Erreur au calcul du maillage.\n" + text
 196     logging.info(text)
 197     raise Exception(text)
 198
 199   _ = bloc1.RemoveOrphanNodes()
 200
 201   skinBlocMeshes = list()
 202   for i_aux, groups_faceCommuneEllipsoideBloc_i in enumerate(groups_faceCommuneEllipsoideBloc):
 203     name = "faceCommuneEllipsoideBloc_{}".format(i_aux)
 204     skinBlocMeshes.append(smesh.CopyMesh(groups_faceCommuneEllipsoideBloc_i, name, 0, 0))
 205   for i_aux, groups_faceExterneBloc_i in enumerate(groups_faceExterneBloc):
 206     name = "faceExterneBloc_{}".format(i_aux)
 207     skinBlocMeshes.append(smesh.CopyMesh(groups_faceExterneBloc_i, name, 0, 0))
 208
 209   meshesBloc = [internalBoundary.GetMesh()]
 210   for skinBlocMeshes_i in skinBlocMeshes:
 211     meshesBloc.append(skinBlocMeshes_i.GetMesh())
 212   blocMesh = smesh.Concatenate(meshesBloc, 1, 1, 1e-05,False)
 213
 214   algo3d = blocMesh.Tetrahedron(algo=smeshBuilder.NETGEN)
 215   putName(algo3d.GetSubMesh(), "bloc", i_pref=nro_cas)
 216   hypo3d = algo3d.MaxElementVolume(1000.0)
 217   putName(hypo3d, "bloc", i_pref=nro_cas)
 218
 219   is_done = blocMesh.Compute()
 220   text = "blocMesh.Compute"
 221   if is_done:
 222     logging.info(text+" OK")
 223   else:
 224     text = "Erreur au calcul du maillage.\n" + text
 225     logging.info(text)
 226     raise Exception(text)
 227
 228   blocComplet = smesh.Concatenate([bloc1.GetMesh(), blocMesh.GetMesh()], 1, 1, 1e-05,False)
 229
 230   return bloc1, blocComplet