Hyperelasticity

Hyperelasticity#

Author: Jørgen S. Dokken and Garth N. Wells

This section shows how to solve the hyperelasticity problem for deformation of a beam.

We will also show how to create a constant boundary condition for a vector function space.

We start by importing DOLFINx and some additional dependencies. Then, we create a slender cantilever consisting of hexahedral elements and create the function space V for our unknown.

from dolfinx import log, default_scalar_type
from dolfinx.fem.petsc import NonlinearProblem
from dolfinx.nls.petsc import NewtonSolver
import pyvista
import numpy as np
import ufl

from mpi4py import MPI
from dolfinx import fem, mesh, plot
L = 20.0
domain = mesh.create_box(MPI.COMM_WORLD, [[0.0, 0.0, 0.0], [L, 1, 1]], [20, 5, 5], mesh.CellType.hexahedron)
V = fem.functionspace(domain, ("Lagrange", 2, (domain.geometry.dim, )))

We create two python functions for determining the facets to apply boundary conditions to

def left(x):
    return np.isclose(x[0], 0)


def right(x):
    return np.isclose(x[0], L)


fdim = domain.topology.dim - 1
left_facets = mesh.locate_entities_boundary(domain, fdim, left)
right_facets = mesh.locate_entities_boundary(domain, fdim, right)

Next, we create a marker based on these two functions

# Concatenate and sort the arrays based on facet indices. Left facets marked with 1, right facets with two
marked_facets = np.hstack([left_facets, right_facets])
marked_values = np.hstack([np.full_like(left_facets, 1), np.full_like(right_facets, 2)])
sorted_facets = np.argsort(marked_facets)
facet_tag = mesh.meshtags(domain, fdim, marked_facets[sorted_facets], marked_values[sorted_facets])

We then create a function for supplying the boundary condition on the left side, which is fixed.

u_bc = np.array((0,) * domain.geometry.dim, dtype=default_scalar_type)

To apply the boundary condition, we identity the dofs located on the facets marked by the MeshTag.

left_dofs = fem.locate_dofs_topological(V, facet_tag.dim, facet_tag.find(1))
bcs = [fem.dirichletbc(u_bc, left_dofs, V)]

Next, we define the body force on the reference configuration (B), and nominal (first Piola-Kirchhoff) traction (T).

B = fem.Constant(domain, default_scalar_type((0, 0, 0)))
T = fem.Constant(domain, default_scalar_type((0, 0, 0)))

Define the test and solution functions on the space \(V\)

v = ufl.TestFunction(V)
u = fem.Function(V)

Define kinematic quantities used in the problem

# Spatial dimension
d = len(u)

# Identity tensor
I = ufl.variable(ufl.Identity(d))

# Deformation gradient
F = ufl.variable(I + ufl.grad(u))

# Right Cauchy-Green tensor
C = ufl.variable(F.T * F)

# Invariants of deformation tensors
Ic = ufl.variable(ufl.tr(C))
J = ufl.variable(ufl.det(F))

Define the elasticity model via a stored strain energy density function \(\psi\), and create the expression for the first Piola-Kirchhoff stress:

# Elasticity parameters
E = default_scalar_type(1.0e4)
nu = default_scalar_type(0.3)
mu = fem.Constant(domain, E / (2 * (1 + nu)))
lmbda = fem.Constant(domain, E * nu / ((1 + nu) * (1 - 2 * nu)))
# Stored strain energy density (compressible neo-Hookean model)
psi = (mu / 2) * (Ic - 3) - mu * ufl.ln(J) + (lmbda / 2) * (ufl.ln(J))**2
# Stress
# Hyper-elasticity
P = ufl.diff(psi, F)

Comparison to linear elasticity

To illustrate the difference between linear and hyperelasticity, the following lines can be uncommented to solve the linear elasticity problem.

# P = 2.0 * mu * ufl.sym(ufl.grad(u)) + lmbda * ufl.tr(ufl.sym(ufl.grad(u))) * I

Define the variational form with traction integral over all facets with value 2. We set the quadrature degree for the integrals to 4.

metadata = {"quadrature_degree": 4}
ds = ufl.Measure('ds', domain=domain, subdomain_data=facet_tag, metadata=metadata)
dx = ufl.Measure("dx", domain=domain, metadata=metadata)
# Define form F (we want to find u such that F(u) = 0)
F = ufl.inner(ufl.grad(v), P) * dx - ufl.inner(v, B) * dx - ufl.inner(v, T) * ds(2)

As the varitional form is non-linear and written on residual form, we use the non-linear problem class from DOLFINx to set up required structures to use a Newton solver.

problem = NonlinearProblem(F, u, bcs)

and then create and customize the Newton solver

solver = NewtonSolver(domain.comm, problem)

# Set Newton solver options
solver.atol = 1e-8
solver.rtol = 1e-8
solver.convergence_criterion = "incremental"

We create a function to plot the solution at each time step.

pyvista.start_xvfb()
plotter = pyvista.Plotter()
plotter.open_gif("deformation.gif", fps=3)

topology, cells, geometry = plot.vtk_mesh(u.function_space)
function_grid = pyvista.UnstructuredGrid(topology, cells, geometry)

values = np.zeros((geometry.shape[0], 3))
values[:, :len(u)] = u.x.array.reshape(geometry.shape[0], len(u))
function_grid["u"] = values
function_grid.set_active_vectors("u")

# Warp mesh by deformation
warped = function_grid.warp_by_vector("u", factor=1)
warped.set_active_vectors("u")

# Add mesh to plotter and visualize
actor = plotter.add_mesh(warped, show_edges=True, lighting=False, clim=[0, 10])

# Compute magnitude of displacement to visualize in GIF
Vs = fem.functionspace(domain, ("Lagrange", 2))
magnitude = fem.Function(Vs)
us = fem.Expression(ufl.sqrt(sum([u[i]**2 for i in range(len(u))])), Vs.element.interpolation_points())
magnitude.interpolate(us)
warped["mag"] = magnitude.x.array

Finally, we solve the problem over several time steps, updating the z-component of the traction

log.set_log_level(log.LogLevel.INFO)
tval0 = -1.5
for n in range(1, 10):
    T.value[2] = n * tval0
    num_its, converged = solver.solve(u)
    assert (converged)
    u.x.scatter_forward()
    print(f"Time step {n}, Number of iterations {num_its}, Load {T.value}")
    function_grid["u"][:, :len(u)] = u.x.array.reshape(geometry.shape[0], len(u))
    magnitude.interpolate(us)
    warped.set_active_scalars("mag")
    warped_n = function_grid.warp_by_vector(factor=1)
    warped.points[:, :] = warped_n.points
    warped.point_data["mag"][:] = magnitude.x.array
    plotter.update_scalar_bar_range([0, 10])
    plotter.write_frame()
plotter.close()
2024-05-20 08:05:19.403 (   6.594s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:05:21.910 (   9.101s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:05:23.907 (  11.098s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 2: r (abs) = 22.2455 (tol = 1e-08) r (rel) = 0.134278(tol = 1e-08)
2024-05-20 08:05:24.239 (  11.429s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:05:26.221 (  13.412s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 3: r (abs) = 2.43261 (tol = 1e-08) r (rel) = 0.0146837(tol = 1e-08)
2024-05-20 08:05:26.556 (  13.746s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:05:28.547 (  15.738s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 4: r (abs) = 4.43158 (tol = 1e-08) r (rel) = 0.0267498(tol = 1e-08)
2024-05-20 08:05:28.872 (  16.063s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:05:30.855 (  18.046s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 5: r (abs) = 0.144189 (tol = 1e-08) r (rel) = 0.000870353(tol = 1e-08)
2024-05-20 08:05:31.189 (  18.380s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:05:33.178 (  20.368s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 6: r (abs) = 0.021424 (tol = 1e-08) r (rel) = 0.000129319(tol = 1e-08)
2024-05-20 08:05:33.502 (  20.693s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:05:35.495 (  22.686s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 7: r (abs) = 4.80065e-06 (tol = 1e-08) r (rel) = 2.89776e-08(tol = 1e-08)
2024-05-20 08:05:35.827 (  23.017s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:05:37.820 (  25.010s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 8: r (abs) = 2.65767e-11 (tol = 1e-08) r (rel) = 1.60422e-13(tol = 1e-08)
2024-05-20 08:05:37.820 (  25.010s) [main            ]       NewtonSolver.cpp:252   INFO| Newton solver finished in 8 iterations and 8 linear solver iterations.
Time step 1, Number of iterations 8, Load [ 0.   0.  -1.5]
2024-05-20 08:05:38.340 (  25.531s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:05:40.657 (  27.848s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:05:42.640 (  29.831s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 2: r (abs) = 17.3254 (tol = 1e-08) r (rel) = 0.117842(tol = 1e-08)
2024-05-20 08:05:42.965 (  30.156s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:05:44.953 (  32.144s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 3: r (abs) = 5.14882 (tol = 1e-08) r (rel) = 0.0350207(tol = 1e-08)
2024-05-20 08:05:45.285 (  32.476s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:05:47.268 (  34.458s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 4: r (abs) = 7.24003 (tol = 1e-08) r (rel) = 0.0492445(tol = 1e-08)
2024-05-20 08:05:47.593 (  34.784s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:05:49.579 (  36.770s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 5: r (abs) = 0.777889 (tol = 1e-08) r (rel) = 0.00529096(tol = 1e-08)
2024-05-20 08:05:49.917 (  37.108s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:05:51.904 (  39.095s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 6: r (abs) = 1.25525 (tol = 1e-08) r (rel) = 0.00853785(tol = 1e-08)
2024-05-20 08:05:52.231 (  39.422s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:05:54.217 (  41.408s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 7: r (abs) = 0.00849512 (tol = 1e-08) r (rel) = 5.77812e-05(tol = 1e-08)
2024-05-20 08:05:54.553 (  41.743s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:05:56.543 (  43.733s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 8: r (abs) = 0.000192107 (tol = 1e-08) r (rel) = 1.30665e-06(tol = 1e-08)
2024-05-20 08:05:56.878 (  44.069s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:05:58.868 (  46.058s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 9: r (abs) = 1.70902e-10 (tol = 1e-08) r (rel) = 1.16242e-12(tol = 1e-08)
2024-05-20 08:05:58.868 (  46.058s) [main            ]       NewtonSolver.cpp:252   INFO| Newton solver finished in 9 iterations and 9 linear solver iterations.
Time step 2, Number of iterations 9, Load [ 0.  0. -3.]
2024-05-20 08:05:59.255 (  46.446s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:01.583 (  48.773s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:03.570 (  50.760s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 2: r (abs) = 10.0011 (tol = 1e-08) r (rel) = 0.0887471(tol = 1e-08)
2024-05-20 08:06:03.904 (  51.095s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:05.894 (  53.085s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 3: r (abs) = 5.33026 (tol = 1e-08) r (rel) = 0.0472992(tol = 1e-08)
2024-05-20 08:06:06.225 (  53.416s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:08.220 (  55.411s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 4: r (abs) = 11.9901 (tol = 1e-08) r (rel) = 0.106397(tol = 1e-08)
2024-05-20 08:06:08.552 (  55.743s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:10.567 (  57.757s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 5: r (abs) = 2.29702 (tol = 1e-08) r (rel) = 0.0203831(tol = 1e-08)
2024-05-20 08:06:10.895 (  58.086s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:12.884 (  60.074s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 6: r (abs) = 3.90234 (tol = 1e-08) r (rel) = 0.0346282(tol = 1e-08)
2024-05-20 08:06:13.219 (  60.409s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:15.206 (  62.397s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 7: r (abs) = 0.236535 (tol = 1e-08) r (rel) = 0.00209895(tol = 1e-08)
2024-05-20 08:06:15.532 (  62.722s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:17.514 (  64.704s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 8: r (abs) = 0.0427142 (tol = 1e-08) r (rel) = 0.000379034(tol = 1e-08)
2024-05-20 08:06:17.845 (  65.036s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:19.833 (  67.024s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 9: r (abs) = 2.87798e-05 (tol = 1e-08) r (rel) = 2.55384e-07(tol = 1e-08)
2024-05-20 08:06:20.162 (  67.353s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:22.149 (  69.340s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 10: r (abs) = 6.09146e-10 (tol = 1e-08) r (rel) = 5.40539e-12(tol = 1e-08)
2024-05-20 08:06:22.149 (  69.340s) [main            ]       NewtonSolver.cpp:252   INFO| Newton solver finished in 10 iterations and 10 linear solver iterations.
Time step 3, Number of iterations 10, Load [ 0.   0.  -4.5]
2024-05-20 08:06:22.534 (  69.725s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:24.857 (  72.048s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:26.842 (  74.033s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 2: r (abs) = 5.50693 (tol = 1e-08) r (rel) = 0.0653918(tol = 1e-08)
2024-05-20 08:06:27.176 (  74.367s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:29.163 (  76.354s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 3: r (abs) = 26.2489 (tol = 1e-08) r (rel) = 0.311692(tol = 1e-08)
2024-05-20 08:06:29.489 (  76.680s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:31.479 (  78.670s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 4: r (abs) = 2.30927 (tol = 1e-08) r (rel) = 0.0274213(tol = 1e-08)
2024-05-20 08:06:31.811 (  79.002s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:33.801 (  80.991s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 5: r (abs) = 14.0562 (tol = 1e-08) r (rel) = 0.16691(tol = 1e-08)
2024-05-20 08:06:34.126 (  81.317s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:36.112 (  83.303s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 6: r (abs) = 0.222774 (tol = 1e-08) r (rel) = 0.00264532(tol = 1e-08)
2024-05-20 08:06:36.447 (  83.638s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:38.434 (  85.625s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 7: r (abs) = 0.286671 (tol = 1e-08) r (rel) = 0.00340406(tol = 1e-08)
2024-05-20 08:06:38.769 (  85.960s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:40.757 (  87.948s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 8: r (abs) = 0.000321869 (tol = 1e-08) r (rel) = 3.82203e-06(tol = 1e-08)
2024-05-20 08:06:41.091 (  88.282s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:43.082 (  90.273s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 9: r (abs) = 2.63796e-07 (tol = 1e-08) r (rel) = 3.13244e-09(tol = 1e-08)
2024-05-20 08:06:43.082 (  90.273s) [main            ]       NewtonSolver.cpp:252   INFO| Newton solver finished in 9 iterations and 9 linear solver iterations.
Time step 4, Number of iterations 9, Load [ 0.  0. -6.]
2024-05-20 08:06:43.468 (  90.659s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:45.789 (  92.980s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:47.777 (  94.967s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 2: r (abs) = 3.19462 (tol = 1e-08) r (rel) = 0.0496479(tol = 1e-08)
2024-05-20 08:06:48.108 (  95.299s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:50.091 (  97.282s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 3: r (abs) = 7.71429 (tol = 1e-08) r (rel) = 0.119888(tol = 1e-08)
2024-05-20 08:06:50.426 (  97.617s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:52.406 (  99.597s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 4: r (abs) = 0.850873 (tol = 1e-08) r (rel) = 0.0132235(tol = 1e-08)
2024-05-20 08:06:52.743 (  99.934s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:54.732 ( 101.922s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 5: r (abs) = 0.371434 (tol = 1e-08) r (rel) = 0.0057725(tol = 1e-08)
2024-05-20 08:06:55.066 ( 102.257s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:57.049 ( 104.239s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 6: r (abs) = 0.00215066 (tol = 1e-08) r (rel) = 3.34236e-05(tol = 1e-08)
2024-05-20 08:06:57.383 ( 104.573s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:06:59.366 ( 106.557s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 7: r (abs) = 2.54607e-06 (tol = 1e-08) r (rel) = 3.95687e-08(tol = 1e-08)
2024-05-20 08:06:59.694 ( 106.885s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:01.673 ( 108.863s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 8: r (abs) = 1.73699e-13 (tol = 1e-08) r (rel) = 2.69948e-15(tol = 1e-08)
2024-05-20 08:07:01.673 ( 108.863s) [main            ]       NewtonSolver.cpp:252   INFO| Newton solver finished in 8 iterations and 8 linear solver iterations.
Time step 5, Number of iterations 8, Load [ 0.   0.  -7.5]
2024-05-20 08:07:02.070 ( 109.261s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:04.389 ( 111.579s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:06.374 ( 113.565s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 2: r (abs) = 2.00649 (tol = 1e-08) r (rel) = 0.0395622(tol = 1e-08)
2024-05-20 08:07:06.709 ( 113.900s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:08.690 ( 115.881s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 3: r (abs) = 4.60977 (tol = 1e-08) r (rel) = 0.0908914(tol = 1e-08)
2024-05-20 08:07:09.018 ( 116.208s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:11.011 ( 118.201s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 4: r (abs) = 0.185372 (tol = 1e-08) r (rel) = 0.00365501(tol = 1e-08)
2024-05-20 08:07:11.344 ( 118.534s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:13.329 ( 120.520s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 5: r (abs) = 0.024688 (tol = 1e-08) r (rel) = 0.000486777(tol = 1e-08)
2024-05-20 08:07:13.667 ( 120.857s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:15.644 ( 122.835s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 6: r (abs) = 5.69255e-06 (tol = 1e-08) r (rel) = 1.12241e-07(tol = 1e-08)
2024-05-20 08:07:15.979 ( 123.169s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:17.964 ( 125.155s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 7: r (abs) = 2.62006e-11 (tol = 1e-08) r (rel) = 5.166e-13(tol = 1e-08)
2024-05-20 08:07:17.964 ( 125.155s) [main            ]       NewtonSolver.cpp:252   INFO| Newton solver finished in 7 iterations and 7 linear solver iterations.
Time step 6, Number of iterations 7, Load [ 0.  0. -9.]
2024-05-20 08:07:18.360 ( 125.551s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:20.670 ( 127.861s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:22.659 ( 129.850s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 2: r (abs) = 1.38506 (tol = 1e-08) r (rel) = 0.0336622(tol = 1e-08)
2024-05-20 08:07:22.996 ( 130.186s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:24.982 ( 132.173s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 3: r (abs) = 3.03739 (tol = 1e-08) r (rel) = 0.07382(tol = 1e-08)
2024-05-20 08:07:25.315 ( 132.505s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:27.293 ( 134.484s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 4: r (abs) = 0.0412386 (tol = 1e-08) r (rel) = 0.00100225(tol = 1e-08)
2024-05-20 08:07:27.631 ( 134.821s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:29.616 ( 136.806s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 5: r (abs) = 0.00205057 (tol = 1e-08) r (rel) = 4.98364e-05(tol = 1e-08)
Time step 7, Number of iterations 6, Load [  0.    0.  -10.5]
2024-05-20 08:07:29.952 ( 137.143s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:31.934 ( 139.125s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 6: r (abs) = 1.78867e-08 (tol = 1e-08) r (rel) = 4.34714e-10(tol = 1e-08)
2024-05-20 08:07:31.934 ( 139.125s) [main            ]       NewtonSolver.cpp:252   INFO| Newton solver finished in 6 iterations and 6 linear solver iterations.
2024-05-20 08:07:32.319 ( 139.510s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:34.645 ( 141.836s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:36.628 ( 143.819s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 2: r (abs) = 1.06336 (tol = 1e-08) r (rel) = 0.031085(tol = 1e-08)
2024-05-20 08:07:36.964 ( 144.155s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:38.951 ( 146.142s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 3: r (abs) = 2.0477 (tol = 1e-08) r (rel) = 0.0598598(tol = 1e-08)
2024-05-20 08:07:39.286 ( 146.477s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:41.270 ( 148.461s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 4: r (abs) = 0.00897719 (tol = 1e-08) r (rel) = 0.000262427(tol = 1e-08)
2024-05-20 08:07:41.601 ( 148.792s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:43.587 ( 150.778s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 5: r (abs) = 0.000167422 (tol = 1e-08) r (rel) = 4.89419e-06(tol = 1e-08)
2024-05-20 08:07:43.923 ( 151.114s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:45.904 ( 153.095s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 6: r (abs) = 3.24488e-11 (tol = 1e-08) r (rel) = 9.48564e-13(tol = 1e-08)
2024-05-20 08:07:45.904 ( 153.095s) [main            ]       NewtonSolver.cpp:252   INFO| Newton solver finished in 6 iterations and 6 linear solver iterations.
Time step 8, Number of iterations 6, Load [  0.   0. -12.]
2024-05-20 08:07:46.299 ( 153.490s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:48.612 ( 155.802s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:50.597 ( 157.787s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 2: r (abs) = 0.898789 (tol = 1e-08) r (rel) = 0.0309666(tol = 1e-08)
2024-05-20 08:07:50.921 ( 158.111s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:52.906 ( 160.097s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 3: r (abs) = 1.38354 (tol = 1e-08) r (rel) = 0.0476679(tol = 1e-08)
2024-05-20 08:07:53.241 ( 160.432s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:55.230 ( 162.421s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 4: r (abs) = 0.00185096 (tol = 1e-08) r (rel) = 6.37724e-05(tol = 1e-08)
2024-05-20 08:07:55.555 ( 162.745s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:57.540 ( 164.731s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 5: r (abs) = 7.87183e-06 (tol = 1e-08) r (rel) = 2.71213e-07(tol = 1e-08)
2024-05-20 08:07:57.865 ( 165.055s) [main            ]              petsc.cpp:700   INFO| PETSc Krylov solver starting to solve system.
2024-05-20 08:07:59.850 ( 167.040s) [main            ]       NewtonSolver.cpp:38    INFO| Newton iteration 6: r (abs) = 2.05992e-13 (tol = 1e-08) r (rel) = 7.09718e-15(tol = 1e-08)
2024-05-20 08:07:59.850 ( 167.040s) [main            ]       NewtonSolver.cpp:252   INFO| Newton solver finished in 6 iterations and 6 linear solver iterations.
Time step 9, Number of iterations 6, Load [  0.    0.  -13.5]
gif