Hyperelastic materials (#30)

* add three hyperelastic materials

Author: lijas

src/FiniteStrain/neohook.jl

@@ -0,0 +1,39 @@
+"""
+    NeoHook(; E, ν)
+
+Hyperelastic material
+#Arguments
+- `λ::Float64`: Lamé parameter
+- `μ::Float64`: Lamé parameter (shear modulus)
+"""
+
+struct NeoHook <: AbstractMaterial
+    λ::Float64
+    μ::Float64
+end
+
+struct NeoHookState <: AbstractMaterialState
+end
+
+function initial_material_state(::NeoHook)
+    return NeoHookState()
+end
+
+function NeoHook(; λ::T, μ::T) where T
+
+    return NeoHook(λ, μ)
+end
+
+function ψ(mp::NeoHook, C::SymmetricTensor{2,3})
+    J = sqrt(det(C))
+    I = tr(C)
+    return mp.μ/2 * (I-3) - mp.μ*log(J) + mp.λ/2 * log(J)^2
+end
+
+function material_response(mp::NeoHook, C::SymmetricTensor{2,3}, state::NeoHookState = NeoHookState(), 
+                           Δt=nothing; cache=nothing, options=nothing)
+    ∂²Ψ∂C², ∂Ψ∂C, _ =  hessian((C) -> ψ(mp, C), C, :all)
+    S = 2.0 * ∂Ψ∂C
+    ∂S∂C = 2.0 * ∂²Ψ∂C²
+    return S, ∂S∂C, NeoHookState()
+end

src/FiniteStrain/stvenant.jl

@@ -0,0 +1,38 @@
+"""
+    StVenant(; E, ν)
+
+Hyperelastic material
+#Arguments
+- `λ::Float64`: Lamé parameter
+- `μ::Float64`: Lamé parameter (shear modulus)
+"""
+
+struct StVenant <: AbstractMaterial
+    λ::Float64
+    μ::Float64
+end
+
+struct StVenantState <: AbstractMaterialState
+end
+
+function initial_material_state(::StVenant)
+    return StVenantState()
+end
+
+function StVenant(; λ::T, μ::T) where T
+    return StVenant(λ, μ)
+end
+
+function _SPK(mp::StVenant, C::SymmetricTensor{2,3})
+    I = one(C)
+    return 0.5*mp.λ*(tr(C) - 3)*I + mp.μ*(C - I)
+end
+
+function material_response(mp::StVenant, C::SymmetricTensor{2,3}, state::StVenantState = StVenantState(), 
+                           Δt=nothing; cache=nothing, options=nothing)
+                           
+    ∂S∂C, S =  gradient((C) -> _SPK(mp, C), C, :all)
+    return S, ∂S∂C, StVenantState()
+end
+
+

src/FiniteStrain/yeoh.jl

@@ -0,0 +1,42 @@
+"""
+    Yeoh(; λ, μ, c₂, c₃)
+
+Hyper elastic material
+#Arguments
+- `λ::Float64`: Lamé parameter
+- `μ::Float64`: Lamé parameter (shear modulus)
+- `c₂::Float64`: Material constant
+- `c₃::Float64`: Material constant
+"""
+struct Yeoh <: AbstractMaterial
+    λ::Float64
+    μ::Float64
+    c₂::Float64
+    c₃::Float64
+end
+
+struct YeohState <: AbstractMaterialState
+end
+
+function initial_material_state(::Yeoh)
+    return YeohState()
+end
+
+function Yeoh(; λ::T, μ::T, c₂::T, c₃::T) where T <: AbstractFloat
+    return Yeoh(λ, μ, c₂, c₃)
+end
+
+function ψ(mp::Yeoh, C::SymmetricTensor{2,3})
+    J = sqrt(det(C))
+    I = tr(C)
+    return mp.μ/2 * (I-3) + mp.c₂*(I-3)^2 + mp.c₃*(I-3)^3 - mp.μ*log(J) + mp.λ/2 * log(J)^2
+end
+
+function material_response(mp::Yeoh, C::SymmetricTensor{2,3}, state::YeohState = YeohState(), 
+                           Δt=nothing; cache=nothing, options=nothing)
+    ∂²Ψ∂C², ∂Ψ∂C, _ =  hessian((C) -> ψ(mp, C), C, :all)
+    S = 2.0 * ∂Ψ∂C
+    ∂S∂C = 2.0 * ∂²Ψ∂C²
+    return S, ∂S∂C, YeohState()
+end
+

src/MaterialModels.jl

@@ -72,6 +72,10 @@ include("CrystalViscoPlastic/slipsystems.jl")
 include("CrystalViscoPlastic/CrystalViscoPlastic.jl")
 include("CrystalViscoPlastic/CrystalViscoPlasticRed.jl")
 
+include("FiniteStrain/neohook.jl")
+include("FiniteStrain/yeoh.jl")
+include("FiniteStrain/stvenant.jl")
+
 include("nonlinear_solver.jl")
 include("wrappers.jl")
 
@@ -83,4 +87,5 @@ export LinearElastic, Plastic
 export LinearElasticState, PlasticState
 export OneD, UniaxialStrain, UniaxialStress, PlaneStrain, PlaneStress
 
+export NeoHook, Yeoh, StVenant
 end

test/jld2_files/NeoHook.jld2

@@ -10,3 +10,6 @@ include("test_plastic.jl")
 include("test_crystal_visco_plastic.jl")
 include("test_crystal_visco_plastic_red.jl")
 include("test_wrappers.jl")
+include("test_neohook.jl")
+include("test_yeoh.jl")
+include("test_stvenant.jl")

test/jld2_files/StVenant.jld2

@@ -0,0 +1,14 @@
+
+
+function get_NeoHook_loading()
+    _C = range(0.0,  0.005, length=3)
+    C = [SymmetricTensor{2,3}((x +1.0, x/10, 0.0, 1.0, 0.0, 1.0)) for x in _C]
+    return C
+end  
+
+@testset "NeoHook" begin
+    m = NeoHook(μ=76.9, λ=115.0)
+
+    loading = get_NeoHook_loading()
+    check_jld2(m, loading, "NeoHook")#, debug_print=false, OVERWRITE_JLD2=true)
+end

test/jld2_files/Yeoh.jld2

@@ -0,0 +1,13 @@
+function get_StVenant_loading()
+    _C = range(0.0,  0.005, length=3)
+    C = [SymmetricTensor{2,3}((x + 1.0, x/10, 0.0, 1.0, 0.0, 1.0)) for x in _C]
+
+    return C
+end  
+
+@testset "StVenant" begin
+    m = StVenant(μ=6.8, λ=62.1)
+
+    loading = get_StVenant_loading()
+    check_jld2(m, loading, "StVenant")#, debug_print=false, OVERWRITE_JLD2=true)
+end