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.

import numpy as np
import ufl

from petsc4py import PETSc
from mpi4py import MPI
from dolfinx import fem
from dolfinx.mesh import CellType, MeshTags, create_box, locate_entities_boundary
L = 20.0
mesh = create_box(MPI.COMM_WORLD,[[0.0,0.0,0.0], [L, 1, 1]], [20, 5, 5], CellType.hexahedron)
V = fem.VectorFunctionSpace(mesh, ("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)

left_facets = locate_entities_boundary(mesh, mesh.topology.dim - 1, left)
right_facets = locate_entities_boundary(mesh, mesh.topology.dim - 1, 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(len(left_facets), 1, dtype=np.int32), np.full(len(right_facets), 2, dtype=np.int32)])
sorted_facets = np.argsort(marked_facets)
facet_tag = MeshTags(mesh, mesh.topology.dim-1, 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,) * mesh.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.indices[facet_tag.values==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(mesh, PETSc.ScalarType((0, 0, 0)))
T = fem.Constant(mesh, 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(mesh, E/(2*(1 + nu)))
lmbda = fem.Constant(mesh, 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', subdomain_data=facet_tag, metadata=metadata)
dx = ufl.Measure("dx", 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.NonlinearProblem(F, u, bcs)

and then create and customize the Newton solver

from dolfinx.nls import NewtonSolver
solver = NewtonSolver(MPI.COMM_WORLD, 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")

from dolfinx.plot import create_vtk_topology
topology, cell_types = create_vtk_topology(mesh, mesh.topology.dim)
grid = pyvista.UnstructuredGrid(topology, cell_types, mesh.geometry.x)

def plot_function(t, uh):
    """
    Create a figure of the concentration uh warped visualized in 3D at timet step t.
    """
    p = pyvista.Plotter()
   
    # Update point values on pyvista grid

    topology, cell_types = create_vtk_topology(V)
     # We create a geometry for our modified mesh using the dof coordinates
    geometry = V.tabulate_dof_coordinates()
    # As we are dealing with a vector field, we reshape the underlying dof array to accommedate for the three dimensional space
    num_dofs = V.dofmap.index_map.size_local + V.dofmap.index_map.num_ghosts
    values = np.zeros((num_dofs, 3), dtype=np.float64)
    values[:, :mesh.geometry.dim] = uh.x.array.real.reshape(num_dofs, V.dofmap.index_map_bs)

    # Create grid defined by the function space for visualization
    function_grid = pyvista.UnstructuredGrid(topology, cell_types, geometry)
    var_name = f"u({t})"
    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)
WARNING:py.warnings:/usr/local/dolfinx-real/lib/python3.8/dist-packages/dolfinx/plot.py:120: UserWarning: Plotting of higher order functions is experimental.
  warnings.warn("Plotting of higher order functions is experimental.")

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

from dolfinx.log import set_log_level, LogLevel
set_log_level(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-01-17 08:23:26.292 (   1.567s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:27.282 (   2.557s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:27.847 (   3.122s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 2: r (abs) = 22.2455 (tol = 1e-08) r (rel) = 0.134278(tol = 1e-08)
2022-01-17 08:23:28.097 (   3.373s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:28.679 (   3.955s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 3: r (abs) = 2.43261 (tol = 1e-08) r (rel) = 0.0146837(tol = 1e-08)
2022-01-17 08:23:28.932 (   4.208s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:29.557 (   4.833s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 4: r (abs) = 4.43158 (tol = 1e-08) r (rel) = 0.0267498(tol = 1e-08)
2022-01-17 08:23:29.809 (   5.084s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:30.387 (   5.663s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 5: r (abs) = 0.144189 (tol = 1e-08) r (rel) = 0.000870353(tol = 1e-08)
2022-01-17 08:23:30.631 (   5.907s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:31.213 (   6.488s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 6: r (abs) = 0.021424 (tol = 1e-08) r (rel) = 0.000129319(tol = 1e-08)
2022-01-17 08:23:31.465 (   6.740s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:32.047 (   7.323s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 7: r (abs) = 4.80068e-06 (tol = 1e-08) r (rel) = 2.89778e-08(tol = 1e-08)
2022-01-17 08:23:32.290 (   7.566s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:32.865 (   8.141s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 8: r (abs) = 2.62988e-11 (tol = 1e-08) r (rel) = 1.58745e-13(tol = 1e-08)
WARNING:py.warnings:/usr/local/dolfinx-real/lib/python3.8/dist-packages/dolfinx/plot.py:120: UserWarning: Plotting of higher order functions is experimental.
  warnings.warn("Plotting of higher order functions is experimental.")
Time step 1, Number of iterations 8, Load [ 0.   0.  -1.5]
2022-01-17 08:23:32.866 (   8.141s) [main            ]       NewtonSolver.cpp:250   INFO| Newton solver finished in 8 iterations and 8 linear solver iterations.
2022-01-17 08:23:33.213 (   8.489s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:34.025 (   9.300s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:34.602 (   9.878s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 2: r (abs) = 17.3254 (tol = 1e-08) r (rel) = 0.117842(tol = 1e-08)
2022-01-17 08:23:34.875 (  10.151s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:35.455 (  10.731s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 3: r (abs) = 5.14882 (tol = 1e-08) r (rel) = 0.0350207(tol = 1e-08)
2022-01-17 08:23:35.699 (  10.975s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:36.289 (  11.564s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 4: r (abs) = 7.24003 (tol = 1e-08) r (rel) = 0.0492445(tol = 1e-08)
2022-01-17 08:23:36.543 (  11.818s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:37.131 (  12.407s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 5: r (abs) = 0.777889 (tol = 1e-08) r (rel) = 0.00529096(tol = 1e-08)
2022-01-17 08:23:37.379 (  12.655s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:37.965 (  13.240s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 6: r (abs) = 1.25525 (tol = 1e-08) r (rel) = 0.00853785(tol = 1e-08)
2022-01-17 08:23:38.212 (  13.487s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:38.781 (  14.056s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 7: r (abs) = 0.00849512 (tol = 1e-08) r (rel) = 5.77813e-05(tol = 1e-08)
2022-01-17 08:23:39.026 (  14.301s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:39.632 (  14.907s) [main            ]       NewtonSolver.cpp:34    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-01-17 08:23:39.878 (  15.154s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:40.457 (  15.733s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 9: r (abs) = 1.70417e-10 (tol = 1e-08) r (rel) = 1.15913e-12(tol = 1e-08)
2022-01-17 08:23:40.458 (  15.733s) [main            ]       NewtonSolver.cpp:250   INFO| Newton solver finished in 9 iterations and 9 linear solver iterations.
2022-01-17 08:23:40.791 (  16.066s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:41.614 (  16.890s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:42.191 (  17.467s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 2: r (abs) = 10.0011 (tol = 1e-08) r (rel) = 0.0887471(tol = 1e-08)
2022-01-17 08:23:42.435 (  17.711s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:43.012 (  18.288s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 3: r (abs) = 5.33026 (tol = 1e-08) r (rel) = 0.0472992(tol = 1e-08)
2022-01-17 08:23:43.254 (  18.529s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:43.830 (  19.106s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 4: r (abs) = 11.9901 (tol = 1e-08) r (rel) = 0.106397(tol = 1e-08)
2022-01-17 08:23:44.076 (  19.351s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:44.647 (  19.923s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 5: r (abs) = 2.29702 (tol = 1e-08) r (rel) = 0.0203831(tol = 1e-08)
2022-01-17 08:23:44.899 (  20.174s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:45.481 (  20.756s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 6: r (abs) = 3.90234 (tol = 1e-08) r (rel) = 0.0346282(tol = 1e-08)
2022-01-17 08:23:45.730 (  21.005s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:46.297 (  21.573s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 7: r (abs) = 0.236535 (tol = 1e-08) r (rel) = 0.00209895(tol = 1e-08)
2022-01-17 08:23:46.548 (  21.823s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:47.129 (  22.404s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 8: r (abs) = 0.0427142 (tol = 1e-08) r (rel) = 0.000379034(tol = 1e-08)
2022-01-17 08:23:47.378 (  22.654s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
Time step 3, Number of iterations 10, Load [ 0.   0.  -4.5]
2022-01-17 08:23:47.962 (  23.237s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 9: r (abs) = 2.87798e-05 (tol = 1e-08) r (rel) = 2.55384e-07(tol = 1e-08)
2022-01-17 08:23:48.227 (  23.502s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:48.797 (  24.072s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 10: r (abs) = 6.0895e-10 (tol = 1e-08) r (rel) = 5.40365e-12(tol = 1e-08)
2022-01-17 08:23:48.797 (  24.072s) [main            ]       NewtonSolver.cpp:250   INFO| Newton solver finished in 10 iterations and 10 linear solver iterations.
2022-01-17 08:23:49.143 (  24.418s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:49.969 (  25.245s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:50.540 (  25.816s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 2: r (abs) = 5.50693 (tol = 1e-08) r (rel) = 0.0653918(tol = 1e-08)
2022-01-17 08:23:50.781 (  26.056s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:51.357 (  26.632s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 3: r (abs) = 26.2489 (tol = 1e-08) r (rel) = 0.311692(tol = 1e-08)
2022-01-17 08:23:51.601 (  26.876s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:52.186 (  27.462s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 4: r (abs) = 2.30927 (tol = 1e-08) r (rel) = 0.0274213(tol = 1e-08)
2022-01-17 08:23:52.442 (  27.717s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:53.028 (  28.304s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 5: r (abs) = 14.0562 (tol = 1e-08) r (rel) = 0.16691(tol = 1e-08)
2022-01-17 08:23:53.287 (  28.563s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:53.867 (  29.143s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 6: r (abs) = 0.222774 (tol = 1e-08) r (rel) = 0.00264532(tol = 1e-08)
2022-01-17 08:23:54.118 (  29.394s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:54.683 (  29.959s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 7: r (abs) = 0.286671 (tol = 1e-08) r (rel) = 0.00340406(tol = 1e-08)
2022-01-17 08:23:54.922 (  30.197s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:55.516 (  30.791s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 8: r (abs) = 0.000321869 (tol = 1e-08) r (rel) = 3.82203e-06(tol = 1e-08)
Time step 4, Number of iterations 9, Load [ 0.  0. -6.]
2022-01-17 08:23:55.766 (  31.042s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:56.334 (  31.610s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 9: r (abs) = 2.63797e-07 (tol = 1e-08) r (rel) = 3.13245e-09(tol = 1e-08)
2022-01-17 08:23:56.334 (  31.610s) [main            ]       NewtonSolver.cpp:250   INFO| Newton solver finished in 9 iterations and 9 linear solver iterations.
2022-01-17 08:23:56.669 (  31.945s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:57.516 (  32.791s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:58.080 (  33.356s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 2: r (abs) = 3.19462 (tol = 1e-08) r (rel) = 0.0496479(tol = 1e-08)
2022-01-17 08:23:58.327 (  33.602s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:58.906 (  34.182s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 3: r (abs) = 7.71429 (tol = 1e-08) r (rel) = 0.119888(tol = 1e-08)
2022-01-17 08:23:59.155 (  34.431s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:23:59.723 (  34.999s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 4: r (abs) = 0.850873 (tol = 1e-08) r (rel) = 0.0132235(tol = 1e-08)
2022-01-17 08:23:59.976 (  35.252s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:00.560 (  35.836s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 5: r (abs) = 0.371434 (tol = 1e-08) r (rel) = 0.0057725(tol = 1e-08)
2022-01-17 08:24:00.809 (  36.084s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:01.391 (  36.667s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 6: r (abs) = 0.00215066 (tol = 1e-08) r (rel) = 3.34236e-05(tol = 1e-08)
2022-01-17 08:24:01.649 (  36.925s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
Time step 5, Number of iterations 8, Load [ 0.   0.  -7.5]
2022-01-17 08:24:02.247 (  37.523s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 7: r (abs) = 2.54607e-06 (tol = 1e-08) r (rel) = 3.95687e-08(tol = 1e-08)
2022-01-17 08:24:02.487 (  37.762s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:03.060 (  38.335s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 8: r (abs) = 3.33871e-13 (tol = 1e-08) r (rel) = 5.18872e-15(tol = 1e-08)
2022-01-17 08:24:03.060 (  38.335s) [main            ]       NewtonSolver.cpp:250   INFO| Newton solver finished in 8 iterations and 8 linear solver iterations.
2022-01-17 08:24:03.411 (  38.686s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:04.251 (  39.526s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:04.832 (  40.107s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 2: r (abs) = 2.00649 (tol = 1e-08) r (rel) = 0.0395622(tol = 1e-08)
2022-01-17 08:24:05.087 (  40.362s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:05.667 (  40.943s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 3: r (abs) = 4.60977 (tol = 1e-08) r (rel) = 0.0908914(tol = 1e-08)
2022-01-17 08:24:05.917 (  41.192s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:06.482 (  41.757s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 4: r (abs) = 0.185372 (tol = 1e-08) r (rel) = 0.00365501(tol = 1e-08)
2022-01-17 08:24:06.729 (  42.005s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:07.309 (  42.585s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 5: r (abs) = 0.024688 (tol = 1e-08) r (rel) = 0.000486777(tol = 1e-08)
2022-01-17 08:24:07.566 (  42.842s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:08.136 (  43.411s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 6: r (abs) = 5.69255e-06 (tol = 1e-08) r (rel) = 1.12241e-07(tol = 1e-08)
Time step 6, Number of iterations 7, Load [ 0.  0. -9.]
2022-01-17 08:24:08.382 (  43.658s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:08.960 (  44.235s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 7: r (abs) = 2.58918e-11 (tol = 1e-08) r (rel) = 5.10512e-13(tol = 1e-08)
2022-01-17 08:24:08.960 (  44.235s) [main            ]       NewtonSolver.cpp:250   INFO| Newton solver finished in 7 iterations and 7 linear solver iterations.
2022-01-17 08:24:09.292 (  44.568s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:10.108 (  45.384s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:10.687 (  45.963s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 2: r (abs) = 1.38506 (tol = 1e-08) r (rel) = 0.0336622(tol = 1e-08)
2022-01-17 08:24:10.937 (  46.213s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:11.522 (  46.798s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 3: r (abs) = 3.03739 (tol = 1e-08) r (rel) = 0.07382(tol = 1e-08)
2022-01-17 08:24:11.770 (  47.046s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:12.347 (  47.623s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 4: r (abs) = 0.0412386 (tol = 1e-08) r (rel) = 0.00100225(tol = 1e-08)
2022-01-17 08:24:12.599 (  47.874s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
Time step 7, Number of iterations 6, Load [  0.    0.  -10.5]
2022-01-17 08:24:13.173 (  48.448s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 5: r (abs) = 0.00205057 (tol = 1e-08) r (rel) = 4.98364e-05(tol = 1e-08)
2022-01-17 08:24:13.419 (  48.695s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:13.996 (  49.271s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 6: r (abs) = 1.78867e-08 (tol = 1e-08) r (rel) = 4.34714e-10(tol = 1e-08)
2022-01-17 08:24:13.996 (  49.271s) [main            ]       NewtonSolver.cpp:250   INFO| Newton solver finished in 6 iterations and 6 linear solver iterations.
2022-01-17 08:24:14.338 (  49.614s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:15.179 (  50.455s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:15.747 (  51.022s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 2: r (abs) = 1.06336 (tol = 1e-08) r (rel) = 0.031085(tol = 1e-08)
2022-01-17 08:24:15.995 (  51.270s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:16.574 (  51.849s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 3: r (abs) = 2.0477 (tol = 1e-08) r (rel) = 0.0598598(tol = 1e-08)
2022-01-17 08:24:16.821 (  52.096s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:17.398 (  52.673s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 4: r (abs) = 0.00897719 (tol = 1e-08) r (rel) = 0.000262427(tol = 1e-08)
2022-01-17 08:24:17.649 (  52.924s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:18.263 (  53.538s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 5: r (abs) = 0.000167422 (tol = 1e-08) r (rel) = 4.89419e-06(tol = 1e-08)
Time step 8, Number of iterations 6, Load [  0.   0. -12.]
2022-01-17 08:24:18.511 (  53.786s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:19.077 (  54.353s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 6: r (abs) = 3.25874e-11 (tol = 1e-08) r (rel) = 9.52617e-13(tol = 1e-08)
2022-01-17 08:24:19.077 (  54.353s) [main            ]       NewtonSolver.cpp:250   INFO| Newton solver finished in 6 iterations and 6 linear solver iterations.
2022-01-17 08:24:19.409 (  54.685s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:20.268 (  55.543s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:20.856 (  56.131s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 2: r (abs) = 0.898789 (tol = 1e-08) r (rel) = 0.0309666(tol = 1e-08)
2022-01-17 08:24:21.101 (  56.377s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:21.673 (  56.949s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 3: r (abs) = 1.38354 (tol = 1e-08) r (rel) = 0.0476679(tol = 1e-08)
2022-01-17 08:24:21.924 (  57.200s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:22.500 (  57.776s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 4: r (abs) = 0.00185096 (tol = 1e-08) r (rel) = 6.37724e-05(tol = 1e-08)
2022-01-17 08:24:22.757 (  58.032s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
Time step 9, Number of iterations 6, Load [  0.    0.  -13.5]
2022-01-17 08:24:23.334 (  58.610s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 5: r (abs) = 7.87183e-06 (tol = 1e-08) r (rel) = 2.71213e-07(tol = 1e-08)
2022-01-17 08:24:23.582 (  58.857s) [main            ]              petsc.cpp:769   INFO| PETSc Krylov solver starting to solve system.
2022-01-17 08:24:24.166 (  59.441s) [main            ]       NewtonSolver.cpp:34    INFO| Newton iteration 6: r (abs) = 7.21799e-13 (tol = 1e-08) r (rel) = 2.48686e-14(tol = 1e-08)
2022-01-17 08:24:24.166 (  59.441s) [main            ]       NewtonSolver.cpp:250   INFO| Newton solver finished in 6 iterations and 6 linear solver iterations.