Integral in Lines

mesh
1

Dear teacher:

  1. I am a learner in FEniCS and Gmsh. Now I want to integral my function in Lines. So I give two domains, this is my Gmsh code:
    //Gmsh Code:
    Geometry.OCCTargetUnit = “M”;

mesh1 = 0.04;
mesh2 = 0.1;
L = 1;
rate = 1/4;

//Subdomain1
p0 = newp; Point(p0) = {0.1,0,0,mesh1};
p1 = newp; Point(p1) = {0.45,rateL/2/2,0,mesh1};
p2 = newp; Point(p2) = {0.1,rateL/2/2,0,mesh1};
p3 = newp; Point(p3) = {0.1,-rateL/2/2,0,mesh1};
p4 = newp; Point(p4) = {0.45,-rateL/2/2,0,mesh1};

l1 = newl; Line(l1) = {p1,p2};
l2 = newl; Line(l2) = {p2,p3};
l3 = newl; Line(l3) = {p3,p4};
l4 = newl; Line(l4) = {p4,p1};
ll1 = newll; Line loop(ll1) = {l1,l2,l3,l4};
s1 = news; Surface(s1) = {ll1};

//Subdomain2
p10 = newp; Point(p10) = {0,0,0,mesh2};
p11 = newp; Point(p11) = {L/2,rateL/2,0,mesh2};
p12 = newp; Point(p12) = {-L/2,rateL/2,0,mesh2};
p13 = newp; Point(p13) = {-L/2,-rateL/2,0,mesh2};
p14 = newp; Point(p14) = {L/2,-rateL/2,0,mesh2};

l11 = newl; Line(l11) = {p11,p12};
l12 = newl; Line(l12) = {p12,p13};
l13 = newl; Line(l13) = {p13,p14};
l14 = newl; Line(l14) = {p14,p11};
ll11 = newll; Line loop(ll11) = {l11,l12,l13,l14};
s11 = news; Surface(s11) = {ll11, ll1};

Physical Surface(1) = {s11,s1};
Physical Line(1) = {l1};
Physical Line(2) = {l11};

  1. Now I want to integral in Lines, So I first put “mesh.geo” into FEniCS,the fenics code is:
    #fenics Code:
    import meshio
    from dolfin import *
    import matplotlib.pyplot as plt
    import numpy

#Mesh
mesh_from_file = meshio.read(“111.msh”)
def create_mesh(mesh, cell_type, prune_z=False):
cells = numpy.vstack([cell.data for cell in mesh.cells if cell.type==cell_type])
cell_data = numpy.hstack([mesh.cell_data_dict[“gmsh:physical”][key]
for key in mesh.cell_data_dict[“gmsh:physical”].keys() if key==cell_type])
# Remove z-coordinates from mesh if we have a 2D cell and all points have the same third coordinate
points= mesh.points
if prune_z:
points = points[:,:2]
mesh_new = meshio.Mesh(points=points, cells={cell_type: cells}, cell_data={“name_to_read”:[cell_data]})
return mesh_new

line_mesh = create_mesh(mesh_from_file, “line”, prune_z=True)
meshio.write(“facet_mesh.xdmf”, line_mesh)
triangle_mesh = create_mesh(mesh_from_file, “triangle”, prune_z=True)
meshio.write(“mesh.xdmf”, triangle_mesh)

mesh = Mesh()
with XDMFFile(“mesh.xdmf”) as infile:
infile.read(mesh)
mvc = MeshValueCollection(“size_t”, mesh, 1)
with XDMFFile(“facet_mesh.xdmf”) as infile:
infile.read(mvc, “name_to_read”)
mf = cpp.mesh.MeshFunctionSizet(mesh, mvc)
ds1 = Measure(“ds”, domain=mesh, subdomain_data=mf)#边界的识别

plot(mesh,title=‘Mesh’)

#Integral
print(‘ds1(1)=’,assemble(1ds1(1)))
print(‘ds1(2)=’,assemble(1ds1(2)))

  1. My Question is:When I want to integral in Inner Line(ds1(1)), the answer is 0;
    When I want to integral in Outer Line(ds1(2)), the answer is 1;It seems that I failed to identify the inner line 1 in fenics. So if I want to integral in ds(1), what can I do?

Thanks for your answer!

2

OK! Just change
ds1 = Measure(“ds”, domain=mesh, subdomain_data=mf)
into
ds1 = Measure(“dS”, domain=mesh, subdomain_data=mf)
Because:
DOLFIN predefines the “measures” dx , ds and dS representing integration over cells, exterior facets (that is, facets on the boundary) and interior facets.
Thanks for your View!