# \[FEniCS 2019.1.0\]

**URL:** <https://fenicsproject.discourse.group/t/fenics-2019-1-0/15993>\
**Category:** General\
**Created:** [October 7, 2024, 6:35pm UTC](https://fenicsproject.discourse.group/t/fenics-2019-1-0/15993 "2024-10-07T18:35:31Z")\
**Posts on this page:** 7\
**Page:** 1

<div class="post-metadata">

**Author:** ![a516095046](https://avatars.discourse-cdn.com/v4/letter/a/e36b37/32.png) [@a516095046](https://fenicsproject.discourse.group/u/a516095046)\
**Post date:** [October 7, 2024, 6:35pm UTC](https://fenicsproject.discourse.group/t/fenics-2019-1-0/15993/1 "2024-10-07T18:35:31Z")

</div>

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](https://fenicsproject.discourse.group/t/minimizing-the-potential-energy-in-a-hyperelasticity-problem/10175/2)

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:

```auto
*** -------------------------------------------------------------------------
*** 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:  
*** -------------------------------------------------------------------------

```

```auto
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()

```

---

<div class="post-metadata">

**Author:** ![dokken](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/dokken/32/1560_2.png) [@dokken](https://fenicsproject.discourse.group/u/dokken)\
**Post date:** [October 7, 2024, 9:25pm UTC](https://fenicsproject.discourse.group/t/fenics-2019-1-0/15993/2 "2024-10-07T21:25:03Z")

</div>

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:

```bash
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](https://fenics.readthedocs.io/projects/ufl/en/latest/manual/form_language.html#variable-derivatives)

---

<div class="post-metadata">

**Author:** ![a516095046](https://avatars.discourse-cdn.com/v4/letter/a/e36b37/32.png) [@a516095046](https://fenicsproject.discourse.group/u/a516095046)\
**Post date:** [October 8, 2024, 12:56am UTC](https://fenicsproject.discourse.group/t/fenics-2019-1-0/15993/3 "2024-10-08T00:56:55Z")

</div>

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

---

<div class="post-metadata">

**Author:** ![a516095046](https://avatars.discourse-cdn.com/v4/letter/a/e36b37/32.png) [@a516095046](https://fenicsproject.discourse.group/u/a516095046)\
**Post date:** [October 8, 2024, 12:59am UTC](https://fenicsproject.discourse.group/t/fenics-2019-1-0/15993/4 "2024-10-08T00:59:49Z")

</div>

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

```auto
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()

```

---

<div class="post-metadata">

**Author:** ![a516095046](https://avatars.discourse-cdn.com/v4/letter/a/e36b37/32.png) [@a516095046](https://fenicsproject.discourse.group/u/a516095046)\
**Post date:** [October 8, 2024, 1:02am UTC](https://fenicsproject.discourse.group/t/fenics-2019-1-0/15993/5 "2024-10-08T01:02:18Z")

</div>

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

---

<div class="post-metadata">

**Author:** ![dokken](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/dokken/32/1560_2.png) [@dokken](https://fenicsproject.discourse.group/u/dokken)\
**Post date:** [October 8, 2024, 6:38am UTC](https://fenicsproject.discourse.group/t/fenics-2019-1-0/15993/6 "2024-10-08T06:38:47Z")

</div>

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.

---

<div class="post-metadata">

**Author:** ![a516095046](https://avatars.discourse-cdn.com/v4/letter/a/e36b37/32.png) [@a516095046](https://fenicsproject.discourse.group/u/a516095046)\
**Post date:** [October 8, 2024, 2:45pm UTC](https://fenicsproject.discourse.group/t/fenics-2019-1-0/15993/7 "2024-10-08T14:45:43Z")

</div>

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!