Hello,
I am trying to define one of my boundaries as a symmetry boundary. I know how to define it mathematically as the dot product of velocity and normal vector are zero along that boundary. However, I have not been able to find any resource on how to actually implement this in FEniCS.
Any help or guidance is greatly appreciated.
Thanks!
Hi,
if V is your velocity function space and your symmetry plane is let’s say x=0 then just impose your DirichletBC on V.sub(0)
My question is more focused on how to actually write out that Dirichlet BC in Fenics, my experience is limited to defining actual velocities on the boundary, i.e. setting velocity to (5,0) at an inlet. I am not sure how I would actually write out an expression with my solution space in it, as a BC.
I would suggest something like this: considering a 2D flow in a channel like the Poiseulle flow through parallel plates - provided four boundaries (left, right, top, bottom), if the bottom boundary is given as symmetry, it would go as follows:
from dolfin import *
from mshr import *
channel = Rectangle(Point(0,0), Point(2, 1))
mesh = generate_mesh(channel, 50)
# Boundaries
left= 'near(x[0], 0)'
right= 'near(x[0], 2)'
bottom = 'near(x[1], 0)'
top= 'near(x[1], 1)'
# bottom is symmetry
symmetry = DirichletBC(V.sub(1), Constant(0.0), bottom)
noslip = DirichletBC(V, Constant((0, 0)), top)
# for a pressure-difference driven flow
P_inlet = DirichletBC(Q, Constant(3.0), left)
P_outlet = DirichletBC(Q, Constant(0.0), right)
# If flow rate is imposed
# Ux = 0.1
# U_inlet = DirichletBC(V, Constant((Ux, 0)), left)
# Collect bcs for pressure-driven flow
bcu = [noslip , symmetry]
bcp = [P_inlet, P_outlet]
# Collect bcs for flow-rate-driven flow
# bcu = [noslip , symmetry, U_inlet]
# bcp = [P_outlet]
Provided the variational form for your problem, this approach should work to implement symmetry boundary. I use it for my applications.
Hello,
In your answer, you indicate imposing Dirichlet directly on the second component of the velocity, V.sub(1), which implies that you are not using a vector valued function for the velocity since you can access each velocity component individually. How do you do that?
I can’t access the individual components of the velocity field to impose Dirichlet only on the transversal velocity… I am very new at Fenics, if you could point me in the right direction I would really appreciate it!
For example, if I am solving Stokes equations, my spaces look like
V = VectorElement("CG", triangle, 2)
Q = FiniteElement("CG", triangle, 1)
W = FunctionSpace(mesh, MixedElement([V, Q]))In my reply I should have defined the function spaces used in the problem. The code snippet from above should include the lines below before defining the boundary condition:
V = VectorFunctionSpace(mesh, 'P', 2)
Q = FunctionSpace(mesh, 'P', 1)
From your code snippet I see that you are using a mixed element formulation, which means W.sub(0) is 2D (i.e. velocity functions will be based there) and W.sub(1) is 1D.
To impose a specific value for the x component of the velocity for DirichletBC using mixed elements I would access the component as below:
from dolfin import *
from mshr import *
channel = Rectangle(Point(0,0), Point(2, 1))
mesh = generate_mesh(channel, 50)
# Boundaries
left= 'near(x[0], 0)'
right= 'near(x[0], 2)'
bottom = 'near(x[1], 0)'
top= 'near(x[1], 1)'
V = VectorElement("CG", triangle, 2)
Q = FiniteElement("CG", triangle, 1)
W = FunctionSpace(mesh, MixedElement([V, Q]))
noslip = DirichletBC(W.sub(0), Constant((0, 0)), top) # no-slip using mixed elements
# addressing a value for the y-component of velocity
symmetry = DirichletBC(W.sub(0).sub(1), Constant(0.0), bottom)
P_inlet = DirichletBC(W.sub(1), Constant(3.0), left)
P_outlet = DirichletBC(W.sub(1), Constant(0.0), right)
# Collect bcs for pressure-driven flow
bcu = [noslip , symmetry]
bcp = [P_inlet, P_outlet]
I hope I addressed your question. I am prone to errors (as everyone), so any suggestions are welcome!
Thank you so much, this is exactly what I wanted! I just need (a lot) more practice with Fenics, I’m still so unfamiliar that I don’t know what we’re allowed to do.
Using W.sub(0).sub(1) makes a lot of sense!
Hello, in the case here I want to implement a symmetric BC on boundary 1. I have the same mixed space as above:
Is this the correct syntax?
symmetry = DirichletBC(W.sub(0).sub(0), Constant(0.0), Boundary1)
To my (very limited) knowledge, yes. if by symmetric you mean symmetry condition, so u*n = 0 at that boundary.
The way I think of it is that W.sub(0) ‘accesses’ the vector space for velocity, and then W.sub(0).sub(0) corresponds to the horizontal component while W.sub(0).sub(1) gets to the vertical component of the velocity. Since your boundary 1 is vertical, the normal direction aligns with (minus) the horizontal component, so setting W.sub(0).sub(0) to zero is the way to go.
