- Thank you received: 24
Unstructed mesh with structured subdomain
5 years 1 month ago - 5 years 1 month ago #1975
by christopher
Replied by christopher on topic Unstructed mesh with structured subdomain
It is, but not straight forward:
create a mesh for the outside, where the boundaries of the top and bottom edge of the enclosed subdomain have copied edges. Then you need to fill the interior with a structured mesh. The following code uses some knowledge of how edges and so on are treated internally and moving code lines may result in unexpected behaviour
Hope this helps
Best
Christopher
create a mesh for the outside, where the boundaries of the top and bottom edge of the enclosed subdomain have copied edges. Then you need to fill the interior with a structured mesh. The following code uses some knowledge of how edges and so on are treated internally and moving code lines may result in unexpected behaviour
Code:
from netgen.geom2d import *
from netgen.meshing import *
geo = SplineGeometry()
geo.AddRectangle([0,0],[1,1])
sd = ((0.33,0.33),(0.66,0.66))
pts = [geo.AppendPoint(*p) for p in [sd[0], (sd[1][0], sd[0][1]), sd[1], (sd[0][1], sd[1][0])]]
index = {}
index["bot"] = geo.Append( ["line", pts[0], pts[1]], leftdomain=0, rightdomain=1, maxh=0.01)
index["top"] = geo.Append( ["line", pts[3], pts[2]], copy=index["bot"], maxh=0.01)
index["right"] = geo.Append( ["line", pts[1], pts[2]], leftdomain=0, rightdomain=1)
index["left"] = geo.Append( ["line", pts[0], pts[3]], copy=index["right"])
geo.SetDomainMaxH(1,0.1)
mesh = geo.GenerateMesh()
# inverse map with index shift
inv_map = {}
for key,val in index.items():
inv_map[val+1] = key
segments = {key : [] for key in index}
for el in mesh.Elements1D():
if el.index in inv_map:
segments[inv_map[el.index]].append(el.vertices)
pnts = {}
for key, val in segments.items():
pnts[key] = [mesh[el[0]] for el in val]
pnts[key].append(mesh[val[-1][1]])
rows = [[seg[0] for seg in segments["bot"]]]
rows[0].append(segments["bot"][-1][1])
for i in range(1,len(segments["left"])):
new_row = [segments["left"][i][0]]
for j in range(1, len(segments["bot"])):
new_row.append(mesh.Add(MeshPoint(Pnt(pnts["bot"][j][0], pnts["left"][i][1], 0))))
new_row.append(segments["right"][i][0])
rows.append(new_row)
rows.append([seg[0] for seg in segments["top"]])
rows[-1].append(segments["top"][-1][1])
mesh.Add(FaceDescriptor(surfnr = 2, bc=2))
for i in range(len(rows)-1):
for j in range(len(rows[0])-1):
el = Element2D(index=2, vertices=[rows[i][j], rows[i][j+1], rows[i+1][j+1], rows[i+1][j]])
mesh.Add(el)
from ngsolve import *
msh = Mesh(mesh)
Draw(msh)
Hope this helps
Best
Christopher
Attachments:
Time to create page: 0.120 seconds