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.

import numpy as np
import ufl

from petsc4py import PETSc
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.VectorFunctionSpace(domain, ("CG", 2))

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=PETSc.ScalarType)

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, PETSc.ScalarType((0, 0, 0)))
T = fem.Constant(domain, PETSc.ScalarType((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 = PETSc.ScalarType(1.0e4)
nu = PETSc.ScalarType(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 = fem.petsc.NonlinearProblem(F, u, bcs)

and then create and customize the Newton solver

from dolfinx import nls
solver = nls.petsc.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.

import pyvista
pyvista.set_jupyter_backend("ipygany")

grid = pyvista.UnstructuredGrid(*plot.create_vtk_mesh(domain, domain.topology.dim))

def plot_function(t, uh):
    """
    Create a figure of the concentration uh warped visualized in 3D at timet step t.
    """
    p = pyvista.Plotter()
    # Create grid defined by the function space for visualization of the function
    topology, cells, geometry = plot.create_vtk_mesh(uh.function_space)
    function_grid = pyvista.UnstructuredGrid(topology, cells, geometry)
    var_name = f"u({t})"
    values = np.zeros((geometry.shape[0], 3))
    values[:, :len(uh)] = uh.x.array.reshape(geometry.shape[0], len(uh))
    function_grid[var_name] = values
    function_grid.set_active_vectors(var_name)
    # Warp mesh by deformation
    warped = function_grid.warp_by_vector(var_name, factor=1)
    
    # Add mesh to plotter and visualize
    actor = p.add_mesh(warped)
    p.show_axes()
    if not pyvista.OFF_SCREEN:
       p.show()
    else:
        pyvista.start_xvfb()
        figure_as_array = p.screenshot(f"diffusion_{t:.2f}.png")
        # Clear plotter for next plot
        p.remove_actor(actor)

plot_function(0, u)

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

from dolfinx import log
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}")
    plot_function(n, u)
2022-07-04 08:25:25.480 (   1.134s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:26.170 (   1.823s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:26.511 (   2.165s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 2: r (abs) = 22.2455 (tol = 1e-08) r (rel) = 0.134278(tol = 1e-08)
2022-07-04 08:25:26.736 (   2.389s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:27.080 (   2.734s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 3: r (abs) = 2.43261 (tol = 1e-08) r (rel) = 0.0146837(tol = 1e-08)
2022-07-04 08:25:27.305 (   2.958s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:27.650 (   3.304s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 4: r (abs) = 4.43158 (tol = 1e-08) r (rel) = 0.0267498(tol = 1e-08)
2022-07-04 08:25:27.875 (   3.528s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:28.218 (   3.872s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 5: r (abs) = 0.144189 (tol = 1e-08) r (rel) = 0.000870353(tol = 1e-08)
2022-07-04 08:25:28.443 (   4.097s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:28.789 (   4.442s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 6: r (abs) = 0.021424 (tol = 1e-08) r (rel) = 0.000129319(tol = 1e-08)
Time step 1, Number of iterations 8, Load [ 0.   0.  -1.5]
2022-07-04 08:25:29.013 (   4.667s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:29.355 (   5.008s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 7: r (abs) = 4.80067e-06 (tol = 1e-08) r (rel) = 2.89777e-08(tol = 1e-08)
2022-07-04 08:25:29.579 (   5.232s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:29.920 (   5.574s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 8: r (abs) = 2.64464e-11 (tol = 1e-08) r (rel) = 1.59635e-13(tol = 1e-08)
2022-07-04 08:25:29.920 (   5.574s) [main            ]       NewtonSolver.cpp:255   INFO| Newton solver finished in 8 iterations and 8 linear solver iterations.
2022-07-04 08:25:30.218 (   5.872s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:30.783 (   6.437s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:31.131 (   6.784s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 2: r (abs) = 17.3254 (tol = 1e-08) r (rel) = 0.117842(tol = 1e-08)
2022-07-04 08:25:31.355 (   7.008s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:31.701 (   7.355s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 3: r (abs) = 5.14882 (tol = 1e-08) r (rel) = 0.0350207(tol = 1e-08)
2022-07-04 08:25:31.926 (   7.579s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:32.273 (   7.926s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 4: r (abs) = 7.24003 (tol = 1e-08) r (rel) = 0.0492445(tol = 1e-08)
2022-07-04 08:25:32.497 (   8.151s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:32.841 (   8.495s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 5: r (abs) = 0.777889 (tol = 1e-08) r (rel) = 0.00529096(tol = 1e-08)
2022-07-04 08:25:33.065 (   8.719s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:33.407 (   9.061s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 6: r (abs) = 1.25525 (tol = 1e-08) r (rel) = 0.00853785(tol = 1e-08)
2022-07-04 08:25:33.631 (   9.285s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:33.976 (   9.630s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 7: r (abs) = 0.00849512 (tol = 1e-08) r (rel) = 5.77813e-05(tol = 1e-08)
2022-07-04 08:25:34.200 (   9.854s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:34.546 (  10.199s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 8: r (abs) = 0.000192107 (tol = 1e-08) r (rel) = 1.30665e-06(tol = 1e-08)
Time step 2, Number of iterations 9, Load [ 0.  0. -3.]
2022-07-04 08:25:34.772 (  10.425s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:35.115 (  10.769s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 9: r (abs) = 1.70757e-10 (tol = 1e-08) r (rel) = 1.16144e-12(tol = 1e-08)
2022-07-04 08:25:35.115 (  10.769s) [main            ]       NewtonSolver.cpp:255   INFO| Newton solver finished in 9 iterations and 9 linear solver iterations.
2022-07-04 08:25:35.406 (  11.060s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:35.973 (  11.627s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:36.318 (  11.972s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 2: r (abs) = 10.0011 (tol = 1e-08) r (rel) = 0.0887471(tol = 1e-08)
2022-07-04 08:25:36.544 (  12.198s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:36.890 (  12.544s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 3: r (abs) = 5.33026 (tol = 1e-08) r (rel) = 0.0472992(tol = 1e-08)
2022-07-04 08:25:37.114 (  12.768s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:37.459 (  13.112s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 4: r (abs) = 11.9901 (tol = 1e-08) r (rel) = 0.106397(tol = 1e-08)
2022-07-04 08:25:37.683 (  13.336s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:38.027 (  13.681s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 5: r (abs) = 2.29702 (tol = 1e-08) r (rel) = 0.0203831(tol = 1e-08)
2022-07-04 08:25:38.252 (  13.905s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:38.596 (  14.250s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 6: r (abs) = 3.90234 (tol = 1e-08) r (rel) = 0.0346282(tol = 1e-08)
2022-07-04 08:25:38.821 (  14.474s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:39.164 (  14.818s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 7: r (abs) = 0.236535 (tol = 1e-08) r (rel) = 0.00209895(tol = 1e-08)
2022-07-04 08:25:39.389 (  15.042s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:39.736 (  15.390s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 8: r (abs) = 0.0427142 (tol = 1e-08) r (rel) = 0.000379034(tol = 1e-08)
2022-07-04 08:25:39.961 (  15.614s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
Time step 3, Number of iterations 10, Load [ 0.   0.  -4.5]
2022-07-04 08:25:40.308 (  15.962s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 9: r (abs) = 2.87798e-05 (tol = 1e-08) r (rel) = 2.55384e-07(tol = 1e-08)
2022-07-04 08:25:40.532 (  16.186s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:40.877 (  16.531s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 10: r (abs) = 6.09481e-10 (tol = 1e-08) r (rel) = 5.40836e-12(tol = 1e-08)
2022-07-04 08:25:40.877 (  16.531s) [main            ]       NewtonSolver.cpp:255   INFO| Newton solver finished in 10 iterations and 10 linear solver iterations.
2022-07-04 08:25:41.168 (  16.821s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:41.734 (  17.388s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:42.077 (  17.731s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 2: r (abs) = 5.50693 (tol = 1e-08) r (rel) = 0.0653918(tol = 1e-08)
2022-07-04 08:25:42.301 (  17.954s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:42.644 (  18.298s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 3: r (abs) = 26.2489 (tol = 1e-08) r (rel) = 0.311692(tol = 1e-08)
2022-07-04 08:25:42.869 (  18.523s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:43.214 (  18.868s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 4: r (abs) = 2.30927 (tol = 1e-08) r (rel) = 0.0274213(tol = 1e-08)
2022-07-04 08:25:43.438 (  19.092s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:43.782 (  19.436s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 5: r (abs) = 14.0562 (tol = 1e-08) r (rel) = 0.16691(tol = 1e-08)
2022-07-04 08:25:44.006 (  19.660s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:44.349 (  20.003s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 6: r (abs) = 0.222774 (tol = 1e-08) r (rel) = 0.00264532(tol = 1e-08)
2022-07-04 08:25:44.574 (  20.228s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:44.919 (  20.573s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 7: r (abs) = 0.286671 (tol = 1e-08) r (rel) = 0.00340406(tol = 1e-08)
2022-07-04 08:25:45.143 (  20.797s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
Time step 4, Number of iterations 9, Load [ 0.  0. -6.]
2022-07-04 08:25:45.489 (  21.143s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 8: r (abs) = 0.000321869 (tol = 1e-08) r (rel) = 3.82203e-06(tol = 1e-08)
2022-07-04 08:25:45.713 (  21.367s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:46.056 (  21.709s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 9: r (abs) = 2.63796e-07 (tol = 1e-08) r (rel) = 3.13244e-09(tol = 1e-08)
2022-07-04 08:25:46.056 (  21.709s) [main            ]       NewtonSolver.cpp:255   INFO| Newton solver finished in 9 iterations and 9 linear solver iterations.
2022-07-04 08:25:46.345 (  21.999s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:46.913 (  22.566s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:47.256 (  22.910s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 2: r (abs) = 3.19462 (tol = 1e-08) r (rel) = 0.0496479(tol = 1e-08)
2022-07-04 08:25:47.480 (  23.133s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:47.826 (  23.480s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 3: r (abs) = 7.71429 (tol = 1e-08) r (rel) = 0.119888(tol = 1e-08)
2022-07-04 08:25:48.051 (  23.704s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:48.395 (  24.048s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 4: r (abs) = 0.850873 (tol = 1e-08) r (rel) = 0.0132235(tol = 1e-08)
2022-07-04 08:25:48.619 (  24.273s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:48.965 (  24.618s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 5: r (abs) = 0.371434 (tol = 1e-08) r (rel) = 0.0057725(tol = 1e-08)
2022-07-04 08:25:49.189 (  24.843s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
Time step 5, Number of iterations 8, Load [ 0.   0.  -7.5]
2022-07-04 08:25:49.535 (  25.189s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 6: r (abs) = 0.00215066 (tol = 1e-08) r (rel) = 3.34236e-05(tol = 1e-08)
2022-07-04 08:25:49.759 (  25.413s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:50.103 (  25.757s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 7: r (abs) = 2.54607e-06 (tol = 1e-08) r (rel) = 3.95687e-08(tol = 1e-08)
2022-07-04 08:25:50.330 (  25.983s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:50.675 (  26.328s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 8: r (abs) = 1.51485e-13 (tol = 1e-08) r (rel) = 2.35425e-15(tol = 1e-08)
2022-07-04 08:25:50.675 (  26.328s) [main            ]       NewtonSolver.cpp:255   INFO| Newton solver finished in 8 iterations and 8 linear solver iterations.
2022-07-04 08:25:50.965 (  26.619s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:51.533 (  27.187s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:51.875 (  27.529s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 2: r (abs) = 2.00649 (tol = 1e-08) r (rel) = 0.0395622(tol = 1e-08)
2022-07-04 08:25:52.100 (  27.754s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:52.443 (  28.097s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 3: r (abs) = 4.60977 (tol = 1e-08) r (rel) = 0.0908914(tol = 1e-08)
2022-07-04 08:25:52.668 (  28.322s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:53.011 (  28.664s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 4: r (abs) = 0.185372 (tol = 1e-08) r (rel) = 0.00365501(tol = 1e-08)
2022-07-04 08:25:53.235 (  28.889s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:53.579 (  29.232s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 5: r (abs) = 0.024688 (tol = 1e-08) r (rel) = 0.000486777(tol = 1e-08)
2022-07-04 08:25:53.804 (  29.458s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:54.150 (  29.804s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 6: r (abs) = 5.69254e-06 (tol = 1e-08) r (rel) = 1.12241e-07(tol = 1e-08)
Time step 6, Number of iterations 7, Load [ 0.  0. -9.]
2022-07-04 08:25:54.375 (  30.028s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:54.721 (  30.374s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 7: r (abs) = 2.67626e-11 (tol = 1e-08) r (rel) = 5.27682e-13(tol = 1e-08)
2022-07-04 08:25:54.721 (  30.375s) [main            ]       NewtonSolver.cpp:255   INFO| Newton solver finished in 7 iterations and 7 linear solver iterations.
2022-07-04 08:25:55.013 (  30.667s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:55.581 (  31.234s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:55.926 (  31.580s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 2: r (abs) = 1.38506 (tol = 1e-08) r (rel) = 0.0336622(tol = 1e-08)
2022-07-04 08:25:56.150 (  31.804s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:56.494 (  32.147s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 3: r (abs) = 3.03739 (tol = 1e-08) r (rel) = 0.07382(tol = 1e-08)
2022-07-04 08:25:56.717 (  32.371s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:57.059 (  32.712s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 4: r (abs) = 0.0412386 (tol = 1e-08) r (rel) = 0.00100225(tol = 1e-08)
2022-07-04 08:25:57.284 (  32.937s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:57.627 (  33.281s) [main            ]       NewtonSolver.cpp:36    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]
2022-07-04 08:25:57.852 (  33.505s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:58.200 (  33.854s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 6: r (abs) = 1.78865e-08 (tol = 1e-08) r (rel) = 4.34709e-10(tol = 1e-08)
2022-07-04 08:25:58.200 (  33.854s) [main            ]       NewtonSolver.cpp:255   INFO| Newton solver finished in 6 iterations and 6 linear solver iterations.
2022-07-04 08:25:58.490 (  34.143s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:59.056 (  34.710s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:59.402 (  35.055s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 2: r (abs) = 1.06336 (tol = 1e-08) r (rel) = 0.031085(tol = 1e-08)
2022-07-04 08:25:59.626 (  35.279s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:25:59.971 (  35.624s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 3: r (abs) = 2.0477 (tol = 1e-08) r (rel) = 0.0598598(tol = 1e-08)
2022-07-04 08:26:00.195 (  35.849s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:26:00.540 (  36.194s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 4: r (abs) = 0.00897719 (tol = 1e-08) r (rel) = 0.000262427(tol = 1e-08)
2022-07-04 08:26:00.765 (  36.419s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
Time step 8, Number of iterations 6, Load [  0.   0. -12.]
2022-07-04 08:26:01.109 (  36.762s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 5: r (abs) = 0.000167422 (tol = 1e-08) r (rel) = 4.89419e-06(tol = 1e-08)
2022-07-04 08:26:01.333 (  36.986s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:26:01.677 (  37.330s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 6: r (abs) = 3.23994e-11 (tol = 1e-08) r (rel) = 9.4712e-13(tol = 1e-08)
2022-07-04 08:26:01.677 (  37.330s) [main            ]       NewtonSolver.cpp:255   INFO| Newton solver finished in 6 iterations and 6 linear solver iterations.
2022-07-04 08:26:01.969 (  37.623s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:26:02.539 (  38.193s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:26:02.884 (  38.538s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 2: r (abs) = 0.898789 (tol = 1e-08) r (rel) = 0.0309666(tol = 1e-08)
2022-07-04 08:26:03.109 (  38.763s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:26:03.456 (  39.109s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 3: r (abs) = 1.38354 (tol = 1e-08) r (rel) = 0.0476679(tol = 1e-08)
2022-07-04 08:26:03.680 (  39.334s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:26:04.024 (  39.677s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 4: r (abs) = 0.00185096 (tol = 1e-08) r (rel) = 6.37724e-05(tol = 1e-08)
2022-07-04 08:26:04.248 (  39.902s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:26:04.595 (  40.248s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 5: r (abs) = 7.87183e-06 (tol = 1e-08) r (rel) = 2.71213e-07(tol = 1e-08)
Time step 9, Number of iterations 6, Load [  0.    0.  -13.5]
2022-07-04 08:26:04.819 (  40.473s) [main            ]              petsc.cpp:677   INFO| PETSc Krylov solver starting to solve system.
2022-07-04 08:26:05.163 (  40.817s) [main            ]       NewtonSolver.cpp:36    INFO| Newton iteration 6: r (abs) = 2.75782e-13 (tol = 1e-08) r (rel) = 9.50169e-15(tol = 1e-08)
2022-07-04 08:26:05.163 (  40.817s) [main            ]       NewtonSolver.cpp:255   INFO| Newton solver finished in 6 iterations and 6 linear solver iterations.