Robin boundary condition for 2D problem

General
1

Hello everyone,
I am trying to implement Robin boundary condition in 2D problem. The schematic for the problem is given below. In the problem 𝜎_12= 𝛼_1 𝑒_1 and 𝜎_22= 𝛼_2 𝑒_2 where 𝑒_1 and 𝑒_2 are components of the vector displacement field u to be solved. Any suggestions or help on how to add the stress terms in the variational form ?

schem_robin
schem_robin599Γ—489 4.26 KB

2

Have you written down the weak form on paper?

Do you know how that would work for the case of a Neumann boundary condition? Can you start from there and see what is the difference in your case?

This is a forum about FEniCS specific support, we are not supposed to answer in detail about all possible background you may be missing on the finite element method.

3

Hello Francesco,

I have the weak form written down as
image
The term on the right side can be written as
image
I know how to implement Neumann boundary condition. So implement this in FEniCSx, I tried the following code:

Coeff = ufl.as_tensor([[alpha1, 0], [0, alpha_2]])
u = ufl.TrialFunction(V)
v = ufl.TestFunction(V)
Tr_top = ufl.dot(Coeff, u)

a = ufl.inner(sigma(u), epsilon(v))*ufl.dx - ufl.dot(Tr_top, v) *ufl.ds

Is this the correct way to implement it ?

4

It follows the same logic as when you write the problem on paper.
So yes, you would include the term

in your variational form.

Please note that since you are using a as a variable name, I am assuming that you are trying to set up a linear system (a(u,v)=L(v)).

Then you should either write everything in a residual form, such as

F = ufl.inner(sigma(u), epsilon(v))*ufl.dx - ufl.dot(Tr_top, v) *ufl.ds

and use

a, L = ufl.system(F)

to get the bi-linear and linear form, or just write them out by hand:

a = ufl.inner(sigma(u), epsilon(v))*ufl.dx 
L = ufl.dot(Tr_top, v) *ufl.ds
1 Like
5

@viveksingh What you have written down is not a Robin condition though. Are you sure that is what you need?

In your context, a Neumann condition would be:
\sigma\cdot n = T_{top} \,\, on \,\,\partial\Omega_{top}
That condition is substituted into the boundary integral that shows up in working out the weak formulation and then moved to the right hand side. And @dokken has provided you the code with how to implement that in FEniCS.

A Robin condition, though, would look like:
\sigma\cdot n + \alpha u = T_{top} \,\, on \,\,\partial\Omega_{top}
Substituting that into the boundary integral gives you an additional term in your left hand side.

Augmenting Dokken’s earlier code snippet:

a = ufl.inner(sigma(u), epsilon(v))*ufl.dx + ufl.dot(alpha * u, v) *ufl.ds
L = ufl.dot(Tr_top, v) *ufl.ds
1 Like