Use of 'derivative' command in nonlinear problems

NB: I’m very new in the FEniCS world, so please have patience and if you can point me out some other tutorials or guides to follow or, even better, strategies to debug my code I would appreciate it very much. Thanks in advance

Dear all,
I would like to ask a clarification about the following command

derivative

I found it in this tutorial about hyperelasticity in the dolfin guide
https://fenicsproject.org/docs/dolfin/latest/python/demos/hyperelasticity/demo_hyperelasticity.py.html
and I would like to understand if I can use it also in my case

My actual problem is that I’m writing the work principle as functional to be minimized
B \int_{\Omega}\chi \:\delta\chi dx + other terms =0
B is bending stifness

where the rotation is
\theta= \arctan{\frac{u_y'(s,t)}{1+u_x'(s,t)}}
and the curvature is
\chi=\frac{d\theta}{ds}
a non linear function of the derivatives of the displacement, u_x and u_y, with respect s \in [0,1] (curvilinear abscissa)

and what I would like to do is to use derivative in order to obtain the following variation of \chi
\delta \chi=f1(u_x',u_y',u_x'',u_y'')\delta u_x'+f2(u_x',u_y',u_x'',u_y'')\delta u_y'+f3(u_x',u_y',u_x'',u_y'')\delta u_x''+f4(u_x',u_y',u_x'',u_y'')\delta u_y''

is that possible?

and I’m coding it in the right way?

from dolfin import *

mesh = UnitIntervalMesh(10) 			# s \in [0,1]
P1 = FiniteElement('P', interval , 2) 	# http://www.femtable.org/ ?hermitiane?
el2 = MixedElement([P1, P1])
V = FunctionSpace(mesh, el2)

u_ = TestFunction(V)
ux_, uy_ = split(u_)

u = Function(V)
ux, uy = split(u)

# Define rotation,curvature, constraints
def theta(ux, uy):
	return atan(div(uy)/(1+div(ux)))

def chi(ux, uy):
	return div(theta(ux, uy))

# Variation through the derivative function
def dchi(ux, uy):
	return derivative(chi(ux, uy), u, u_)

Thanks in advance to all the community
Please, if something is not clear, let me know.

It looks like your understanding of derivative is correct, although, assuming your original energy functional is

E = \frac{1}{2}B\int_\Omega\chi^2\,d\mathbf{x} + \cdots \text{ ,}

you could alternatively just apply deriative to that directly instead of doing the first step of the derivative yourself, since

\delta E = B\int_\Omega \chi\,\delta\chi\, d\mathbf{x} + \cdots\text{ .}

However, a bigger problem here is that your energy functional is not well-defined when \mathbf{u} is in a standard C^0 finite element space. The function space defined by FiniteElement('P', interval , 2) is not a C^1 Hermite spline spline as suggested in the comments. It is only C^0 at element boundaries, so the second derivatives of \mathbf{u} in \chi include Dirac deltas at element boundaries, so \chi^2 is not well-defined.

To use this formulation as-is, you would need a C^1 function space, which is currently not implemented in the standard FEniCS toolchain. I believe there is some support for C^1 elements in Firedrake, which has a similar Python API, and I also wrote a small library myself, tIGAr, which allows you to compute with smoother splines in FEniCS, albeit with slightly different usage.

1 Like

@kamensky thank you very much for your answer and for your useful comments.

I’m definitely going to check your paper but I have already seen that is way too advanced and technical for the competence I have in the field.

To be honest, I have actually tried to install tIGAr but, unfortunately, I encountered some issue during the installation process :exploding_head:

I will let you know.
Thanks again,
Andrea