Modeling Elastohydrodynamic Lubrication in Fenics- Coupled PDEs

variational formulation
1

Hi,
I am trying to solve a problem in elastohydrodynamic lubrication using Fenics, and I am trying to figure out if the following system can be implemented in Fenics.

image
The configuration is shown in figure above has the following formulation.

Equation 1 :Elasticity equation

\nabla \cdot \sigma=0 ,\text{on} \, \Omega

with \sigma=C:\epsilon where C is the elastic stiffness. In this the unknown is the displacement field U

Equation 2 : Advection Diffusion type equation on {\Omega_c}.This is a 1D formulation

\nabla \cdot (\bar{\varepsilon} \nabla P)-\frac{\partial (\bar{\rho} H)}{\partial X}=0 \, ,\text{on} \, \Omega_c

This equation is also a free boundary problem with BCs P=0 on the extents of \Omega_c and laso a constraint P>=0 on f \Omega_c which needs a penalty formulation.The parameters P , H refer to pressure and the film thickness, which are both unknown apriori. The film thickness typically uses the trial solution

H=H_0 +X^2/2+W

where H_0 is unknown(we need to take an initial assumption) W is the deflection of the top boundary that comes from the solution of Eqn 1.

Equation 3 : This is a force balance constraint which can be written as

\int_{\Omega_C} PdX =\frac{\pi}{2}

this is the equation that allows that needs to be solved for H_0.

For the other bounds of the domain, the traction free BC is imposed on \partial \Omega_4 and
\partial \Omega_5 , and the boundaries \partial \Omega_1 ,\partial \Omega_2 and \partial \Omega_3 are all fixed

I was able to look at the literature in this area, and one popular apporoach to solve this found in
the thesis attached,(https://laur.lau.edu.lb:8443/xmlui/handle/10725/6860), where the authors pose this to be a set of weak forms in Galerkin:
Find (P,U,H_0) such that one has

Equation\, 1:\int_{\Omega} - C{\epsilon}(U) \cdot {\epsilon}(W_U) d\Omega + \int_{\Omega_c} -PW_U d\Omega_c =0 \\ Equation\, 2:\int_{\Omega_c} - \bar{\varepsilon} \nabla P\cdot \nabla W_p d\Omega_c + \int_{\Omega_c} \bar{\rho} H\frac{\partial W_p }{\partial X} d\Omega_c =0 \\ Equation\, 3:\int_{\Omega_C} Pd\Omega_c -\frac{\pi}{2}=0

I am curious as to know if such a system can be solved in Fenics, instead of writing my own solver in MATLAB or python. Could anyone help me a bit as to how to go about modeling this?

2

I dont know if you managed this but i tried Habchi’s formulations in Comsol and it was very difficult to achieve. I failed a lot then somehow got there with some fiddling.
I would also like to use an open source solution for the EHL problem, would be interested if you got there in the end.

3

Hi delta_dubo! I am currently simulating EHL problems in Comsol. Can you guide me on how to implement Habachi formulation in Comsol? You can also reply to me at imran.sajid@hitecuni.edu.pk