1 # -*- coding: iso-8859-1 -*-
2 # This illustrates how to make a mesh partition using the value of a
3 # field defined on this mesh (for example to extract the cells where
4 # the field takes a value greater than a threshold L.
5 # (Anthony Geay, nov. 2012)
7 # WRN: this use case does not require a med input file because the
8 # data (mesh and field) to work with are created from scratch.
10 from MEDCoupling import *
12 # =======================================================
13 # Creation of the input data (mesh and field)
14 # =======================================================
16 # We prepare the input field from scratch instead of load it from a
17 # file, but there is no difference
18 from MEDCouplingDataForTest import MEDCouplingDataForTest
19 m3D=MEDCouplingDataForTest.build3DTargetMesh_1()
21 a=DataArrayDouble.New([1.,1.,1.,-1.,-1.,-1.,-1.,-1.])
22 field=MEDCouplingFieldDouble.New(ON_CELLS,ONE_TIME)
25 field.checkCoherency()
28 # Save the field (and associated mesh)
29 from MEDLoader import MEDLoader
30 MEDLoader.WriteField("partition_input.med",field,True)
32 # =======================================================
33 # Determine the border skin mesh
34 # =======================================================
36 # We have to determine the 2D mesh that delimits the volume where the
37 # field is greater than a threshold L from the volume where the field
38 # is lower than this threshold (in this example L=0).
40 # WRN: This works in SALOME V660 only
44 arr = field.getArray()
45 ids = arr.getIdsInRange(L,1e300)
46 m3DSub = field.getMesh()[ids]
47 skin = m3DSub.computeSkin()
48 MEDLoader.WriteUMesh("partition_skin.med",skin,True);
51 # =======================================================
52 # Compare to the result in SALOME V650
53 # =======================================================
54 # SALOME V650 requires a more complicated syntax.
55 m2D,desc,descI,revDesc,revDescI=m3DSub.buildDescendingConnectivity()
56 numberOf3DVolSharing=revDescI.deltaShiftIndex()
57 ids2D=numberOf3DVolSharing.getIdsEqual(1)
59 # We can check if the two skins are identical
60 print "Are two meshes equal between V660 and V650 ?",skin.isEqual(skin_V650,1e-12)