Dynamic Hyperelasticity Problem

Hello
I am trying to implement a dynamic hyperelasticity problem (Neohookean Model). A static demo is available here
My question is all about implementation of the variational formulation. I looked at this document (page.9) and here is what I did for a time dependent hyperelasticity problem:

1
discretization in time:
2
Lets say \rho=1 and Stress is defined as:
3
Where J = det(F) and b is same as C in the FEniCS demo (Right Cauchy-Green tensor)
Finally after multiplying a test function (v) and integration by part, here is the variational form :

Here is a minimal code:

from dolfin import *
import numpy as np

A1 = [0,2,4,6,8]
A2 = [0,3,6,9,12]

mesh1 = Mesh("mesh.xml")
facets1 = MeshFunction("size_t", mesh1, "mesh_facet_region.xml")
domains1 = MeshFunction("size_t", mesh1, "mesh_physical_region.xml")

#function space
V = VectorFunctionSpace(mesh1, "Lagrange", 1)

W_const = Constant((0.0, 0.0))
u1= interpolate(W_const, V)
u0= interpolate(W_const, V)

dt = 0.5
t = 0
T = 2

# Define Dirichlet boundary conditions
U_const = Constant((0.0, 0.0))
bc_top = DirichletBC(V, U_const, facets1,35)

bcs = [bc_top]

# Define functions
du = TrialFunction(V)
v  = TestFunction(V)
u  = Function(V)

# Kinematics
d = u.geometric_dimension()
I = Identity(d)             # Identity tensor
F = I + grad(u)             # Deformation gradient
C = F.T*F                   # Right Cauchy-Green tensor

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

# Elasticity parameters
E, nu = 2E6, 0.37
mu  = Constant(E/(2*(1 + nu)))
lmbdaa= Constant(E*nu/((1 + nu)*(1 - 2*nu)))

sigma =  inv(J)*(lmbdaa*ln(J)*I+mu*(C_1-I))

# Define vectorial body force
delta_A2  = Expression(' alpha*t',degree=1, alpha=1, t=0)
delta_A1 = Expression(' beta*t',degree=1, beta=1, t=0)

force =  (delta_A2 - delta_A1)

def BF():
    return as_vector([force] + [0.0])

# Variational Form. Is that right?!
F =  inner(u,v)*dx + dt*dt*inner(sigma, grad(v)) * dx -dt*dt*inner(v, BF()) * dx - 2*inner(u1,v)*dx+inner(u0,v)*dx

counter = 0
while t <= T:

    delta_A1.t = A1[counter]
    delta_A2.t = A2[counter]

    J = derivative(F, u, du)

    problem = NonlinearVariationalProblem(F, u, bcs, J)
    solver = NonlinearVariationalSolver(problem)

    # Solve variational problem
    solver.solve()

    u0.assign(u1)
    u1.assign(u)

    t += dt
    counter = counter + 1

As I said this code works without error but I am not sure about the variational formulation. Is this variational form correct or I should take a new approach?

The form in your code looks reasonable. The formula you gave above, however, is missing a factor \Delta t^2 in front of the volume force term.

2 Likes

Thanks for your response.

Dear Leo , i am working on dynamic hyperelasticity problem as well. The FEniCS mechanics implementation offers one probably option to implement dynamic hyperelasticity using generalized alpha method. Are you able to use your mentioned variational formulation or did you face numerical issues. Please do share if you find feasible.

Thanks

Hi
I used the variational formulation you can see above and it worked.