Salome HOME
Update of CheckDone
[modules/smesh.git] / src / SMESH_SWIG / smeshBuilder.py
index 25a947352351a2e561162e588452f021f2a864e9..2401ed74cfcc7be98c4e92706c120e297250a1c7 100644 (file)
@@ -1,4 +1,4 @@
-# Copyright (C) 2007-2022  CEA/DEN, EDF R&D, OPEN CASCADE
+# Copyright (C) 2007-2024  CEA, EDF, OPEN CASCADE
 #
 # This library is free software; you can redistribute it and/or
 # modify it under the terms of the GNU Lesser General Public
@@ -462,6 +462,22 @@ class smeshBuilder( SMESH._objref_SMESH_Gen, object ):
             obj,name = name,obj
         return Mesh(self, self.geompyD, obj, name)
 
+    def ParallelMesh(self, obj, name=0, split_geom=True):
+        """
+        Create a parallel mesh.
+
+        Parameters:
+            obj: geometrical object for meshing
+            name: the name for the new mesh.
+            split_geom: If True split the geometry and create the assoicated
+            sub meshes
+
+        Returns:
+            an instance of class :class:`ParallelMesh`.
+        """
+        return ParallelMesh(self, self.geompyD, obj,
+                            split_geom=split_geom, name=name)
+
     def RemoveMesh( self, mesh ):
         """
         Delete a mesh
@@ -789,7 +805,6 @@ class smeshBuilder( SMESH._objref_SMESH_Gen, object ):
         """
         if isinstance( mesh, Mesh ):
             mesh = mesh.GetMesh()
-        print("calling createdualmesh from Python")
         dualMesh = SMESH._objref_SMESH_Gen.CreateDualMesh(self, mesh, meshName, adaptToShape)
         return Mesh(self, self.geompyD, dualMesh)
 
@@ -1188,6 +1203,8 @@ class smeshBuilder( SMESH._objref_SMESH_Gen, object ):
             functor = aFilterMgr.CreateAspectRatio3D()
         elif theCriterion == FT_Warping:
             functor = aFilterMgr.CreateWarping()
+        elif theCriterion == FT_Warping3D:
+            functor = aFilterMgr.CreateWarping3D()
         elif theCriterion == FT_MinimumAngle:
             functor = aFilterMgr.CreateMinimumAngle()
         elif theCriterion == FT_Taper:
@@ -1218,6 +1235,8 @@ class smeshBuilder( SMESH._objref_SMESH_Gen, object ):
             functor = aFilterMgr.CreateNodeConnectivityNumber()
         elif theCriterion == FT_BallDiameter:
             functor = aFilterMgr.CreateBallDiameter()
+        elif theCriterion == FT_ScaledJacobian:
+            functor = aFilterMgr.CreateScaledJacobian()
         else:
             print("Error: given parameter is not numerical functor type.")
         aFilterMgr.UnRegister()
@@ -1578,7 +1597,7 @@ class Mesh(metaclass = MeshMeta):
     mesh = 0
     editor = 0
 
-    def __init__(self, smeshpyD, geompyD, obj=0, name=0):
+    def __init__(self, smeshpyD, geompyD, obj=0, name=0, parallel=False):
 
         """
         Constructor
@@ -1611,7 +1630,12 @@ class Mesh(metaclass = MeshMeta):
                     else:
                         geo_name = "%s_%s to mesh"%(self.geom.GetShapeType(), id(self.geom)%100)
                     geompyD.addToStudy( self.geom, geo_name )
-                self.SetMesh( self.smeshpyD.CreateMesh(self.geom) )
+                if parallel and isinstance(self, ParallelMesh):
+                    mymesh = self.smeshpyD.CreateParallelMesh(self.geom)
+                    mymesh2 = mymesh._narrow(SMESH._objref_SMESH_Mesh)
+                    self.SetMesh( mymesh )
+                else:
+                    self.SetMesh( self.smeshpyD.CreateMesh(self.geom) )
 
             elif isinstance(obj, SMESH._objref_SMESH_Mesh):
                 self.SetMesh(obj)
@@ -1642,7 +1666,7 @@ class Mesh(metaclass = MeshMeta):
         Destructor. Clean-up resources
         """
         if self.mesh:
-            self.mesh.UnRegister()
+            #self.mesh.UnRegister()
             pass
         pass
 
@@ -1863,31 +1887,6 @@ class Mesh(metaclass = MeshMeta):
                 geom = self.geom
         return self.smeshpyD.Evaluate(self.mesh, geom)
 
-    def ParallelCompute(self, nbThreads, mesherNbThreads=1, geom=0, discardModifs=False, refresh=False):
-        """
-        Parallel computation of the mesh and return the status of the computation
-        The mesh must contains have be constructed using create_parallel_mesh
-
-        Parameters:
-                nbThreads: Number of threads to use for a parallel computation
-                geom: geomtrical shape on which mesh data should be computed
-                discardModifs: if True and the mesh has been edited since
-                        a last total re-compute and that may prevent successful partial re-compute,
-                        then the mesh is cleaned before Compute()
-                refresh: if *True*, Object Browser is automatically updated (when running in GUI)
-
-        Returns:
-                True or False
-        """
-        if (nbThreads <= 1):
-            raise ValueError("nbThreads must be strictly greater than 1")
-        if (mesherNbThreads < 1):
-            raise ValueError("nbThreads must be greater than 1")
-
-        self.mesh.SetMesherNbThreads(mesherNbThreads)
-        self.mesh.SetNbThreads(nbThreads)
-        return self.Compute(geom=geom, discardModifs=discardModifs, refresh=refresh)
-
     def Compute(self, geom=0, discardModifs=False, refresh=False):
         """
         Compute the mesh and return the status of the computation
@@ -1898,7 +1897,6 @@ class Mesh(metaclass = MeshMeta):
                         a last total re-compute and that may prevent successful partial re-compute,
                         then the mesh is cleaned before Compute()
                 refresh: if *True*, Object Browser is automatically updated (when running in GUI)
-                nbThreads: Number of threads to use for a parallel computation
 
         Returns:
                 True or False
@@ -2003,6 +2001,13 @@ class Mesh(metaclass = MeshMeta):
 
         return ok
 
+    def CheckCompute(self):
+        """
+        Check if the mesh was properly compute
+        """
+        if not self.mesh.IsComputedOK():
+            raise Exception("Could not compute {}".format(self.GetName()))
+
     def GetComputeErrors(self, shape=0 ):
         """
         Return a list of error messages (:class:`SMESH.ComputeError`) of the last :meth:`Compute`
@@ -5471,6 +5476,32 @@ class Mesh(metaclass = MeshMeta):
         if mesh: mesh = self.smeshpyD.Mesh(mesh)
         return mesh, group
 
+    def MakeBoundaryOfEachElement(self, groupName="", meshName="", toCopyAll=False, groups=[] ):
+        """
+        Create boundary elements around the whole mesh or groups of elements
+
+        Parameters:
+                groupName: a name of group to store all boundary elements in,
+                        "" means not to create the group
+                meshName: a name of a new mesh, which is a copy of the initial
+                        mesh + created boundary elements; "" means not to create the new mesh
+                toCopyAll: if True, the whole initial mesh will be copied into
+                        the new mesh else only boundary elements will be copied into the new mesh
+                groups: list of :class:`sub-meshes, groups or filters <SMESH.SMESH_IDSource>` of elements to make boundary around
+
+        Returns:
+                tuple( long, mesh, group )
+                       - long - number of added boundary elements
+                       - mesh - the :class:`Mesh` where elements were added to
+                       - group - the :class:`group <SMESH.SMESH_Group>` of boundary elements or None
+        """
+        dimension=SMESH.BND_2DFROM3D
+        toCreateAllElements = True # create all boundary elements in the mesh
+        nb, mesh, group = self.editor.MakeBoundaryElements( dimension,groupName,meshName,
+                                                           toCopyAll,toCreateAllElements,groups)
+        if mesh: mesh = self.smeshpyD.Mesh(mesh)
+        return nb, mesh, group
+
     def MakeBoundaryElements(self, dimension=SMESH.BND_2DFROM3D, groupName="", meshName="",
                              toCopyAll=False, groups=[]):
         """
@@ -5494,9 +5525,9 @@ class Mesh(metaclass = MeshMeta):
                        - mesh - the :class:`Mesh` where elements were added to
                        - group - the :class:`group <SMESH.SMESH_Group>` of boundary elements or None
         """
-
+        toCreateAllElements = False # create only elements in the boundary of the solid
         nb, mesh, group = self.editor.MakeBoundaryElements(dimension,groupName,meshName,
-                                                           toCopyAll,groups)
+                                                           toCopyAll,toCreateAllElements,groups)
         if mesh: mesh = self.smeshpyD.Mesh(mesh)
         return nb, mesh, group
 
@@ -7438,6 +7469,19 @@ class Mesh(metaclass = MeshMeta):
 
         return self.FunctorValue(SMESH.FT_Warping, elemId)
 
+    def GetWarping3D(self, elemId):
+        """
+        Get warping angle of faces element of 3D elements.
+
+        Parameters:
+            elemId: mesh element ID
+
+        Returns:
+            element's warping angle value
+        """
+
+        return self.FunctorValue(SMESH.FT_Warping3D, elemId)
+
     def GetMinimumAngle(self, elemId):
         """
         Get minimum angle of 2D element.
@@ -7477,6 +7521,19 @@ class Mesh(metaclass = MeshMeta):
 
         return self.FunctorValue(SMESH.FT_Skew, elemId)
 
+    def GetScaledJacobian(self, elemId):
+        """
+        Get the scaled jacobian of 3D element id
+
+        Parameters:
+            elemId: mesh element ID
+
+        Returns:
+            the scaled jacobian
+        """
+
+        return self.FunctorValue(SMESH.FT_ScaledJacobian, elemId)
+
     def GetMinMax(self, funType, meshPart=None):
         """
         Return minimal and maximal value of a given functor.
@@ -7511,6 +7568,447 @@ class Mesh(metaclass = MeshMeta):
 
     pass # end of Mesh class
 
+def _copy_gmsh_param(dim, local_param, global_param):
+    if dim==3:
+        local_param.SetMaxSize(global_param.GetMaxSize())
+        local_param.SetMinSize(global_param.GetMinSize())
+        local_param.Set3DAlgo(global_param.Get3DAlgo())
+        local_param.SetRecombineAll(global_param.GetRecombineAll())
+        local_param.SetSubdivAlgo(global_param.GetSubdivAlgo())
+        local_param.SetRemeshAlgo(global_param.GetRemeshAlgo())
+        local_param.SetRemeshPara(global_param.GetRemeshPara())
+        local_param.SetSmouthSteps(global_param.GetSmouthSteps())
+        local_param.SetSizeFactor(global_param.GetSizeFactor())
+        local_param.SetUseIncomplElem(global_param.GetUseIncomplElem())
+        local_param.SetMeshCurvatureSize(global_param.GetMeshCurvatureSize())
+        local_param.SetSecondOrder(global_param.GetSecondOrder())
+        local_param.SetIs2d(global_param.GetIs2d())
+    elif dim==2:
+        local_param.SetMaxSize(global_param.GetMaxSize())
+        local_param.SetMinSize(global_param.GetMinSize())
+        local_param.Set2DAlgo(global_param.Get2DAlgo())
+        local_param.SetRecomb2DAlgo(global_param.GetRecomb2DAlgo())
+        local_param.SetRecombineAll(global_param.GetRecombineAll())
+        local_param.SetSubdivAlgo(global_param.GetSubdivAlgo())
+        local_param.SetRemeshAlgo(global_param.GetRemeshAlgo())
+        local_param.SetRemeshPara(global_param.GetRemeshPara())
+        local_param.SetSmouthSteps(global_param.GetSmouthSteps())
+        local_param.SetSizeFactor(global_param.GetSizeFactor())
+        local_param.SetUseIncomplElem(global_param.GetUseIncomplElem())
+        local_param.SetMeshCurvatureSize(global_param.GetMeshCurvatureSize())
+        local_param.SetSecondOrder(global_param.GetSecondOrder())
+        local_param.SetIs2d(global_param.GetIs2d())
+
+def _copy_netgen_param(dim, local_param, global_param):
+    """
+    Create 1D/2D/3D netgen parameters from a NETGEN 1D2D3D parameter
+    """
+    if dim==1:
+        #TODO: More global conversion ? or let user define it
+        local_param.NumberOfSegments(int(global_param.GetMaxSize()))
+    elif dim==2:
+        local_param.SetMaxSize(global_param.GetMaxSize())
+        local_param.SetMinSize(global_param.GetMinSize())
+        local_param.SetOptimize(global_param.GetOptimize())
+        local_param.SetFineness(global_param.GetFineness())
+        local_param.SetNbSegPerEdge(global_param.GetNbSegPerEdge())
+        local_param.SetNbSegPerRadius(global_param.GetNbSegPerRadius())
+        #TODO: Why the 0.9 to have same results
+        local_param.SetGrowthRate(global_param.GetGrowthRate()*0.9)
+        local_param.SetChordalError(global_param.GetChordalError())
+        local_param.SetChordalErrorEnabled(global_param.GetChordalErrorEnabled())
+        local_param.SetUseSurfaceCurvature(global_param.GetUseSurfaceCurvature())
+        local_param.SetUseDelauney(global_param.GetUseDelauney())
+        local_param.SetQuadAllowed(global_param.GetQuadAllowed())
+        local_param.SetWorstElemMeasure(global_param.GetWorstElemMeasure())
+        local_param.SetCheckChartBoundary(global_param.GetCheckChartBoundary())
+        local_param.SetNbThreads(global_param.GetNbThreads())
+    else:
+        local_param.SetMaxSize(global_param.GetMaxSize())
+        local_param.SetMinSize(global_param.GetMinSize())
+        local_param.SetOptimize(global_param.GetOptimize())
+        local_param.SetCheckOverlapping(global_param.GetCheckOverlapping())
+        local_param.SetCheckChartBoundary(global_param.GetCheckChartBoundary())
+        local_param.SetFineness(global_param.GetFineness())
+        local_param.SetNbSegPerEdge(global_param.GetNbSegPerEdge())
+        local_param.SetNbSegPerRadius(global_param.GetNbSegPerRadius())
+        local_param.SetGrowthRate(global_param.GetGrowthRate())
+        local_param.SetNbThreads(global_param.GetNbThreads())
+
+
+def _shaperstudy2geom(geompyD, shaper_obj):
+    """
+    Convertion of shaper object to geom object
+
+    Parameters:
+        geompyD: geomBuilder instance
+        shaper_obj: Shaper study object
+
+    Returns:
+        geom object
+
+    """
+    import tempfile
+    #Writing shaperstudy object into a brep file
+    fid, tmp_file = tempfile.mkstemp(suffix='.brep')
+    with open(fid, 'wb') as f:
+        f.write(shaper_obj.GetShapeStream())
+    # Reimporting brep file into geom
+    real_geom = geompyD.ImportBREP(tmp_file)
+    os.remove(tmp_file)
+
+    return real_geom
+
+
+def _split_geom(geompyD, geom):
+    """
+    Splitting geometry into n solids and a 2D/1D compound
+
+    Parameters:
+        geompyD: geomBuilder instance
+        geom: geometrical object for meshing
+
+    Returns:
+        compound containing all the 1D,2D elements
+        list of solids
+    """
+
+    # Splitting geometry into 3D elements and all the 2D/1D into one compound
+    object_solids = geompyD.ExtractShapes(geom, geompyD.ShapeType["SOLID"],
+                                          True)
+
+    solids = []
+    isolid = 0
+    for solid in object_solids:
+        isolid += 1
+        geompyD.addToStudyInFather( geom, solid, 'Solid_{}'.format(isolid) )
+        solids.append(solid)
+    # If geom is a solid ExtractShapes will return nothin in that case geom is the solids
+    if isolid == 0:
+       solids = [geom]
+
+    faces = []
+    iface = 0
+    solid_faces = geompyD.ExtractShapes(geom, geompyD.ShapeType["FACE"],
+                                        True)
+    for face in solid_faces:
+        faces.append(face)
+        iface += 1
+        geompyD.addToStudyInFather(geom, face,
+                                   'Face_{}'.format(iface))
+
+    return faces, solids
+
+
+MULTITHREAD, MULTINODE = range(2)
+class ParallelismSettings:
+    """
+    Defines the parameters for the parallelism of ParallelMesh
+    """
+    def __init__(self, mesh):
+        """
+        Construsctor
+
+        Parameters:
+        mesh: Instance of ParallelMesh
+        """
+        if not(isinstance(mesh, ParallelMesh)):
+            raise ValueError("mesh should be a ParallelMesh")
+
+        self._mesh = mesh
+
+
+class MTParallelismSettings(ParallelismSettings):
+    """
+    Defines the parameters for the parallelism of ParallelMesh using MultiThreading
+    """
+    def __init__(self, mesh):
+        ParallelismSettings.__init__(self, mesh)
+
+    # Multithreading methods
+    def SetNbThreads(self, nbThreads):
+        """ Set the number of threads for multithread """
+        if nbThreads < 1:
+            raise ValueError("Number of threads must be stricly greater than 1")
+
+        self._mesh.mesh.SetNbThreads(nbThreads)
+
+    def GetNbThreads(self):
+        """ Get Number of threads """
+        return self._mesh.mesh.GetNbThreads()
+
+    def __str__(self):
+        """ str conversion """
+        string = "\nParameter for MultiThreading parallelism:\n"
+        string += "NbThreads: {}\n".format(self.GetNbThreads())
+
+        return string
+
+
+class MNParallelismSettings(ParallelismSettings):
+    """
+    Defines the parameters for the parallelism of ParallelMesh using MultiNodal
+    """
+    def __init__(self, mesh):
+        ParallelismSettings.__init__(self, mesh)
+
+    def GetResource(self):
+        """ Get the resource on which to run """
+        return self._mesh.mesh.GetResource()
+
+    def SetResource(self, resource):
+        """ Set the resource on which to run """
+        self._mesh.mesh.SetResource(resource)
+
+    def SetNbProc(self, nbProc):
+        """ Set the number of Processor for multinode """
+        if nbProc < 1:
+            raise ValueError("Number of Proc must be stricly greater than 1")
+        self._mesh.mesh.SetNbProc(nbProc)
+
+    def GetNbProc(self):
+        """ Get Number of Processor """
+        return self._mesh.mesh.GetNbProc()
+
+    def SetNbProcPerNode(self, nbProcPerNode):
+        """ Set the number of Processor Per Node for multinode """
+        if nbProcPerNode < 1:
+            raise ValueError("Number of Processor Per Node must be stricly greater than 1")
+
+        self._mesh.mesh.SetNbProcPerNode(nbProcPerNode)
+
+    def GetNbProcPerNode(self):
+        """ Get Number of Processor Per Node """
+        return self._mesh.mesh.GetNbProcPerNode()
+
+    def SetNbNode(self, nbNode):
+        """ Set the number of Node for multinode """
+        if nbNode < 1:
+            raise ValueError("Number of Node must be stricly greater than 1")
+        self._mesh.mesh.SetNbNode(nbNode)
+
+    def GetNbNode(self):
+        """ Get Number of Node """
+        return self._mesh.mesh.GetNbNode()
+
+    def SetWcKey(self, wcKey):
+        """ Set the number of Node for multinode """
+        self._mesh.mesh.SetWcKey(wcKey)
+
+    def GetWcKey(self):
+        """ Get Number of Node """
+        return self._mesh.mesh.GetWcKey()
+
+    def SetWalltime(self, walltime):
+        """ Set the number of Node for multinode """
+        self._mesh.mesh.SetWalltime(walltime)
+
+    def GetWalltime(self):
+        """ Get Number of Node """
+        return self._mesh.mesh.GetWalltime()
+
+    def __str__(self):
+        """ str conversion """
+        string = "\nParameter for MultiNode parallelism:\n"
+        string += "Reource: {}\n".format(self.GetResource())
+        string += "NbProc: {}\n".format(self.GetNbProc())
+        string += "NbProcPerNode: {}\n".format(self.GetNbProcPerNode())
+        string += "NbNode: {}\n".format(self.GetNbNode())
+        string += "WcKey: {}\n".format(self.GetWcKey())
+        string += "Walltime: {}\n".format(self.GetWalltime())
+
+        return string
+
+
+class ParallelMesh(Mesh):
+    """
+    Surcharge on Mesh for parallel computation of a mesh
+    """
+    def __init__(self, smeshpyD, geompyD, geom, split_geom=True, name=0):
+        """
+        Create a parallel mesh.
+
+        Parameters:
+            smeshpyD: instance of smeshBuilder
+            geompyD: instance of geomBuilder
+            geom: geometrical object for meshing
+            split_geom: If true will divide geometry on solids and 1D/2D
+            coumpound and create the associated submeshes
+            name: the name for the new mesh.
+
+        Returns:
+            an instance of class :class:`ParallelMesh`.
+        """
+
+        if not isinstance(geom, geomBuilder.GEOM._objref_GEOM_Object):
+            raise ValueError("geom argument must be a geometry")
+
+        try:
+            import SHAPERSTUDY
+            shaper_object = SHAPERSTUDY.SHAPERSTUDY_ORB._objref_SHAPER_Object
+            has_shaper = True
+        except ImportError:
+            shaper_object = int
+            has_shaper = False
+
+        # If we have a shaper object converting it into geom (temporary solution)
+        if isinstance(geom, shaper_object):
+            self._geom_obj = _shaperstudy2geom(geompyD, geom)
+        elif isinstance(geom, geomBuilder.GEOM._objref_GEOM_Object):
+            self._geom_obj = geom
+        else:
+            msg= ""
+            if not has_shaper:
+                msg = "\nShaper was not compiled"
+            raise Exception("Could not handle geom format {}.{} ".format(type(geom), msg))
+
+        # Splitting geometry into one geom containing 1D and 2D elements and a
+        # list of 3D elements
+        super(ParallelMesh, self).__init__(smeshpyD, geompyD, self._geom_obj, name, parallel=True)
+
+        if split_geom:
+            self._faces, self._solids = _split_geom(geompyD, self._geom_obj)
+
+        self._param = None
+
+    def _build_submeshes(self, mesher2D, mesher3D):
+        """
+        Contruct the submeshes for a parallel use of smesh
+
+        Parameters:
+            mesher2D: name of 2D mesher for 2D parallel compute (NETGEN)
+            mesher3D: name of 3D mesher for 3D parallel compute (NETGEN or
+            GMSH)
+        """
+
+        # Building global 2D mesher
+        if mesher3D:
+            if mesher3D == "NETGEN":
+                algo2D = "NETGEN_2D"
+            elif mesher3D == "GMSH":
+                algo2D = "GMSH_2D"
+            else:
+                raise ValueError("mesher3D should be either NETGEN or GMSH")
+
+            self._algo2d = self.Triangle(geom=self._geom_obj, algo=algo2D)
+
+        # Parallel 2D
+        if mesher2D:
+            #Means that we want to mesh face of solids in parallel and not
+            #the volume
+            self._algo2d = []
+            #For the moment use AutomaticLength based on finesse
+            # TODO: replace by input hypothesis
+            self._algo1d = self.Segment(geom=self._geom_obj)
+
+            for face_id, face in enumerate(self._faces):
+                name = "face_{}".format(face_id)
+                algo2d = self.Triangle(geom=face, algo="NETGEN_2D_Remote")
+                self._algo2d.append(algo2d)
+
+        if mesher3D:
+            self._algo3d = []
+            for solid_id, solid in enumerate(self._solids):
+                name = "Solid_{}".format(solid_id)
+                if ( mesher3D == "NETGEN" ):
+                    algo3d = self.Tetrahedron(geom=solid, algo="NETGEN_3D_Remote")
+                    self._algo3d.append(algo3d)
+                elif ( mesher3D == "GMSH" ):
+                    algo3d = self.Tetrahedron(geom=solid, algo="GMSH_3D_Remote")
+                    self._algo3d.append(algo3d)
+
+    def GetNbSolids(self):
+        """
+        Return the number of 3D solids
+        """
+        return len(self._solids)
+
+    def GetNbFaces(self):
+        """
+        Return the number of 2D faces
+        """
+        return len(self._faces)
+
+    def GetParallelismMethod(self):
+        """ Get the parallelims method """
+        return self.mesh.GetParallelismMethod()
+
+    def SetParallelismMethod(self, method):
+        """ Set the parallelims method """
+        if method not in [MULTITHREAD , MULTINODE]:
+            raise ValueError("Parallelism method can only be 0:MultiThread or 1:MultiNode")
+
+        self.mesh.SetParallelismMethod(method)
+
+        if method == MULTITHREAD:
+            self._param = MTParallelismSettings(self)
+        else:
+            self._param = MNParallelismSettings(self)
+
+    def GetParallelismSettings(self):
+        """
+        Return class to set parameters for the parallelism
+        """
+        if self._param is None:
+            raise Exception("You need to set Parallelism method first (SetParallelismMethod)")
+        return self._param
+
+    def AddGlobalHypothesis(self, hyp):
+        """
+        Split hypothesis to apply it to all the submeshes:
+        - the 1D+2D
+        - each of the 3D solids
+
+        Parameters:
+                hyp: a hypothesis to assign
+
+        """
+        if isinstance(hyp, NETGENPlugin._objref_NETGENPlugin_Hypothesis):
+            copy_param = _copy_netgen_param
+            mesher3D = "NETGEN"
+        elif isinstance(hyp, GMSHPlugin._objref_GMSHPlugin_Hypothesis):
+            copy_param = _copy_gmsh_param
+            mesher3D = "GMSH"
+        else:
+            raise ValueError("param must come from NETGENPlugin or GMSHPlugin")
+
+        self.mesh.SetParallelismDimension(3)
+        self._build_submeshes(None, mesher3D)
+
+        param2d = self._algo2d.Parameters()
+        copy_param(2, param2d, hyp)
+
+        for algo3d in self._algo3d:
+            param3d = algo3d.Parameters()
+            copy_param(3, param3d, hyp)
+
+    def Add2DGlobalHypothesis(self, hyp):
+        """
+        Split hypothesis to apply it to all the submeshes:
+        - the 1D
+        - each of the 2D faces
+
+        Parameters:
+                hyp: a hypothesis to assign
+
+        """
+        if isinstance(hyp, NETGENPlugin._objref_NETGENPlugin_Hypothesis):
+            copy_param = _copy_netgen_param
+            mesher2D = "NETGEN"
+        else:
+            raise ValueError("param must come from NETGENPlugin")
+
+        self.mesh.SetParallelismDimension(2)
+        self._build_submeshes(mesher2D, None)
+
+        param1d = self._algo1d
+        copy_param(1, param1d, hyp)
+
+        for algo2d in self._algo2d:
+            param2d = algo2d.Parameters()
+            copy_param(2, param2d, hyp)
+
+    pass # End of ParallelMesh
 
 class meshProxy(SMESH._objref_SMESH_Mesh):
     """
@@ -7554,10 +8052,20 @@ class meshProxy(SMESH._objref_SMESH_Mesh):
         if len( args ) == 1:
             args += True,
         return SMESH._objref_SMESH_Mesh.ExportDAT(self, *args)
-    pass
+
 omniORB.registerObjref(SMESH._objref_SMESH_Mesh._NP_RepositoryId, meshProxy)
 
 
+class parallelMeshProxy(SMESH._objref_SMESH_ParallelMesh):
+    def __init__(self,*args):
+        SMESH._objref_SMESH_ParallelMesh.__init__(self,*args)
+    def __deepcopy__(self, memo=None):
+        new = self.__class__(self)
+        return new
+omniORB.registerObjref(SMESH._objref_SMESH_ParallelMesh._NP_RepositoryId, parallelMeshProxy)
+
+
+
 class submeshProxy(SMESH._objref_SMESH_subMesh):
 
     """