Spuriuos oscillations with curvature related forces

Question

Hallo
I have to deal with problems including surface tension. This forcing is proportional to the curvature of the interface. The treatment is as follows.
$\int_{\Gamma} w^{t} k \gamma d\Gamma = \int_{Gamma}\gamma t \cdot (t \cdot \nabla)w d\Gamma + b.t.$
Where k is the local vector curvature, w is the test function and $\gamma$ is the surface tension coefficient, t is the tangent vector to the interface.
The ufl code I wrote is as follows.

P0 = FiniteElement("DG",triangle,0)
P1 = FiniteElement("DG",triangle,0)
P2 = VectorElement("Lagrange",triangle,2)

V = MixedElement([P2,P1])

(s, ps) = TrialFunctions(V)
(w, qs) = TestFunctions(V)

gamma = Constant(triangle)
kmesh = Constant(triangle)

n = FacetNormal(triangle)

ss = Coefficient(P2)
pps = Coefficient(P1)

t=as_vector([n[1], -n[0]])

F = -(\
((kmesh*w[j].dx(k)*ss[j].dx(k)\
+kmesh*w[j].dx(k)*ss[k].dx(j)\
-w[i].dx(i)*pps\
-qs*ss[i].dx(i)))*dx\
-gamma*t[k]('+')*avg(w[j]).dx(k)*t[j]('+')*dS(2)\
)

J = -derivative(F,(ss, pps), (s, ps))

forms=[F, J]

The exact solution for a circular domain should be a jump in pressure that equals the surface tension coefficient $\gamma$, but I find spurious overshoots in pressure just in the cells close to the interface. Unfortunately I am not able to put her an image of them.

Thanks for any suggestion.
Hallo
Stefano

Answer 1

I have done a similar thing in FEniCS using a variational approach and automatic differentiation for an hyperelastic material with surface tension.

You can find the paper with the FEniCS code in the supplemental material here

http://journals.aps.org/prl/abstract/10.1103/PhysRevLett.111.114301#supplemental

If you do not have access to the review I can send you the paper and the code

Corrado

Comments

Thank you very much. I read the article. It is very interesting and quite similar to what I'm doing. I'm searching the steady state solution of a bubble, so I have to solve both the internal and external fields. The boundary force acts so on a boundary immersed in the domain. What I posted is just a test that reproduces the error. I have no access to the code I would be very grate to you if you could send it to me. stefano.pedraglio@libero.it.
Thank you very much.