Hello everyone,
I want to extract markers from gmsh model. so I added two physical group tagged 101,102 as following:
gmsh.initialize()
gmsh.option.setNumber("General.Terminal", 1)
gmsh.model.add("modelo_6")
geom = gmsh.model.geo # Line to abbreviate geometry commands
L = 1 # Length from the center of the inner mesh to the edges
# %% Spline control points are created:
tm=4
geom.addPoint(1, 0, 0,tm)
geom.addPoint(0.8, 0.027, 0,tm)
geom.addPoint(0.55, 0.058, 0,tm)
geom.addPoint(0.3, 0.08, 0,tm)
geom.addPoint(0.15, 0.078, 0,tm)
geom.addPoint(0, 0.052, 0,tm)
geom.addPoint(0, 0, 0,tm)
geom.addPoint(0, -0.01, 0,tm)
geom.addPoint(0.15, -0.052, 0,tm)
geom.addPoint(0.3, -0.05, 0,tm)
geom.addPoint(0.55, -0.043, 0,tm)
geom.addPoint(0.8, -0.017, 0,tm)
# The B-Spline is created by taking the previous points as control points:
af = geom.addBSpline ([1 , 2 , 3 , 4 , 5 , 6 , 7 , 8 , 9 , 10 , 11 , 12 , 1 ])
N=151
geom.mesh.setTransfiniteCurve(af, N)
cint = geom.addCurveLoop([af]) # Curve loop interior
# %% The outer surface border is created:
p1 = geom.addPoint(0.5-L, -L, 0)
p2 = geom.addPoint(0.5+L, -L, 0)
p3 = geom.addPoint ( 0.5 + L , L , 0 )
p4 = geom.addPoint ( 0.5 - L , L , 0 )
l1 = geom.addLine(p1, p2)
l2 = geom.addLine(p2, p3)
l3 = geom.addLine(p3, p4)
l4 = geom.addLine(p4, p1)
cext = geom.addCurveLoop([l1, l2, l3, l4]) # Curve loop exterior
s1 = geom.addPlaneSurface ([ cext , cint ]) # Create the outer surface
s2 = geom.addPlaneSurface ([cint ]) # Create the outer surface
# Create physical surface
gmsh.model.addPhysicalGroup(2, [s1], tag=101)
gmsh.model.setPhysicalName(2, 101, "out")
gmsh.model.addPhysicalGroup(2, [s2], tag=102)
gmsh.model.setPhysicalName(2, 102, "in")
# %% Finally the mesh is created and saved:
geom.synchronize()
gmsh.model.mesh.generate(2)
Then I got the markers using io.extract_gmsh_topology_and_markers:
gdim = 2
element_types=2
geometry_data = io.extract_gmsh_geometry(gmsh.model)
topology_data=io.extract_gmsh_topology_and_markers(gmsh.model)
ufl_domain = io.ufl_mesh_from_gmsh(element_types, gdim)
cells= topology_data[2]['topology']
marker=topology_data[2]['cell_data']
# Create distributed mesh
domain = mesh.create_mesh(MPI.COMM_WORLD,cells, geometry_data[:, :gdim], ufl_domain)
and assigned different values to the different regions:
V = fem.FunctionSpace(domain, ("DG", 0))
ub = fem.Function(V)
ub.vector.array[marker==102]=1
but, as it can be seen below, the markers end up with wrong regions (the region tagged 102 is supposed to be inside of airfoil). What is the correct way to separate different regions?
