KratosMultiphysics / Kratos

Kratos Multiphysics (A.K.A Kratos) is a framework for building parallel multi-disciplinary simulation software. Modularity, extensibility and HPC are the main objectives. Kratos has BSD license and is written in C++ with extensive Python interface.
https://kratosmultiphysics.github.io/Kratos/
Other
1k stars 244 forks source link

[GeoMechanicsApplication] Investigate and test Hencky ( = natural = logarithmic ) strains #11643

Closed indigocoral closed 10 months ago

indigocoral commented 10 months ago

Related issue: #11642

WPK4FEM commented 10 months ago

Looked up Wikipedia ( https://en.wikipedia.org/wiki/Strain_(mechanics) ) to get the names of strain measures correctly: Engineering strain ( = Cauchy strain ) = (L1 -L0)/L0 Logarithm strain ( = Hencky strain = true strain ) = ln ( L1 / L0 ) Green strain = (1/2) * ( L1^2 - L0^2) / l0^2 or more usefull the finite strain derivations from the deformation gradient F ( at https://en.wikipedia.org/wiki/Finite_strain_theory#Finite_strain_tensors ): Cauchy strain tensor (right Cauchy–Green deformation tensor) C = F^T F Green-Lagrangian strain tensor: E = (1/2) ( C - I ) Logarithmic strain, Natural strain, True strain, or Hencky strain E = (1/2) ln( C )

For comparison with old D Settlement the Hencky strain is implemented. 2 copies exist ( in small_strain_U_Pw_diff_order_element and in U_Pw_small_strain_element. The updated Lagrangian elements inherit from these elements. )

Tests in applications/GeoMechanicsApplication/tests found so far: 0.

WPK4FEM commented 10 months ago

Separated the functionality such that the geomechanics application now only holds it once. Added 3 unit tests to check the behavior of CalculateHenckyStrain:

  1. Incompressible elongation: F is a diagonal tensor
  2. Rigid rotation F^T =F^-1, resulting strains should be zero
  3. Pure shear. Eigenvectors are at 45 and 135 deg, The principal Hencky strains should be 1 compressive, 1 tensile.

The set of these 3 should be sufficient as any deformation can be decomposed in contributions of these 3.

AlejandroCornejo commented 10 months ago

just for you to know, in the ConstitutiveLawsApp we compute the Henky strain as:

template<SizeType TVoigtSize>
void AdvancedConstitutiveLawUtilities<TVoigtSize>::CalculateHenckyStrain(
    const MatrixType& rCauchyTensor,
    Vector& rStrainVector
    )
{
    // Doing resize in case is needed
    if (rStrainVector.size() != VoigtSize)
        rStrainVector.resize(VoigtSize, false);

    // Declare the different matrix
    BoundedMatrixType eigen_values_matrix, eigen_vectors_matrix;

    // Decompose matrix
    MathUtils<double>::GaussSeidelEigenSystem(rCauchyTensor, eigen_vectors_matrix, eigen_values_matrix, 1.0e-16, 20);

    // Calculate the eigenvalues of the E matrix
    for (IndexType i = 0; i < Dimension; ++i) {
        eigen_values_matrix(i, i) = 0.5 * std::log(eigen_values_matrix(i, i));
    }

    // Calculate E matrix
    BoundedMatrixType E_matrix;
    MathUtils<double>::BDBtProductOperation(E_matrix, eigen_values_matrix, eigen_vectors_matrix);

    // Hencky Strain Calculation
    rStrainVector = MathUtils<double>::StrainTensorToVector(E_matrix, TVoigtSize);
}
WPK4FEM commented 10 months ago

Dear Alejandro Cornejo,

Thank you for your comment. I am aware of the existence of Hencky strain in the ConstitutiveLawApp. In fact I think that a former colleague copied it from there and placed two copies in the GeoMechanicsApplication. I would like these implementations to exist only once and be fully tested with nightly tests at its single place of existance. However I'm not aware of a possibility to do that without including the complete ConstitutiveLawApp into GeoMechanicsApplication. Therefore, I have reduced the number of copies in GeoMechanicsApplication to one and made some unit tests ( stretch, rigid rotation and pure shear ) for it such that at least there is test coverage there.

Regards, Wijtze Pieter Kikstra ( Deltares )