Hi all,
I have created a 2d mesh with triangular elements in GMSH and have imported it in FEniCS.
There are 4 boundaries in my mesh i.e. free, fixed, top and bottom. These names are given in the .geo file.
I want to specify Dirichlet BC for fixed, Robin BC for free, and Neuman BC for top and bottom boundary.
I don’t want to specify BCs as it given in tutorials eg.
# Defining Dirichlet boundary condition
tol = 1E-14 # tolerance for coordinate comparisons
def Dirichlet_boundary0(x, on_boundary):
return on_boundary and abs(x[0]) < tol
bc = DirichletBC(V, Constant(0), Dirichlet_boundary0)
Because this way will be limited to simple problems, it won’t work if a boundary is a circle or any other shape.
I want to specify the BC by using the physical groups created in the GMSH i.e. fixed, free, top and bottom. For this, I will have to create subdomains for each boundary.
I am attaching the .geo file and the code here.
// Gmsh project created on Wed Jun 16 11:22:49 2021
SetFactory("OpenCASCADE");
//+
Point(1) = {0, 0, 0, 3.0};
//+
Point(2) = {2, 0, 0, 3.0};
//+
Point(3) = {2, 1, 0, 3.0};
//+
Point(4) = {0, 1, 0, 3.0};
//+
Line(1) = {4, 1};
//+
Line(2) = {1, 2};
//+
Line(3) = {2, 3};
//+
Line(4) = {3, 4};
//+
Curve Loop(1) = {1, 2, 3, 4};
//+
Plane Surface(1) = {1};
//+
Physical Curve("fixed", 5) = {1};
//+
Physical Curve("free", 6) = {3};
//+
Physical Curve("top", 7) = {4};
//+
Physical Curve("bottom", 8) = {2};
//+
Physical Surface("fluid", 9) = {1};
This is the code.
from fenics import *
import matplotlib.pyplot as plt
pip install --user meshio
pip install --user h5py
import h5py
import meshio
mesh_from_file = meshio.read("loaded_bar.msh")
import numpy
def create_mesh(mesh, cell_type, prune_z=False):
cells = mesh.get_cells_type(cell_type)
cell_data = mesh.get_cell_data("gmsh:physical", cell_type)
out_mesh = meshio.Mesh(points=mesh.points, cells={cell_type: cells}, cell_data={"name_to_read":[cell_data]})
if prune_z:
out_mesh.prune_z_0()
return out_mesh
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, 2)
Thanks