[FEniCS 2019.1.0]

General
1

Hello,
I am trying to rewrite the hyperplasticity demo in terms of the second Piola–Kirchhoff stress sigma (partial of phi wrt partial of Eps ) and the Lagrangian Green strain Eps = 0.5*(F.T*F-I),

referring to this post Minimizing the potential energy in a hyperelasticity problem - #2 by kamensky

The main term on the weak form here would be sigma contraction with dEps, which I computed with dEps = fe.derivative(Eps, u, u_test)

However, I got a error message after just one Newton iteration:

*** -------------------------------------------------------------------------
*** Error:   Unable to successfully call PETSc function 'MatSetValuesLocal'.
*** Reason:  PETSc error code is: 63 (Argument out of range).
*** Where:   This error was encountered inside /Users/runner/miniforge3/conda-bld/fenics-pkgs_1719292218781/work/dolfin/la/PETScMatrix.cpp.
*** Process: 0
*** 
*** DOLFIN version: 2019.1.0
*** Git changeset:  
*** -------------------------------------------------------------------------



import fenics as fe
import ufl
def bottom(x, on_boundary):
    return (on_boundary and fe.near(x[1], 0.0))
def top(x, on_boundary):
    return (on_boundary and fe.near(x[1], 1.0))

E = 20.0e6
mu = 0.5*E
Bulk = E*mu / (3*(3*mu-E))
mesh = fe.UnitSquareMesh(5,5)
V = fe.VectorFunctionSpace(mesh, "CG", 1)
u_tr = fe.TrialFunction(V)
u_test = fe.TestFunction(V)
u = fe.Function(V)
bc = [ fe.DirichletBC(V, fe.Constant((0.0, 0.0)), bottom),  fe.DirichletBC(V, fe.Constant((0.0, 0.1)), top)]

I = fe.Identity(2)
F = I + fe.grad(u)  
C = F.T*F  
J = fe.det(F)  

## Gent
I_1 = fe.tr(C)
phi=( mu * (I_1 - 2 - 3 * fe.ln(J)) + Bulk * fe.ln(J) )


Eps = 0.5*(F.T*F-I)
dEps  = fe.derivative(Eps, u, u_test)
Eps = fe.variable(Eps)
sigma = fe.derivative(phi, Eps)

Form = fe.inner(sigma, dEps) * fe.dx      

Jac = fe.derivative(Form, u, u_tr)

problem = fe.NonlinearVariationalProblem(Form, u, bc, Jac)
solver = fe.NonlinearVariationalSolver(problem)
solver.solve()
2

Please provide a more descriptive title to your question.

Also, please ensure that your problem is runnable, as I get the following error by running your code:

Calling FFC just-in-time (JIT) compiler, this may take some time.
Calling FFC just-in-time (JIT) compiler, this may take some time.
Traceback (most recent call last):
  File "mwe.py", line 30, in <module>
    sigma = fe.derivative(phi, Eps)
  File "/usr/local/lib/python3.6/dist-packages/dolfin/fem/formmanipulations.py", line 61, in derivative
    form_arguments = form.arguments()
AttributeError: 'Sum' object has no attribute 'arguments'

because you have not multiplied phi by an integration measure.#
Are you sure you don’t want to use ufl.diff?
See: Form language — Unified Form Language (UFL) 2021.1.0 documentation

3

Hi dokken,

Sorry it was a typo, sigma = fe.diff(phi, Eps) is correct, it wouldn’t let me edit this post to correct that though.
*it also wouldn’t let me update the title

4

Here is the corrected MWE code that produces the above error:

import fenics as fe
import ufl
def bottom(x, on_boundary):
    return (on_boundary and fe.near(x[1], 0.0))
def top(x, on_boundary):
    return (on_boundary and fe.near(x[1], 1.0))

E = 20.0e6
mu = 0.5*E
Bulk = E*mu / (3*(3*mu-E))
mesh = fe.UnitSquareMesh(5,5)
V = fe.VectorFunctionSpace(mesh, "CG", 1)
u_tr = fe.TrialFunction(V)
u_test = fe.TestFunction(V)
u = fe.Function(V)
bc = [ fe.DirichletBC(V, fe.Constant((0.0, 0.0)), bottom),  fe.DirichletBC(V, fe.Constant((0.0, 0.1)), top)]

I = fe.Identity(2)
F = I + fe.grad(u)  
C = F.T*F  
J = fe.det(F)  

## Gent
I_1 = fe.tr(C)
phi=( mu * (I_1 - 2 - 3 * fe.ln(J)) + Bulk * fe.ln(J) )


Eps = 0.5*(F.T*F-I)
dEps  = fe.derivative(Eps, u, u_test)
Eps = fe.variable(Eps)
sigma = fe.diff(phi, Eps)

Form = fe.inner(sigma, dEps) * fe.dx      

Jac = fe.derivative(Form, u, u_tr)

problem = fe.NonlinearVariationalProblem(Form, u, bc, Jac)
solver = fe.NonlinearVariationalSolver(problem)
solver.solve()
5

I think the error lies in sigma = fe.diff(phi, Eps), phi was not explicitly defined from Eps, but I can’t think of any work around at the moment

6

Well, it has to if you want to differentiate through it.
Otherwise, how should one do the differentiation, if the variable you are differentiating with respect to isn’t formed by the expression above.

7

Thanks dokken!

I figured out a way to reformulate phi explicitly in terms of Eps and it works now. I get mixed up with different autodiff sometimes. for example autograd in PyTorch is able to compute the implicit derivatives of phi wrt Eps as long as they were both defined from the same u. Thank you very much for helping though!