FrontISTR  5.7.0
Large-scale structural analysis program with finit element method
mhyperelastic Module Reference

This module provides functions for hyperelastic calculation. More...

Functions/Subroutines

subroutine cderiv (matl, sectType, ctn, itn, inv1b, inv2b, inv3b, dibdc, d2ibdc2, strain, direction, inv4b, dibdc_ani, d2ibdc2_ani)
 This subroutine calculates derivative of the invariant with respect to Cauchy-Green tensor. More...
 
subroutine cderiv_aniso (ctn, inv3, didc_3, d2idc2_3, direction, inv4b, dibdc_ani, d2ibdc2_ani)
 This subroutine calculates derivative of the 4th reduced invariant I4 with respect to Cauchy-Green tensor where I4 = I_3^{-1/3}*a_i*C_{ij}*a_j and a_i = direc(i). More...
 
subroutine calelasticarrudaboyce (matl, sectType, cijkl, strain)
 This subroutine provides elastic tangent coefficient for Arruda-Boyce hyperelastic material. More...
 
subroutine calupdateelasticarrudaboyce (matl, sectType, dstrain, dstress)
 This subroutine provides to update stress and strain for Arrude-Royce material. More...
 
subroutine calelasticmooneyrivlin (matl, sectType, cijkl, strain, hdflag)
 This subroutine provides elastic tangent coefficient for Mooney-Rivlin hyperelastic material. More...
 
subroutine calupdateelasticmooneyrivlin (matl, sectType, strain, stress, strain_energy, hdflag)
 This subroutine provides to update stress and strain for Mooney-Rivlin material. More...
 
subroutine calelasticmooneyrivlinaniso (matl, sectType, cijkl, strain, cdsys, hdflag)
 This subroutine provides elastic tangent coefficient for anisotropic Mooney-Rivlin hyperelastic material. More...
 
subroutine calupdateelasticmooneyrivlinaniso (matl, sectType, strain, stress, cdsys, hdflag)
 This subroutine provides to update stress and strain for anisotropic Mooney-Rivlin material. More...
 

Detailed Description

This module provides functions for hyperelastic calculation.

Function/Subroutine Documentation

◆ calelasticarrudaboyce()

subroutine mhyperelastic::calelasticarrudaboyce ( type( tmaterial ), intent(in)  matl,
integer, intent(in)  sectType,
real(kind=kreal), dimension(3,3,3,3), intent(out)  cijkl,
real(kind=kreal), dimension(6), intent(in)  strain 
)

This subroutine provides elastic tangent coefficient for Arruda-Boyce hyperelastic material.

Parameters
[in]matlmaterial rpoperties
[in]secttypenot used currently
[out]cijklconstitutive relation
[in]strainCauchy-Lagrange strain tensor

Definition at line 198 of file Hyperelastic.f90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ calelasticmooneyrivlin()

subroutine mhyperelastic::calelasticmooneyrivlin ( type( tmaterial ), intent(in)  matl,
integer, intent(in)  sectType,
real(kind=kreal), dimension(3,3,3,3), intent(out)  cijkl,
real(kind=kreal), dimension(6), intent(in)  strain,
integer(kind=kint), intent(in), optional  hdflag 
)

This subroutine provides elastic tangent coefficient for Mooney-Rivlin hyperelastic material.

Parameters
[in]matlmaterial rpoperties
[in]secttypenot used curr
[out]cijklconstitutive relation
[in]strainCauchy-Lagrange strain tensor

Definition at line 287 of file Hyperelastic.f90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ calelasticmooneyrivlinaniso()

subroutine mhyperelastic::calelasticmooneyrivlinaniso ( type( tmaterial ), intent(in)  matl,
integer, intent(in)  sectType,
real(kind=kreal), dimension(3,3,3,3), intent(out)  cijkl,
real(kind=kreal), dimension(6), intent(in)  strain,
real(kind=kreal), dimension(3,3), intent(in)  cdsys,
integer(kind=kint), intent(in), optional  hdflag 
)

This subroutine provides elastic tangent coefficient for anisotropic Mooney-Rivlin hyperelastic material.

Parameters
[in]matlmaterial rpoperties
[in]secttypenot used curr
[out]cijklconstitutive relation
[in]strainCauchy-Lagrange strain tensor

Definition at line 387 of file Hyperelastic.f90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ calupdateelasticarrudaboyce()

subroutine mhyperelastic::calupdateelasticarrudaboyce ( type( tmaterial ), intent(in)  matl,
integer, intent(in)  sectType,
real(kind=kreal), dimension(6), intent(in)  dstrain,
real(kind=kreal), dimension(6), intent(out)  dstress 
)

This subroutine provides to update stress and strain for Arrude-Royce material.

Definition at line 242 of file Hyperelastic.f90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ calupdateelasticmooneyrivlin()

subroutine mhyperelastic::calupdateelasticmooneyrivlin ( type( tmaterial ), intent(in)  matl,
integer, intent(in)  sectType,
real(kind=kreal), dimension(6), intent(in)  strain,
real(kind=kreal), dimension(6), intent(out)  stress,
real(kind=kreal), intent(out)  strain_energy,
integer(kind=kint), intent(in), optional  hdflag 
)

This subroutine provides to update stress and strain for Mooney-Rivlin material.

Parameters
[in]matlmaterial properties
[in]secttypenot used currently
[out]stress2nd Piola-Kirchhoff stress
[in]strainGreen-Lagrangen strain
[out]strain_energystrain energy

Definition at line 334 of file Hyperelastic.f90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ calupdateelasticmooneyrivlinaniso()

subroutine mhyperelastic::calupdateelasticmooneyrivlinaniso ( type( tmaterial ), intent(in)  matl,
integer, intent(in)  sectType,
real(kind=kreal), dimension(6), intent(in)  strain,
real(kind=kreal), dimension(6), intent(out)  stress,
real(kind=kreal), dimension(3,3), intent(in)  cdsys,
integer(kind=kint), intent(in), optional  hdflag 
)

This subroutine provides to update stress and strain for anisotropic Mooney-Rivlin material.

Parameters
[in]matlmaterial properties
[in]secttypenot used currently
[out]stress2nd Piola-Kirchhoff stress
[in]strainGreen-Lagrangen strain

Definition at line 440 of file Hyperelastic.f90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ cderiv()

subroutine mhyperelastic::cderiv ( type( tmaterial ), intent(in)  matl,
integer, intent(in)  sectType,
real(kind=kreal), dimension(3,3), intent(out)  ctn,
real(kind=kreal), dimension(3,3), intent(out)  itn,
real(kind=kreal), intent(out)  inv1b,
real(kind=kreal), intent(out)  inv2b,
real(kind=kreal), intent(out)  inv3b,
real(kind=kreal), dimension(3,3,3), intent(out)  dibdc,
real(kind=kreal), dimension(3,3,3,3,3), intent(out)  d2ibdc2,
real(kind=kreal), dimension(6), intent(in)  strain,
real(kind=kreal), dimension(3), intent(in), optional  direction,
real(kind=kreal), intent(out), optional  inv4b,
real(kind=kreal), dimension(3,3), intent(out), optional  dibdc_ani,
real(kind=kreal), dimension(3,3,3,3), intent(out), optional  d2ibdc2_ani 
)

This subroutine calculates derivative of the invariant with respect to Cauchy-Green tensor.

Parameters
[in]matlmaterial rpoperties
[in]secttypenot used currently
[out]inv1binvariants
[out]inv2binvariants
[out]inv3binvariants
[out]dibdcderivative of the invariant with respect to c(i,j)
[out]d2ibdc2derivative of the invariant with respect to c(i,j)
[in]strainCauchy-Lagrange strain tensor
[out]ctnright Cauchy-Green deformation tensor
[out]itnidentity tensor
[in]directiondirection vector of anisotropic materials
[out]inv4binvariants
[out]dibdc_aniderivative of the 4th anisotropic invariant
[out]d2ibdc2_ani2nd derivative of the 4th anisotropic invariant

Definition at line 16 of file Hyperelastic.f90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ cderiv_aniso()

subroutine mhyperelastic::cderiv_aniso ( real(kind=kreal), dimension(3,3), intent(in)  ctn,
real(kind=kreal), intent(in)  inv3,
real(kind=kreal), dimension(3,3), intent(in)  didc_3,
real(kind=kreal), dimension(3,3,3,3), intent(in)  d2idc2_3,
real(kind=kreal), dimension(3), intent(in)  direction,
real(kind=kreal), intent(out)  inv4b,
real(kind=kreal), dimension(3,3), intent(out)  dibdc_ani,
real(kind=kreal), dimension(3,3,3,3), intent(out)  d2ibdc2_ani 
)

This subroutine calculates derivative of the 4th reduced invariant I4 with respect to Cauchy-Green tensor where I4 = I_3^{-1/3}*a_i*C_{ij}*a_j and a_i = direc(i).

Parameters
[in]ctnright Cauchy-Green deformation tensor
[in]inv33rd invariants
[in]didc_3derivative of the 3rd invariant with respect to c(i,j)
[in]d2idc2_3derivative of the 3rd invariant with respect to c(i,j)
[in]directiondirection vector of anisotropic materials
[out]inv4b4th reduced invariant I4
[out]dibdc_aniderivative of the 4th anisotropic invariant
[out]d2ibdc2_ani2nd derivative of the 4th anisotropic invariant

Definition at line 143 of file Hyperelastic.f90.

Here is the caller graph for this function: