### (Deleted) Instability of nonlinear elastodynamics

68

views

0

I played around with the undocumented elastodynamics demo and had some numerical issues. I changed the definition of the stress (for the hyperelastic Saint-Venant Kirchhoff model) but kept almost everything else the same. With this material model, the equations are nonlinear, so instead of forms 'a' and 'L', now I have 'F' and 'J' with 'u' being a Coefficient. Here are the changed parts in the ufl:

In the cpp part, the appropriate values are assigned to the forms 'F' and 'J'. Following the Cahn-Hilliard demo, I use the Newton method in each timestep to solve the problem. These are the modified lines of the cpp code:

The results seem to be accurate for some timesteps but after some time the Newton iteration fails. It seems to me, that the problem occurs when the deformation changes direction. I tried to reduce the timestep but the problem remained. What can be the reason for this instability? As far as I am concerned, the generalized-$\alpha$

Here is the deformation of a cantilever beam demonstrating the problem:

I post the Newton iteration before the error and the error:

Solving nonlinear variational problem.

Newton iteration 0: r (abs) = 3.049e-02 (tol = 1.000e-10) r (rel) = 1.000e+00 (tol = 1.000e-10)

Newton iteration 1: r (abs) = 6.863e-03 (tol = 1.000e-10) r (rel) = 2.251e-01 (tol = 1.000e-10)

Newton iteration 2: r (abs) = 5.976e-04 (tol = 1.000e-10) r (rel) = 1.960e-02 (tol = 1.000e-10)

Newton iteration 3: r (abs) = 1.934e-05 (tol = 1.000e-10) r (rel) = 6.342e-04 (tol = 1.000e-10)

Newton iteration 4: r (abs) = 4.702e-08 (tol = 1.000e-10) r (rel) = 1.542e-06 (tol = 1.000e-10)

Newton iteration 5: r (abs) = 3.160e-13 (tol = 1.000e-10) r (rel) = 1.036e-11 (tol = 1.000e-10)

Newton solver finished in 5 iterations and 5 linear solver iterations.

Time: 6.825

Solving nonlinear variational problem.

Newton iteration 0: r (abs) = 3.174e-02 (tol = 1.000e-10) r (rel) = 1.000e+00 (tol = 1.000e-10)

Newton iteration 1: r (abs) = 8.671e-03 (tol = 1.000e-10) r (rel) = 2.732e-01 (tol = 1.000e-10)

Newton iteration 2: r (abs) = 1.097e-03 (tol = 1.000e-10) r (rel) = 3.456e-02 (tol = 1.000e-10)

Newton iteration 3: r (abs) = 2.184e-04 (tol = 1.000e-10) r (rel) = 6.883e-03 (tol = 1.000e-10)

Newton iteration 4: r (abs) = 9.723e-03 (tol = 1.000e-10) r (rel) = 3.064e-01 (tol = 1.000e-10)

Newton iteration 5: r (abs) = 2.484e-03 (tol = 1.000e-10) r (rel) = 7.826e-02 (tol = 1.000e-10)

Newton iteration 6: r (abs) = 6.795e-04 (tol = 1.000e-10) r (rel) = 2.141e-02 (tol = 1.000e-10)

Newton iteration 7: r (abs) = 2.479e-04 (tol = 1.000e-10) r (rel) = 7.809e-03 (tol = 1.000e-10)

Newton iteration 8: r (abs) = 4.503e-04 (tol = 1.000e-10) r (rel) = 1.419e-02 (tol = 1.000e-10)

Newton iteration 9: r (abs) = 2.143e-04 (tol = 1.000e-10) r (rel) = 6.752e-03 (tol = 1.000e-10)

Newton iteration 10: r (abs) = 2.080e-02 (tol = 1.000e-10) r (rel) = 6.552e-01 (tol = 1.000e-10)

Newton iteration 11: r (abs) = 5.253e-03 (tol = 1.000e-10) r (rel) = 1.655e-01 (tol = 1.000e-10)

Newton iteration 12: r (abs) = 1.369e-03 (tol = 1.000e-10) r (rel) = 4.313e-02 (tol = 1.000e-10)

Newton iteration 13: r (abs) = 4.055e-04 (tol = 1.000e-10) r (rel) = 1.278e-02 (tol = 1.000e-10)

Newton iteration 14: r (abs) = 2.144e-04 (tol = 1.000e-10) r (rel) = 6.754e-03 (tol = 1.000e-10)

Newton iteration 15: r (abs) = 1.833e-02 (tol = 1.000e-10) r (rel) = 5.775e-01 (tol = 1.000e-10)

Newton iteration 16: r (abs) = 4.636e-03 (tol = 1.000e-10) r (rel) = 1.461e-01 (tol = 1.000e-10)

Newton iteration 17: r (abs) = 1.215e-03 (tol = 1.000e-10) r (rel) = 3.828e-02 (tol = 1.000e-10)

Newton iteration 18: r (abs) = 3.686e-04 (tol = 1.000e-10) r (rel) = 1.161e-02 (tol = 1.000e-10)

Newton iteration 19: r (abs) = 2.193e-04 (tol = 1.000e-10) r (rel) = 6.911e-03 (tol = 1.000e-10)

Newton iteration 20: r (abs) = 2.149e-03 (tol = 1.000e-10) r (rel) = 6.770e-02 (tol = 1.000e-10)

Newton iteration 21: r (abs) = 5.965e-04 (tol = 1.000e-10) r (rel) = 1.879e-02 (tol = 1.000e-10)

Newton iteration 22: r (abs) = 2.324e-04 (tol = 1.000e-10) r (rel) = 7.322e-03 (tol = 1.000e-10)

Newton iteration 23: r (abs) = 7.238e-04 (tol = 1.000e-10) r (rel) = 2.281e-02 (tol = 1.000e-10)

Newton iteration 24: r (abs) = 2.568e-04 (tol = 1.000e-10) r (rel) = 8.091e-03 (tol = 1.000e-10)

Newton iteration 25: r (abs) = 3.830e-04 (tol = 1.000e-10) r (rel) = 1.207e-02 (tol = 1.000e-10)

Newton iteration 26: r (abs) = 2.167e-04 (tol = 1.000e-10) r (rel) = 6.827e-03 (tol = 1.000e-10)

Newton iteration 27: r (abs) = 4.022e-03 (tol = 1.000e-10) r (rel) = 1.267e-01 (tol = 1.000e-10)

Newton iteration 28: r (abs) = 1.062e-03 (tol = 1.000e-10) r (rel) = 3.346e-02 (tol = 1.000e-10)

Newton iteration 29: r (abs) = 3.324e-04 (tol = 1.000e-10) r (rel) = 1.047e-02 (tol = 1.000e-10)

Newton iteration 30: r (abs) = 2.328e-04 (tol = 1.000e-10) r (rel) = 7.335e-03 (tol = 1.000e-10)

Newton iteration 31: r (abs) = 7.113e-04 (tol = 1.000e-10) r (rel) = 2.241e-02 (tol = 1.000e-10)

Newton iteration 32: r (abs) = 2.542e-04 (tol = 1.000e-10) r (rel) = 8.010e-03 (tol = 1.000e-10)

Newton iteration 33: r (abs) = 3.992e-04 (tol = 1.000e-10) r (rel) = 1.258e-02 (tol = 1.000e-10)

Newton iteration 34: r (abs) = 2.148e-04 (tol = 1.000e-10) r (rel) = 6.769e-03 (tol = 1.000e-10)

Newton iteration 35: r (abs) = 1.068e-02 (tol = 1.000e-10) r (rel) = 3.367e-01 (tol = 1.000e-10)

Newton iteration 36: r (abs) = 2.726e-03 (tol = 1.000e-10) r (rel) = 8.588e-02 (tol = 1.000e-10)

Newton iteration 37: r (abs) = 7.394e-04 (tol = 1.000e-10) r (rel) = 2.330e-02 (tol = 1.000e-10)

Newton iteration 38: r (abs) = 2.600e-04 (tol = 1.000e-10) r (rel) = 8.193e-03 (tol = 1.000e-10)

Newton iteration 39: r (abs) = 3.653e-04 (tol = 1.000e-10) r (rel) = 1.151e-02 (tol = 1.000e-10)

Newton iteration 40: r (abs) = 2.201e-04 (tol = 1.000e-10) r (rel) = 6.936e-03 (tol = 1.000e-10)

Newton iteration 41: r (abs) = 1.896e-03 (tol = 1.000e-10) r (rel) = 5.975e-02 (tol = 1.000e-10)

Newton iteration 42: r (abs) = 5.343e-04 (tol = 1.000e-10) r (rel) = 1.683e-02 (tol = 1.000e-10)

Newton iteration 43: r (abs) = 2.226e-04 (tol = 1.000e-10) r (rel) = 7.015e-03 (tol = 1.000e-10)

Newton iteration 44: r (abs) = 1.393e-03 (tol = 1.000e-10) r (rel) = 4.390e-02 (tol = 1.000e-10)

Newton iteration 45: r (abs) = 4.114e-04 (tol = 1.000e-10) r (rel) = 1.296e-02 (tol = 1.000e-10)

Newton iteration 46: r (abs) = 2.141e-04 (tol = 1.000e-10) r (rel) = 6.745e-03 (tol = 1.000e-10)

Newton iteration 47: r (abs) = 3.521e-02 (tol = 1.000e-10) r (rel) = 1.109e+00 (tol = 1.000e-10)

Newton iteration 48: r (abs) = 8.856e-03 (tol = 1.000e-10) r (rel) = 2.791e-01 (tol = 1.000e-10)

Newton iteration 49: r (abs) = 2.269e-03 (tol = 1.000e-10) r (rel) = 7.149e-02 (tol = 1.000e-10)

Newton iteration 50: r (abs) = 6.262e-04 (tol = 1.000e-10) r (rel) = 1.973e-02 (tol = 1.000e-10)

terminate called after throwing an instance of 'std::runtime_error'

what():

*** -------------------------------------------------------------------------

*** DOLFIN encountered an error. If you are not able to resolve this issue

*** using the information listed below, you can ask for help at

***

*** fenics-support@googlegroups.com

***

*** Remember to include the error message listed below and, if possible,

*** include a *minimal* running example to reproduce the error.

***

*** -------------------------------------------------------------------------

*** Error: Unable to solve nonlinear system with NewtonSolver.

*** Reason: Newton solver did not converge because maximum number of iterations reached.

*** Where: This error was encountered inside NewtonSolver.cpp.

*** Process: 0

***

*** DOLFIN version: 2017.2.0

*** Git changeset: unknown

*** -------------------------------------------------------------------------

```
# Stress tensor
def sigma(r):
I = Identity(r.geometric_dimension())
F0 = I + grad(r)
E = 0.5*(F0.T*F0-I)
return lmbda*tr(E)*I + 2.0*mu*E
# Linear form
E1 = factor_m1*inner(r, u)*dx \
+factor_d1*inner(r, u)*dx \
+(1.0-alpha_f)*inner(grad(r), sigma(u))*dx
E2 = factor_m1*inner(r, u0)*dx + factor_m2*inner(r, v0)*dx + factor_m3*inner(r, a0)*dx \
+ factor_d1*inner(r, u0)*dx + factor_d2*inner(r, v0)*dx + factor_d3*inner(r, a0)*dx \
- alpha_f*inner(grad(r), sigma(u0))*dx \
+ inner(r, f)*dx + (1.0-alpha_f)*inner(r, p)*ds(3) + alpha_f*inner(r, p0)*ds(3)
F = E1 - E2
# Jacobian
du = TrialFunction(element)
J = derivative(F, u, du)
```

In the cpp part, the appropriate values are assigned to the forms 'F' and 'J'. Following the Cahn-Hilliard demo, I use the Newton method in each timestep to solve the problem. These are the modified lines of the cpp code:

```
// Create solver Parameters
Parameters params("nonlinear_variational_solver");
Parameters newton_params("newton_solver");
newton_params.add("relative_tolerance", 1e-10);
params.add(newton_params);
// Start time stepping
std::size_t step = 0;
while (t < T)
{
// Update for next time step
t += *dt;
cout << "Time: " << t << endl;
// Solve nonlinear variational problem
solve(F == 0, *u, bc, J, params);
// Update velocity and acceleration
update_a(a, *u, *a0, *v0, *u0, *beta, *dt);
update_v(v, a, *a0, *v0, *gamma, *dt);
*u0 = *u; *v0 = v; *a0 = a;
Save solutions to file
if (step % 1 0 == 0)
{
file_u << *u;
}
++step;
}
```

The results seem to be accurate for some timesteps but after some time the Newton iteration fails. It seems to me, that the problem occurs when the deformation changes direction. I tried to reduce the timestep but the problem remained. What can be the reason for this instability? As far as I am concerned, the generalized-$\alpha$

`α`method should work for nonlinear problems as well.Here is the deformation of a cantilever beam demonstrating the problem:

I post the Newton iteration before the error and the error:

Solving nonlinear variational problem.

Newton iteration 0: r (abs) = 3.049e-02 (tol = 1.000e-10) r (rel) = 1.000e+00 (tol = 1.000e-10)

Newton iteration 1: r (abs) = 6.863e-03 (tol = 1.000e-10) r (rel) = 2.251e-01 (tol = 1.000e-10)

Newton iteration 2: r (abs) = 5.976e-04 (tol = 1.000e-10) r (rel) = 1.960e-02 (tol = 1.000e-10)

Newton iteration 3: r (abs) = 1.934e-05 (tol = 1.000e-10) r (rel) = 6.342e-04 (tol = 1.000e-10)

Newton iteration 4: r (abs) = 4.702e-08 (tol = 1.000e-10) r (rel) = 1.542e-06 (tol = 1.000e-10)

Newton iteration 5: r (abs) = 3.160e-13 (tol = 1.000e-10) r (rel) = 1.036e-11 (tol = 1.000e-10)

Newton solver finished in 5 iterations and 5 linear solver iterations.

Time: 6.825

Solving nonlinear variational problem.

Newton iteration 0: r (abs) = 3.174e-02 (tol = 1.000e-10) r (rel) = 1.000e+00 (tol = 1.000e-10)

Newton iteration 1: r (abs) = 8.671e-03 (tol = 1.000e-10) r (rel) = 2.732e-01 (tol = 1.000e-10)

Newton iteration 2: r (abs) = 1.097e-03 (tol = 1.000e-10) r (rel) = 3.456e-02 (tol = 1.000e-10)

Newton iteration 3: r (abs) = 2.184e-04 (tol = 1.000e-10) r (rel) = 6.883e-03 (tol = 1.000e-10)

Newton iteration 4: r (abs) = 9.723e-03 (tol = 1.000e-10) r (rel) = 3.064e-01 (tol = 1.000e-10)

Newton iteration 5: r (abs) = 2.484e-03 (tol = 1.000e-10) r (rel) = 7.826e-02 (tol = 1.000e-10)

Newton iteration 6: r (abs) = 6.795e-04 (tol = 1.000e-10) r (rel) = 2.141e-02 (tol = 1.000e-10)

Newton iteration 7: r (abs) = 2.479e-04 (tol = 1.000e-10) r (rel) = 7.809e-03 (tol = 1.000e-10)

Newton iteration 8: r (abs) = 4.503e-04 (tol = 1.000e-10) r (rel) = 1.419e-02 (tol = 1.000e-10)

Newton iteration 9: r (abs) = 2.143e-04 (tol = 1.000e-10) r (rel) = 6.752e-03 (tol = 1.000e-10)

Newton iteration 10: r (abs) = 2.080e-02 (tol = 1.000e-10) r (rel) = 6.552e-01 (tol = 1.000e-10)

Newton iteration 11: r (abs) = 5.253e-03 (tol = 1.000e-10) r (rel) = 1.655e-01 (tol = 1.000e-10)

Newton iteration 12: r (abs) = 1.369e-03 (tol = 1.000e-10) r (rel) = 4.313e-02 (tol = 1.000e-10)

Newton iteration 13: r (abs) = 4.055e-04 (tol = 1.000e-10) r (rel) = 1.278e-02 (tol = 1.000e-10)

Newton iteration 14: r (abs) = 2.144e-04 (tol = 1.000e-10) r (rel) = 6.754e-03 (tol = 1.000e-10)

Newton iteration 15: r (abs) = 1.833e-02 (tol = 1.000e-10) r (rel) = 5.775e-01 (tol = 1.000e-10)

Newton iteration 16: r (abs) = 4.636e-03 (tol = 1.000e-10) r (rel) = 1.461e-01 (tol = 1.000e-10)

Newton iteration 17: r (abs) = 1.215e-03 (tol = 1.000e-10) r (rel) = 3.828e-02 (tol = 1.000e-10)

Newton iteration 18: r (abs) = 3.686e-04 (tol = 1.000e-10) r (rel) = 1.161e-02 (tol = 1.000e-10)

Newton iteration 19: r (abs) = 2.193e-04 (tol = 1.000e-10) r (rel) = 6.911e-03 (tol = 1.000e-10)

Newton iteration 20: r (abs) = 2.149e-03 (tol = 1.000e-10) r (rel) = 6.770e-02 (tol = 1.000e-10)

Newton iteration 21: r (abs) = 5.965e-04 (tol = 1.000e-10) r (rel) = 1.879e-02 (tol = 1.000e-10)

Newton iteration 22: r (abs) = 2.324e-04 (tol = 1.000e-10) r (rel) = 7.322e-03 (tol = 1.000e-10)

Newton iteration 23: r (abs) = 7.238e-04 (tol = 1.000e-10) r (rel) = 2.281e-02 (tol = 1.000e-10)

Newton iteration 24: r (abs) = 2.568e-04 (tol = 1.000e-10) r (rel) = 8.091e-03 (tol = 1.000e-10)

Newton iteration 25: r (abs) = 3.830e-04 (tol = 1.000e-10) r (rel) = 1.207e-02 (tol = 1.000e-10)

Newton iteration 26: r (abs) = 2.167e-04 (tol = 1.000e-10) r (rel) = 6.827e-03 (tol = 1.000e-10)

Newton iteration 27: r (abs) = 4.022e-03 (tol = 1.000e-10) r (rel) = 1.267e-01 (tol = 1.000e-10)

Newton iteration 28: r (abs) = 1.062e-03 (tol = 1.000e-10) r (rel) = 3.346e-02 (tol = 1.000e-10)

Newton iteration 29: r (abs) = 3.324e-04 (tol = 1.000e-10) r (rel) = 1.047e-02 (tol = 1.000e-10)

Newton iteration 30: r (abs) = 2.328e-04 (tol = 1.000e-10) r (rel) = 7.335e-03 (tol = 1.000e-10)

Newton iteration 31: r (abs) = 7.113e-04 (tol = 1.000e-10) r (rel) = 2.241e-02 (tol = 1.000e-10)

Newton iteration 32: r (abs) = 2.542e-04 (tol = 1.000e-10) r (rel) = 8.010e-03 (tol = 1.000e-10)

Newton iteration 33: r (abs) = 3.992e-04 (tol = 1.000e-10) r (rel) = 1.258e-02 (tol = 1.000e-10)

Newton iteration 34: r (abs) = 2.148e-04 (tol = 1.000e-10) r (rel) = 6.769e-03 (tol = 1.000e-10)

Newton iteration 35: r (abs) = 1.068e-02 (tol = 1.000e-10) r (rel) = 3.367e-01 (tol = 1.000e-10)

Newton iteration 36: r (abs) = 2.726e-03 (tol = 1.000e-10) r (rel) = 8.588e-02 (tol = 1.000e-10)

Newton iteration 37: r (abs) = 7.394e-04 (tol = 1.000e-10) r (rel) = 2.330e-02 (tol = 1.000e-10)

Newton iteration 38: r (abs) = 2.600e-04 (tol = 1.000e-10) r (rel) = 8.193e-03 (tol = 1.000e-10)

Newton iteration 39: r (abs) = 3.653e-04 (tol = 1.000e-10) r (rel) = 1.151e-02 (tol = 1.000e-10)

Newton iteration 40: r (abs) = 2.201e-04 (tol = 1.000e-10) r (rel) = 6.936e-03 (tol = 1.000e-10)

Newton iteration 41: r (abs) = 1.896e-03 (tol = 1.000e-10) r (rel) = 5.975e-02 (tol = 1.000e-10)

Newton iteration 42: r (abs) = 5.343e-04 (tol = 1.000e-10) r (rel) = 1.683e-02 (tol = 1.000e-10)

Newton iteration 43: r (abs) = 2.226e-04 (tol = 1.000e-10) r (rel) = 7.015e-03 (tol = 1.000e-10)

Newton iteration 44: r (abs) = 1.393e-03 (tol = 1.000e-10) r (rel) = 4.390e-02 (tol = 1.000e-10)

Newton iteration 45: r (abs) = 4.114e-04 (tol = 1.000e-10) r (rel) = 1.296e-02 (tol = 1.000e-10)

Newton iteration 46: r (abs) = 2.141e-04 (tol = 1.000e-10) r (rel) = 6.745e-03 (tol = 1.000e-10)

Newton iteration 47: r (abs) = 3.521e-02 (tol = 1.000e-10) r (rel) = 1.109e+00 (tol = 1.000e-10)

Newton iteration 48: r (abs) = 8.856e-03 (tol = 1.000e-10) r (rel) = 2.791e-01 (tol = 1.000e-10)

Newton iteration 49: r (abs) = 2.269e-03 (tol = 1.000e-10) r (rel) = 7.149e-02 (tol = 1.000e-10)

Newton iteration 50: r (abs) = 6.262e-04 (tol = 1.000e-10) r (rel) = 1.973e-02 (tol = 1.000e-10)

terminate called after throwing an instance of 'std::runtime_error'

what():

*** -------------------------------------------------------------------------

*** DOLFIN encountered an error. If you are not able to resolve this issue

*** using the information listed below, you can ask for help at

***

*** fenics-support@googlegroups.com

***

*** Remember to include the error message listed below and, if possible,

*** include a *minimal* running example to reproduce the error.

***

*** -------------------------------------------------------------------------

*** Error: Unable to solve nonlinear system with NewtonSolver.

*** Reason: Newton solver did not converge because maximum number of iterations reached.

*** Where: This error was encountered inside NewtonSolver.cpp.

*** Process: 0

***

*** DOLFIN version: 2017.2.0

*** Git changeset: unknown

*** -------------------------------------------------------------------------

Community: FEniCS Project

### 1 Answer

0

Try using mumps. Also post complete error that your are facing. Are you using default solver ?
I use the default solver. Trying to set mums as a linear solver did not help unfortunately. I added the error into my question.

written
29 days ago by
fehere

Please login to add an answer/comment or follow this question.

The thread is closed. No new answer/comment may be added.