Merge branch 'master' into core/nonconvereged_solutions_for_nr_strategy

KratosMultiphysics · May 15, 2024 · bcdcb04 · bcdcb04
2 parents 949f214 + 085a2af
commit bcdcb04
Show file tree

Hide file tree

Showing 338 changed files with 130,061 additions and 119,499 deletions.
diff --git a/.gitignore b/.gitignore
@@ -80,6 +80,7 @@ kratos/scripts/wheels/README.md
 
 # Most used editor files
 *.vscode
+*.idx
 remote_*.md
 
 /compile_commands.json

diff --git a/CMakeLists.txt b/CMakeLists.txt
@@ -27,8 +27,8 @@ endif(POLICY CMP0057)
 
 # Set here the version number **** only update upon tagging a release!
 set (KratosMultiphysics_MAJOR_VERSION 9)
-set (KratosMultiphysics_MINOR_VERSION 4)
-set (KratosMultiphysics_PATCH_VERSION 6)
+set (KratosMultiphysics_MINOR_VERSION 5)
+set (KratosMultiphysics_PATCH_VERSION 0)
 
 # If KRATOS_SOURCE_DIR is not defined use the CMAKE_SOURCE_DIR
 if(NOT DEFINED KRATOS_SOURCE_DIR)

diff --git a/README.md b/README.md
@@ -4,13 +4,13 @@
 
 [![Release][release-image]][releases]
 <a href="https://github.com/KratosMultiphysics/Kratos/releases/latest"><img src="https://img.shields.io/github/release-date/KratosMultiphysics/Kratos?label="></a>
-<a href="https://github.com/KratosMultiphysics/Kratos/compare/Release-8.1...master"><img src="https://img.shields.io/github/commits-since/KratosMultiphysics/Kratos/latest?label=commits%20since"></a>
+<a href="https://github.com/KratosMultiphysics/Kratos/compare/Release-9.5...master"><img src="https://img.shields.io/github/commits-since/KratosMultiphysics/Kratos/latest?label=commits%20since"></a>
 <a href="https://github.com/KratosMultiphysics/Kratos/commit/master"><img src="https://img.shields.io/github/last-commit/KratosMultiphysics/Kratos?label=latest%20commit"></a>
 
 [![PyPI pyversions](https://img.shields.io/pypi/pyversions/KratosMultiphysics.svg)](https://pypi.org/project/KratosMultiphysics/)
 [![Downloads](https://pepy.tech/badge/KratosMultiphysics/month)](https://pepy.tech/project/KratosMultiphysics)
 
-[release-image]: https://img.shields.io/badge/release-9.4-green.svg?style=flat
+[release-image]: https://img.shields.io/badge/release-9.5-green.svg?style=flat
 [releases]: https://github.com/KratosMultiphysics/Kratos/releases
 
 [license-image]: https://img.shields.io/badge/license-BSD-green.svg?style=flat
@@ -157,6 +157,6 @@ Some users of the technologies developed in Kratos are:
 
 # How to cite Kratos?
 Please, use the following references when citing Kratos in your work.
-- Dadvand, P., Rossi, R. & Oñate, E. An Object-oriented Environment for Developing Finite Element Codes for Multi-disciplinary Applications. Arch Computat Methods Eng 17, 253–297 (2010). https://doi.org/10.1007/s11831-010-9045-2
-- Dadvand, P., Rossi, R., Gil, M., Martorell, X., Cotela, J., Juanpere, E., Idelsohn, S., Oñate, E. (2013). Migration of a generic multi-physics framework to HPC environments. Computers & Fluids. 80. 301–309. 10.1016/j.compfluid.2012.02.004.
-- Vicente Mataix Ferrándiz, Philipp Bucher, Rubén Zorrilla, Riccardo Rossi, Jordi Cotela, Alejandro Cornejo Velázquez, Miguel Angel Celigueta, Josep Maria, Tobias Teschemacher, Carlos Roig, Miguel Maso, Guillermo Casas, Suneth Warnakulasuriya, Marc Núñez, Pooyan Dadvand, Salva Latorre, Ignasi de Pouplana, Joaquín Irazábal González, Ferran Arrufat, … Javi Gárate. (2022). KratosMultiphysics/Kratos: Release 9.2 (v9.2). Zenodo. https://doi.org/10.5281/zenodo.3234644
+- [Dadvand, P., Rossi, R. & Oñate, E. An Object-oriented Environment for Developing Finite Element Codes for Multi-disciplinary Applications. Arch Computat Methods Eng 17, 253–297 (2010). https://doi.org/10.1007/s11831-010-9045-2](https://doi.org/10.1007/s11831-010-9045-2)
+- [Dadvand, P., Rossi, R., Gil, M., Martorell, X., Cotela, J., Juanpere, E., Idelsohn, S., Oñate, E. (2013). Migration of a generic multi-physics framework to HPC environments. Computers & Fluids. 80. 301–309. 10.1016/j.compfluid.2012.02.004.](10.1016/j.compfluid.2012.02.004)
+- [Vicente Mataix Ferrándiz, Philipp Bucher, Rubén Zorrilla, Suneth Warnakulasuriya, Riccardo Rossi, Alejandro Cornejo, jcotela, Carlos Roig, Josep Maria, tteschemacher, Miguel Masó, Guillermo Casas, Marc Núñez, Pooyan Dadvand, Salva Latorre, Ignasi de Pouplana, Joaquín Irazábal González, AFranci, Ferran Arrufat, riccardotosi, Aditya Ghantasala, Klaus Bernd Sautter, Peter Wilson, dbaumgaertner, Bodhinanda Chandra, Armin Geiser, Inigo Lopez, lluís, jgonzalezusua, Javi Gárate. (2024). KratosMultiphysics/Kratos: Release 9.5 (v9.5). Zenodo. https://doi.org/10.5281/zenodo.3234644](https://zenodo.org/records/6926179)
diff --git a/applications/CoSimulationApplication/custom_external_libraries/CoSimIO/CMakeLists.txt b/applications/CoSimulationApplication/custom_external_libraries/CoSimIO/CMakeLists.txt
@@ -85,7 +85,7 @@ if(MSVC)
     endif()
 
 elseif(${CMAKE_COMPILER_IS_GNUCXX})
-    set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++11")
+    set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++17")
     set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -std=c89")
 
     if (CO_SIM_IO_STRICT_COMPILER)
@@ -101,7 +101,7 @@ elseif(${CMAKE_COMPILER_IS_GNUCXX})
     # Note: This command makes sure that this option comes pretty late on the cmdline.
     link_libraries("$<$<AND:$<CXX_COMPILER_ID:GNU>,$<VERSION_GREATER:$<CXX_COMPILER_VERSION>,7.0>,$<VERSION_LESS:$<CXX_COMPILER_VERSION>,9.0>>:-lstdc++fs>")
 elseif(CMAKE_CXX_COMPILER_ID MATCHES "Clang")
-    set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++11")
+    set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++17")
     set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -std=c89")
 
     if (CO_SIM_IO_STRICT_COMPILER)
@@ -110,7 +110,7 @@ elseif(CMAKE_CXX_COMPILER_ID MATCHES "Clang")
     endif()
 
 elseif(CMAKE_CXX_COMPILER_ID STREQUAL "Intel")
-    set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++11")
+    set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++17")
     set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -std=c89")
 
     if (CO_SIM_IO_STRICT_COMPILER)
@@ -119,7 +119,7 @@ elseif(CMAKE_CXX_COMPILER_ID STREQUAL "Intel")
     endif()
 
 else()
-    set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++11 -Wall -Wpedantic")
+    set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++17 -Wall -Wpedantic")
     set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -std=c89 -Wall -Wpedantic")
 endif()
 

diff --git a/applications/ConstitutiveLawsApplication/constitutive_laws_application_variables.h b/applications/ConstitutiveLawsApplication/constitutive_laws_application_variables.h
@@ -26,7 +26,8 @@ namespace Kratos
 {
     enum class SofteningType {Linear = 0, Exponential = 1, HardeningDamage = 2, CurveFittingDamage = 3};
 
-    enum class TangentOperatorEstimation {Analytic = 0, FirstOrderPerturbation = 1, SecondOrderPerturbation = 2, Secant = 3, SecondOrderPerturbationV2 = 4, InitialStiffness = 5};
+    enum class TangentOperatorEstimation {Analytic = 0, FirstOrderPerturbation = 1, SecondOrderPerturbation = 2,
+        Secant = 3, SecondOrderPerturbationV2 = 4, InitialStiffness = 5, OrthogonalSecant = 6};
 
 ///@name  Functions
 ///@{

diff --git a/..._constitutive/auxiliary_files/cl_integrators/generic_cl_integrator_kinematic_plasticity.h b/..._constitutive/auxiliary_files/cl_integrators/generic_cl_integrator_kinematic_plasticity.h
@@ -243,7 +243,7 @@ class GenericConstitutiveLawIntegratorKinematicPlasticity
         const double Denominator
         )
     {
-        rTangent = rElasticMatrix - outer_prod(Vector(prod(rElasticMatrix, rGFluxVector)), Vector(prod(rElasticMatrix, rFFluxVector))) * Denominator;
+        noalias(rTangent) = rElasticMatrix - outer_prod(Vector(prod(rElasticMatrix, rGFluxVector)), Vector(prod(rElasticMatrix, rFFluxVector))) * Denominator;
     }
 
 
@@ -281,10 +281,10 @@ class GenericConstitutiveLawIntegratorKinematicPlasticity
     {
         BoundedArrayType deviator = ZeroVector(VoigtSize);
         BoundedArrayType h_capa = ZeroVector(VoigtSize);
-        double J2, tensile_indicator_factor, compression_indicator_factor, slope, hardening_parameter, equivalent_plastic_strain;
+        double J2, tensile_indicator_factor, compression_indicator_factor, slope, hardening_parameter, equivalent_plastic_strain, I1;
 
         YieldSurfaceType::CalculateEquivalentStress( rPredictiveStressVector, rStrainVector, rUniaxialStress, rValues);
-        const double I1 = rPredictiveStressVector[0] + rPredictiveStressVector[1] + rPredictiveStressVector[2];
+        AdvancedConstitutiveLawUtilities<VoigtSize>::CalculateI1Invariant(rPredictiveStressVector, I1);
         AdvancedConstitutiveLawUtilities<VoigtSize>::CalculateJ2Invariant(rPredictiveStressVector, I1, deviator, J2);
         CalculateDerivativeYieldSurface(rPredictiveStressVector, deviator, J2, rYieldSurfaceDerivative, rValues);
         CalculateDerivativePlasticPotential(rPredictiveStressVector, deviator, J2, rDerivativePlasticPotential, rValues);

diff --git a/...ion/custom_constitutive/auxiliary_files/cl_integrators/generic_cl_integrator_plasticity.h b/...ion/custom_constitutive/auxiliary_files/cl_integrators/generic_cl_integrator_plasticity.h
@@ -245,10 +245,10 @@ class GenericConstitutiveLawIntegratorPlasticity
     {
         array_1d<double, VoigtSize> deviator = ZeroVector(VoigtSize);
         array_1d<double, VoigtSize> h_capa = ZeroVector(VoigtSize);
-        double J2, tensile_indicator_factor, compression_indicator_factor, slope, hardening_parameter, equivalent_plastic_strain;
+        double J2, tensile_indicator_factor, compression_indicator_factor, slope, hardening_parameter, equivalent_plastic_strain, I1;
 
         YieldSurfaceType::CalculateEquivalentStress( rPredictiveStressVector, rStrainVector, rUniaxialStress, rValues);
-        const double I1 = rPredictiveStressVector[0] + rPredictiveStressVector[1] + rPredictiveStressVector[2];
+        AdvancedConstitutiveLawUtilities<VoigtSize>::CalculateI1Invariant(rPredictiveStressVector, I1);
         AdvancedConstitutiveLawUtilities<VoigtSize>::CalculateJ2Invariant(rPredictiveStressVector, I1, deviator, J2);
         CalculateFFluxVector(rPredictiveStressVector, deviator, J2, rFflux, rValues);
         CalculateGFluxVector(rPredictiveStressVector, deviator, J2, rGflux, rValues);
@@ -278,7 +278,7 @@ class GenericConstitutiveLawIntegratorPlasticity
         const double Denominator
         )
     {
-        rTangent = rElasticMatrix - outer_prod(Vector(prod(rElasticMatrix, rGFluxVector)), Vector(prod(rElasticMatrix, rFFluxVector))) * Denominator;
+        noalias(rTangent) = rElasticMatrix - outer_prod(Vector(prod(rElasticMatrix, rGFluxVector)), Vector(prod(rElasticMatrix, rFFluxVector))) * Denominator;
     }
 
 

diff --git a/...tion/custom_constitutive/auxiliary_files/plastic_potentials/von_mises_plastic_potential.h b/...tion/custom_constitutive/auxiliary_files/plastic_potentials/von_mises_plastic_potential.h
@@ -116,16 +116,10 @@ class VonMisesPlasticPotential
         )
     {
         array_1d<double, VoigtSize> first_vector, second_vector, third_vector;
-
-        AdvancedConstitutiveLawUtilities<VoigtSize>::CalculateFirstVector(first_vector);
         AdvancedConstitutiveLawUtilities<VoigtSize>::CalculateSecondVector(rDeviator, J2, second_vector);
-        AdvancedConstitutiveLawUtilities<VoigtSize>::CalculateThirdVector(rDeviator, J2, third_vector);
-
-        const double c1 = 0.0;
         const double c2 = std::sqrt(3.0);
-        const double c3 = 0.0;
 
-        noalias(rGFlux) = c1 * first_vector + c2 * second_vector + c3 * third_vector;
+        noalias(rGFlux) = c2 * second_vector;
     }
 
     /**

diff --git a/...ication/custom_constitutive/auxiliary_files/yield_surfaces/drucker_prager_yield_surface.h b/...ication/custom_constitutive/auxiliary_files/yield_surfaces/drucker_prager_yield_surface.h
@@ -234,7 +234,7 @@ class DruckerPragerYieldSurface
         AdvancedConstitutiveLawUtilities<VoigtSize>::CalculateFirstVector(first_vector);
         AdvancedConstitutiveLawUtilities<VoigtSize>::CalculateSecondVector(rDeviator, J2, second_vector);
 
-        const double friction_angle = r_material_properties[FRICTION_ANGLE] * Globals::Pi / 180.0;;
+        const double friction_angle = r_material_properties[FRICTION_ANGLE] * Globals::Pi / 180.0;
         const double sin_phi = std::sin(friction_angle);
         const double Root3 = std::sqrt(3.0);
 

diff --git a/...tom_constitutive/finite_strains/plasticity/generic_finite_strain_isotropic_plasticity.cpp b/...tom_constitutive/finite_strains/plasticity/generic_finite_strain_isotropic_plasticity.cpp
@@ -164,9 +164,7 @@ void GenericFiniteStrainIsotropicPlasticity<TConstLawIntegratorType>::
                 noalias(r_integrated_stress_vector) = predictive_stress_vector;
 
                 if (r_constitutive_law_options.Is(ConstitutiveLaw::COMPUTE_CONSTITUTIVE_TENSOR)) {
-                    this->CalculateTangentTensor(rValues, ConstitutiveLaw::StressMeasure_Kirchhoff);
-                } else {
-                    BaseType::CalculateElasticMatrix( r_constitutive_matrix, rValues);
+                    this->CalculateTangentTensor(rValues, ConstitutiveLaw::StressMeasure_Kirchhoff, plastic_strain);
                 }
             }
         }
@@ -301,7 +299,8 @@ template<class TConstLawIntegratorType>
 void GenericFiniteStrainIsotropicPlasticity<TConstLawIntegratorType>::
     CalculateTangentTensor(
         ConstitutiveLaw::Parameters& rValues,
-        const ConstitutiveLaw::StressMeasure& rStressMeasure
+        const ConstitutiveLaw::StressMeasure& rStressMeasure,
+        const Vector& rPlasticStrain
         )
 {
     const Properties& r_material_properties = rValues.GetMaterialProperties();
@@ -317,6 +316,10 @@ void GenericFiniteStrainIsotropicPlasticity<TConstLawIntegratorType>::
     } else if (tangent_operator_estimation == TangentOperatorEstimation::SecondOrderPerturbation) {
         // Calculates the Tangent Constitutive Tensor by perturbation (second order)
         TangentOperatorCalculatorUtility::CalculateTangentTensor(rValues, this, rStressMeasure, consider_perturbation_threshold, 2);
+    } else if (tangent_operator_estimation == TangentOperatorEstimation::Secant) {
+        const Vector num = prod(rValues.GetConstitutiveMatrix(), rPlasticStrain);
+        const double denom = inner_prod(rValues.GetStrainVector(), num);
+        noalias(rValues.GetConstitutiveMatrix()) -= outer_prod(num, num) / denom;
     }
 }
 

diff --git a/...ustom_constitutive/finite_strains/plasticity/generic_finite_strain_isotropic_plasticity.h b/...ustom_constitutive/finite_strains/plasticity/generic_finite_strain_isotropic_plasticity.h
@@ -315,7 +315,8 @@ class KRATOS_API(CONSTITUTIVE_LAWS_APPLICATION) GenericFiniteStrainIsotropicPlas
      */
     void CalculateTangentTensor(
         ConstitutiveLaw::Parameters &rValues,
-        const ConstitutiveLaw::StressMeasure& rStressMeasure = ConstitutiveLaw::StressMeasure_Cauchy
+        const ConstitutiveLaw::StressMeasure& rStressMeasure,
+        const Vector& rPlasticStrain
         );
 
     ///@}

diff --git a/...tom_constitutive/finite_strains/plasticity/generic_finite_strain_kinematic_plasticity.cpp b/...tom_constitutive/finite_strains/plasticity/generic_finite_strain_kinematic_plasticity.cpp
@@ -170,9 +170,7 @@ void GenericFiniteStrainKinematicPlasticity<TConstLawIntegratorType>::
                 noalias(r_integrated_stress_vector) = predictive_stress_vector;
 
                 if (r_constitutive_law_options.Is(ConstitutiveLaw::COMPUTE_CONSTITUTIVE_TENSOR)) {
-                    this->CalculateTangentTensor(rValues, ConstitutiveLaw::StressMeasure_Kirchhoff);
-                } else {
-                    BaseType::CalculateElasticMatrix( r_constitutive_matrix, rValues);
+                    this->CalculateTangentTensor(rValues, ConstitutiveLaw::StressMeasure_Kirchhoff, plastic_strain);
                 }
             }
         }
@@ -314,7 +312,8 @@ template<class TConstLawIntegratorType>
 void GenericFiniteStrainKinematicPlasticity<TConstLawIntegratorType>::
     CalculateTangentTensor(
         ConstitutiveLaw::Parameters& rValues,
-        const ConstitutiveLaw::StressMeasure& rStressMeasure
+        const ConstitutiveLaw::StressMeasure& rStressMeasure,
+        const Vector& rPlasticStrain
         )
 {
     const Properties& r_material_properties = rValues.GetMaterialProperties();
@@ -330,6 +329,10 @@ void GenericFiniteStrainKinematicPlasticity<TConstLawIntegratorType>::
     } else if (tangent_operator_estimation == TangentOperatorEstimation::SecondOrderPerturbation) {
         // Calculates the Tangent Constitutive Tensor by perturbation (second order)
         TangentOperatorCalculatorUtility::CalculateTangentTensor(rValues, this, ConstitutiveLaw::StressMeasure_Cauchy, consider_perturbation_threshold, 2);
+    } else if (tangent_operator_estimation == TangentOperatorEstimation::Secant) {
+        const Vector num = prod(rValues.GetConstitutiveMatrix(), rPlasticStrain);
+        const double denom = inner_prod(rValues.GetStrainVector(), num);
+        noalias(rValues.GetConstitutiveMatrix()) -= outer_prod(num, num) / denom;
     }
 }
 

diff --git a/...ustom_constitutive/finite_strains/plasticity/generic_finite_strain_kinematic_plasticity.h b/...ustom_constitutive/finite_strains/plasticity/generic_finite_strain_kinematic_plasticity.h
@@ -315,7 +315,8 @@ class KRATOS_API(CONSTITUTIVE_LAWS_APPLICATION) GenericFiniteStrainKinematicPlas
      */
     void CalculateTangentTensor(
         ConstitutiveLaw::Parameters &rValues,
-        const ConstitutiveLaw::StressMeasure& rStressMeasure = ConstitutiveLaw::StressMeasure_Cauchy
+        const ConstitutiveLaw::StressMeasure& rStressMeasure,
+        const Vector& rPlasticStrain
         );
 
     ///@}

diff --git a/...cation/custom_constitutive/small_strains/damage/generic_small_strain_isotropic_damage.cpp b/...cation/custom_constitutive/small_strains/damage/generic_small_strain_isotropic_damage.cpp
@@ -129,6 +129,7 @@ void GenericSmallStrainIsotropicDamage<TConstLawIntegratorType>::CalculateMateri
             noalias(r_integrated_stress_vector) = predictive_stress_vector;
 
             if (r_constitutive_law_options.Is(ConstitutiveLaw::COMPUTE_CONSTITUTIVE_TENSOR)) {
+                r_constitutive_matrix *= (1.0 - damage);
                 this->CalculateTangentTensor(rValues);
             }
         }
@@ -163,7 +164,7 @@ void GenericSmallStrainIsotropicDamage<TConstLawIntegratorType>::CalculateTangen
         // Calculates the Tangent Constitutive Tensor by perturbation (second order)
         TangentOperatorCalculatorUtility::CalculateTangentTensor(rValues, this, ConstitutiveLaw::StressMeasure_Cauchy, consider_perturbation_threshold, 2);
     } else if (tangent_operator_estimation == TangentOperatorEstimation::Secant) {
-        rValues.GetConstitutiveMatrix() *= (1.0 - mDamage);
+        return;
     } else if (tangent_operator_estimation == TangentOperatorEstimation::InitialStiffness) {
         return;
     }

diff --git a/...ustom_constitutive/small_strains/plasticity/generic_small_strain_isotropic_plasticity.cpp b/...ustom_constitutive/small_strains/plasticity/generic_small_strain_isotropic_plasticity.cpp
@@ -114,8 +114,8 @@ void GenericSmallStrainIsotropicPlasticity<TConstLawIntegratorType>::CalculateMa
         this->template AddInitialStrainVectorContribution<Vector>(r_strain_vector);
 
         // We compute the stress or the constitutive matrix
-        if (r_constitutive_law_options.Is( ConstitutiveLaw::COMPUTE_STRESS) ||
-            r_constitutive_law_options.Is( ConstitutiveLaw::COMPUTE_CONSTITUTIVE_TENSOR)) {
+        if (r_constitutive_law_options.Is(ConstitutiveLaw::COMPUTE_STRESS) ||
+            r_constitutive_law_options.Is(ConstitutiveLaw::COMPUTE_CONSTITUTIVE_TENSOR)) {
 
             // We get some variables
             double threshold = this->GetThreshold();
@@ -167,9 +167,7 @@ void GenericSmallStrainIsotropicPlasticity<TConstLawIntegratorType>::CalculateMa
                 noalias(r_integrated_stress_vector) = predictive_stress_vector;
 
                 if (r_constitutive_law_options.Is(ConstitutiveLaw::COMPUTE_CONSTITUTIVE_TENSOR)) {
-                    this->CalculateTangentTensor(rValues); // this modifies the ConstitutiveMatrix
-                } else {
-                    BaseType::CalculateElasticMatrix( r_constitutive_matrix, rValues);
+                    this->CalculateTangentTensor(rValues, plastic_strain); // this modifies the ConstitutiveMatrix
                 }
             }
         }
@@ -180,7 +178,9 @@ void GenericSmallStrainIsotropicPlasticity<TConstLawIntegratorType>::CalculateMa
 /***********************************************************************************/
 
 template <class TConstLawIntegratorType>
-void GenericSmallStrainIsotropicPlasticity<TConstLawIntegratorType>::CalculateTangentTensor(ConstitutiveLaw::Parameters& rValues)
+void GenericSmallStrainIsotropicPlasticity<TConstLawIntegratorType>::CalculateTangentTensor(
+    ConstitutiveLaw::Parameters& rValues,
+    const Vector& rPlasticStrain)
 {
     const Properties& r_material_properties = rValues.GetMaterialProperties();
 
@@ -195,11 +195,17 @@ void GenericSmallStrainIsotropicPlasticity<TConstLawIntegratorType>::CalculateTa
     } else if (tangent_operator_estimation == TangentOperatorEstimation::SecondOrderPerturbation) {
         // Calculates the Tangent Constitutive Tensor by perturbation (second order)
         TangentOperatorCalculatorUtility::CalculateTangentTensor(rValues, this, ConstitutiveLaw::StressMeasure_Cauchy, consider_perturbation_threshold, 2);
+    } else if (tangent_operator_estimation == TangentOperatorEstimation::Secant) {
+        const Vector num = prod(rValues.GetConstitutiveMatrix(), rPlasticStrain);
+        const double denom = inner_prod(rValues.GetStrainVector(), num);
+        noalias(rValues.GetConstitutiveMatrix()) -= outer_prod(num, num) / denom;
     } else if (tangent_operator_estimation == TangentOperatorEstimation::SecondOrderPerturbationV2) {
         // Calculates the Tangent Constitutive Tensor by perturbation (second order)
         TangentOperatorCalculatorUtility::CalculateTangentTensor(rValues, this, ConstitutiveLaw::StressMeasure_Cauchy, consider_perturbation_threshold, 4);
-    }  else if (tangent_operator_estimation == TangentOperatorEstimation::InitialStiffness) {
-        BaseType::CalculateElasticMatrix( rValues.GetConstitutiveMatrix(), rValues);
+    } else if (tangent_operator_estimation == TangentOperatorEstimation::InitialStiffness) {
+        BaseType::CalculateElasticMatrix(rValues.GetConstitutiveMatrix(), rValues);
+    } else if (tangent_operator_estimation == TangentOperatorEstimation::OrthogonalSecant) {
+        TangentOperatorCalculatorUtility::CalculateOrthogonalSecantTensor(rValues);
     }
 }
 
@@ -358,7 +364,6 @@ void GenericSmallStrainIsotropicPlasticity<TConstLawIntegratorType>::FinalizeMat
             plastic_dissipation, plastic_strain_increment,
             r_constitutive_matrix, plastic_strain, rValues,
             characteristic_length);
-        BaseType::CalculateElasticMatrix(r_constitutive_matrix, rValues);
     }
 
     mPlasticDissipation = plastic_dissipation;