From 8accfb939eaca69351307f5f6c460f88aaf88264 Mon Sep 17 00:00:00 2001 From: Nathan Miller Date: Sun, 7 Jan 2024 21:19:37 -0700 Subject: [PATCH 01/15] MAINT: Initial commit of the cohesion-based drucker prager second order yield surface declarations --- ..._hydraMicromorphicDruckerPragerPlasticity.h | 18 ++++++++++++++++++ 1 file changed, 18 insertions(+) diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h index 016df50a..2c880352 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h @@ -73,6 +73,24 @@ namespace tardigradeHydra{ typedef std::vector< constantType > constantVector; //!< Define a vector of constants typedef std::vector< std::vector< constantType > > constantMatrix; //!< Define a matrix of constants + virtual void computeSecondOrderDruckerPragerYieldSurface( const floatVector &stressMeasure, const floatType &cohesion, + const floatVector &preceedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + floatType &yieldValue ); + + virtual void computeSecondOrderDruckerPragerYieldSurface( const floatVector &stressMeasure, const floatType &cohesion, + const floatVector &preceedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + floatType &yieldValue, floatVector &dFdStress, floatType &dFdc, + floatVector &dFdPreceedingF, double tol = 1e-9 ); + + virtual void computeSecondOrderDruckerPragerYieldSurface( const floatVector &stressMeasure, const floatType &cohesion, + const floatVector &preceedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + floatType &yieldValue, floatVector &dFdStress, floatType &dFdc, + floatVector &dFdPreceedingF, floatMatrix &d2FdStress2, + floatMatrix &d2FdStressdPreceedingF, double tol = 1e-9 ); + /*! * The residual for a micromorphic Drucker Prager plasticity model */ From 93b7b1cf7ac1debe75f886791f23590e0751ccfc Mon Sep 17 00:00:00 2001 From: Nathan Miller Date: Sun, 7 Jan 2024 21:23:22 -0700 Subject: [PATCH 02/15] MAINT: Initial commit of the cohesion-based drucker prager higher order yield surface declarations --- ...hydraMicromorphicDruckerPragerPlasticity.h | 21 +++++++++++++++++++ 1 file changed, 21 insertions(+) diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h index 2c880352..044feaac 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h @@ -91,6 +91,27 @@ namespace tardigradeHydra{ floatVector &dFdPreceedingF, floatMatrix &d2FdStress2, floatMatrix &d2FdStressdPreceedingF, double tol = 1e-9 ); + virtual void computeHigherOrderDruckerPragerYieldEquation( const variableVector &stressMeasure, + const variableVector &cohesion, + const variableVector &preceedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + variableVector &yieldValue ); + + virtual void computeHigherOrderDruckerPragerYieldEquation( const variableVector &stressMeasure, + const variableVector &cohesion, + const variableVector &preceedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + variableVector &yieldValue, variableMatrix &dFdStress, variableMatrix &dFdc, + variableMatrix &dFdPreceedingF ); + + virtual void computeHigherOrderDruckerPragerYieldEquation( const variableVector &stressMEasure, + const variableVector &cohesion, + const variableVector &preceedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + variableVector &yieldValue, variableMatrix &dFdStress, variableMatrix &dFdc, + variableMatrix &dFdPreceedingF, variableMatrix &d2FdStress2, + variableMatrix &d2FdStressdPreceedingF ); + /*! * The residual for a micromorphic Drucker Prager plasticity model */ From edd59b4f1395414366fad2abcee0f827c1b7cc73 Mon Sep 17 00:00:00 2001 From: Nathan Miller Date: Sun, 7 Jan 2024 21:25:08 -0700 Subject: [PATCH 03/15] MAINT: Added declaration of the drucker prager parameter calculation --- src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h | 3 +++ 1 file changed, 3 insertions(+) diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h index 044feaac..118fa0c6 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h @@ -73,6 +73,9 @@ namespace tardigradeHydra{ typedef std::vector< constantType > constantVector; //!< Define a vector of constants typedef std::vector< std::vector< constantType > > constantMatrix; //!< Define a matrix of constants + virtual void computeDruckerPragerInternalParameters( const parameterType &frictionAngle, const parametertype &beta, + parameterType &A, parameterType &B ); + virtual void computeSecondOrderDruckerPragerYieldSurface( const floatVector &stressMeasure, const floatType &cohesion, const floatVector &preceedingDeformationGradient, const parameterType &frictionAngle, const parameterType &beta, From c975a5998fc17d5724636d9cbf1fb6bf96f8e7cc Mon Sep 17 00:00:00 2001 From: Nathan Miller Date: Sun, 7 Jan 2024 21:31:47 -0700 Subject: [PATCH 04/15] FEAT: Added the calculation of the Drucker-Prager internal parameters --- ...draMicromorphicDruckerPragerPlasticity.cpp | 33 ++++++++ ...hydraMicromorphicDruckerPragerPlasticity.h | 78 +++++++++---------- ...draMicromorphicDruckerPragerPlasticity.cpp | 23 ++++++ 3 files changed, 95 insertions(+), 39 deletions(-) diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index b429abf9..89e445a0 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -15,6 +15,39 @@ namespace tardigradeHydra{ namespace micromorphicDruckerPragerPlasticity{ + void computeDruckerPragerInternalParameters( const parameterType &frictionAngle, const parameterType &beta, + parameterType &A, parameterType &B ){ + /*! + * Compute the Drucker-Prager internal parameters + * + * :param const parameterType &frictionAngle: The material friction angle ( 0 < frictionAngle < pi / 2 ); + * :param const parameterType &beta: The beta parameter. + * :param parameterType &A: The A parameter. + * :param parameterType &B: The B parameter. + */ + + //Make sure the parameters are within bounds + TARDIGRADE_ERROR_TOOLS_CATCH( + if ( ( 0 > frictionAngle ) || ( frictionAngle > 1.570796 ) ){ + throw std::runtime_error( "The friction angle must be betwen 0 and pi / 2 not " + std::to_string( frictionAngle ) ); + } + ) + + TARDIGRADE_ERROR_TOOLS_CATCH( + if ( abs( beta ) > 1 ){ + throw std::runtime_error( "Beta must be between -1 and 1 not " + std::to_string( beta ) ); + } + ) + + //Compute the parameters + parameterType betaAngle = 2. * std::sqrt(6.) / ( 3. + beta * std::sin( frictionAngle ) ); + + A = betaAngle * std::cos( frictionAngle ); + + B = betaAngle * std::sin( frictionAngle ); + + } + void residual::setMacroDrivingStress( ){ /*! * Set the macro (i.e. the stress associated with the Cauchy stress) driving stress (stress in current configuration of plastic configuration) diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h index 118fa0c6..77b6b93a 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h @@ -73,47 +73,47 @@ namespace tardigradeHydra{ typedef std::vector< constantType > constantVector; //!< Define a vector of constants typedef std::vector< std::vector< constantType > > constantMatrix; //!< Define a matrix of constants - virtual void computeDruckerPragerInternalParameters( const parameterType &frictionAngle, const parametertype &beta, - parameterType &A, parameterType &B ); - - virtual void computeSecondOrderDruckerPragerYieldSurface( const floatVector &stressMeasure, const floatType &cohesion, - const floatVector &preceedingDeformationGradient, - const parameterType &frictionAngle, const parameterType &beta, - floatType &yieldValue ); - - virtual void computeSecondOrderDruckerPragerYieldSurface( const floatVector &stressMeasure, const floatType &cohesion, - const floatVector &preceedingDeformationGradient, - const parameterType &frictionAngle, const parameterType &beta, - floatType &yieldValue, floatVector &dFdStress, floatType &dFdc, - floatVector &dFdPreceedingF, double tol = 1e-9 ); - - virtual void computeSecondOrderDruckerPragerYieldSurface( const floatVector &stressMeasure, const floatType &cohesion, - const floatVector &preceedingDeformationGradient, - const parameterType &frictionAngle, const parameterType &beta, - floatType &yieldValue, floatVector &dFdStress, floatType &dFdc, - floatVector &dFdPreceedingF, floatMatrix &d2FdStress2, - floatMatrix &d2FdStressdPreceedingF, double tol = 1e-9 ); - - virtual void computeHigherOrderDruckerPragerYieldEquation( const variableVector &stressMeasure, - const variableVector &cohesion, - const variableVector &preceedingDeformationGradient, - const parameterType &frictionAngle, const parameterType &beta, - variableVector &yieldValue ); + void computeDruckerPragerInternalParameters( const parameterType &frictionAngle, const parameterType &beta, + parameterType &A, parameterType &B ); + + void computeSecondOrderDruckerPragerYieldSurface( const floatVector &stressMeasure, const floatType &cohesion, + const floatVector &preceedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + floatType &yieldValue ); + + void computeSecondOrderDruckerPragerYieldSurface( const floatVector &stressMeasure, const floatType &cohesion, + const floatVector &preceedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + floatType &yieldValue, floatVector &dFdStress, floatType &dFdc, + floatVector &dFdPreceedingF, double tol = 1e-9 ); + + void computeSecondOrderDruckerPragerYieldSurface( const floatVector &stressMeasure, const floatType &cohesion, + const floatVector &preceedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + floatType &yieldValue, floatVector &dFdStress, floatType &dFdc, + floatVector &dFdPreceedingF, floatMatrix &d2FdStress2, + floatMatrix &d2FdStressdPreceedingF, double tol = 1e-9 ); + + void computeHigherOrderDruckerPragerYieldEquation( const variableVector &stressMeasure, + const variableVector &cohesion, + const variableVector &preceedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + variableVector &yieldValue ); - virtual void computeHigherOrderDruckerPragerYieldEquation( const variableVector &stressMeasure, - const variableVector &cohesion, - const variableVector &preceedingDeformationGradient, - const parameterType &frictionAngle, const parameterType &beta, - variableVector &yieldValue, variableMatrix &dFdStress, variableMatrix &dFdc, - variableMatrix &dFdPreceedingF ); + void computeHigherOrderDruckerPragerYieldEquation( const variableVector &stressMeasure, + const variableVector &cohesion, + const variableVector &preceedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + variableVector &yieldValue, variableMatrix &dFdStress, variableMatrix &dFdc, + variableMatrix &dFdPreceedingF ); - virtual void computeHigherOrderDruckerPragerYieldEquation( const variableVector &stressMEasure, - const variableVector &cohesion, - const variableVector &preceedingDeformationGradient, - const parameterType &frictionAngle, const parameterType &beta, - variableVector &yieldValue, variableMatrix &dFdStress, variableMatrix &dFdc, - variableMatrix &dFdPreceedingF, variableMatrix &d2FdStress2, - variableMatrix &d2FdStressdPreceedingF ); + void computeHigherOrderDruckerPragerYieldEquation( const variableVector &stressMEasure, + const variableVector &cohesion, + const variableVector &preceedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + variableVector &yieldValue, variableMatrix &dFdStress, variableMatrix &dFdc, + variableMatrix &dFdPreceedingF, variableMatrix &d2FdStress2, + variableMatrix &d2FdStressdPreceedingF ); /*! * The residual for a micromorphic Drucker Prager plasticity model diff --git a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index 08f71542..80823780 100644 --- a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -1889,3 +1889,26 @@ BOOST_AUTO_TEST_CASE( test_setCohesion ){ BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousdMicroGradientCohesiondXAssembled, previousdMicroGradientCohesiondX ) ); } + +BOOST_AUTO_TEST_CASE( testComputeDruckerPragerInternalParameters ){ + /*! + * Test the computation of the internal parameters for the Drucker-Prager plasticity. + * + */ + + parameterType frictionAngle = 0.25; + parameterType beta = .423; + + parameterType Aanswer = 1.528893501990677; + parameterType Banswer = 0.39039060414065774; + + parameterType Aresult, Bresult; + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeDruckerPragerInternalParameters( frictionAngle, beta, Aresult, Bresult ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( Aresult, Aanswer ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( Bresult, Banswer ) ); + +} + From c9e3aac75559b456d73debee67c1465f205b547f Mon Sep 17 00:00:00 2001 From: Nathan Miller Date: Sun, 7 Jan 2024 21:59:16 -0700 Subject: [PATCH 05/15] WIP: Added second order Drucker Prager yield equation and tests --- ...draMicromorphicDruckerPragerPlasticity.cpp | 48 ++++++ ...hydraMicromorphicDruckerPragerPlasticity.h | 36 ++-- ...draMicromorphicDruckerPragerPlasticity.cpp | 162 ++++++++++++++++++ 3 files changed, 228 insertions(+), 18 deletions(-) diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index 89e445a0..1c3bc982 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -48,6 +48,54 @@ namespace tardigradeHydra{ } + void computeSecondOrderDruckerPragerYieldEquation( const variableVector &stressMeasure, const variableType &cohesion, + const variableVector &preceedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + variableType &yieldValue ){ + /*! + * Compute the second-order Drucker Prager Yield equation + * + * \f$ F = ||dev ( stressMeasure ) || - \left( A^{\phi} \bar{c} - B^{\phi} \bar{p} \right) \leq 0 + * + * || dev ( stressMeasure ) || = \sqrt{ dev( referenceStressMeasure ) : dev( referenceStressMeasure ) } + * dev( referenceStressMeasure ) : dev( referenceStressMeasure ) = dev( referenceStressMeasure )_{IJ} dev( referenceStressMeasure )_{IJ} + * dev( referenceStressMeasure )_{IJ} = referenceStressMeasure_{IJ} - \bar{p} elasticRightCauchyGreen_{IJ}^{-1} + * \bar{p} = \frac{1}{3} elasticRightCauchyGreen_{IJ} referenceStressMeasure_{IJ} + * A^{angle} = \beta^{angle} \cos( frictionAngle ) + * B^{angle} = \beta^{angle} \sin( frictionAngle ) + * \beta^{angle} = \frac{2 \sqrt{6} }{3 + \beta \sin( frictionAngle ) } + * \f$ + * + * :param const variableVector &stressMeasure: The stress measure + * :param const variableType &cohesion: The cohesion measure. + * :param const variableVector &preceedingDeformationGradient: The deformation gradients preceeding the configuration of the stress measure + * :param const parameterType &frictionAngle: The friction angle + * :param const parameterType &beta: The beta parameter + * :param variableType &yieldValue: The yield value. + */ + + parameterType AAngle, BAngle; + TARDIGRADE_ERROR_TOOLS_CATCH( tardigradeHydra::micromorphicDruckerPragerPlasticity::computeDruckerPragerInternalParameters( frictionAngle, beta, AAngle, BAngle ) ); + + //Compute the right Cauchy-Green deformation tensor + floatVector rightCauchyGreen; + TARDIGRADE_ERROR_TOOLS_CATCH_NODE_POINTER( tardigradeConstitutiveTools::computeRightCauchyGreen( preceedingDeformationGradient, rightCauchyGreen ) ); + + //Compute the decomposition of the stress + variableType pressure; + variableVector deviatoricReferenceStress; + + TARDIGRADE_ERROR_TOOLS_CATCH_NODE_POINTER( tardigradeMicromorphicTools::computeSecondOrderReferenceStressDecomposition( stressMeasure, + rightCauchyGreen, deviatoricReferenceStress, pressure ) ); + + //Compute the l2norm of the deviatoric stress + variableType normDevStress = tardigradeVectorTools::l2norm( deviatoricReferenceStress ); + + //Evaluate the yield equation + yieldValue = normDevStress - ( AAngle * cohesion - BAngle * pressure ); + + } + void residual::setMacroDrivingStress( ){ /*! * Set the macro (i.e. the stress associated with the Cauchy stress) driving stress (stress in current configuration of plastic configuration) diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h index 77b6b93a..eba5c89e 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h @@ -76,23 +76,23 @@ namespace tardigradeHydra{ void computeDruckerPragerInternalParameters( const parameterType &frictionAngle, const parameterType &beta, parameterType &A, parameterType &B ); - void computeSecondOrderDruckerPragerYieldSurface( const floatVector &stressMeasure, const floatType &cohesion, - const floatVector &preceedingDeformationGradient, - const parameterType &frictionAngle, const parameterType &beta, - floatType &yieldValue ); - - void computeSecondOrderDruckerPragerYieldSurface( const floatVector &stressMeasure, const floatType &cohesion, - const floatVector &preceedingDeformationGradient, - const parameterType &frictionAngle, const parameterType &beta, - floatType &yieldValue, floatVector &dFdStress, floatType &dFdc, - floatVector &dFdPreceedingF, double tol = 1e-9 ); - - void computeSecondOrderDruckerPragerYieldSurface( const floatVector &stressMeasure, const floatType &cohesion, - const floatVector &preceedingDeformationGradient, - const parameterType &frictionAngle, const parameterType &beta, - floatType &yieldValue, floatVector &dFdStress, floatType &dFdc, - floatVector &dFdPreceedingF, floatMatrix &d2FdStress2, - floatMatrix &d2FdStressdPreceedingF, double tol = 1e-9 ); + void computeSecondOrderDruckerPragerYieldEquation( const floatVector &stressMeasure, const floatType &cohesion, + const floatVector &preceedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + floatType &yieldValue ); + + void computeSecondOrderDruckerPragerYieldEquation( const floatVector &stressMeasure, const floatType &cohesion, + const floatVector &preceedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + floatType &yieldValue, floatVector &dFdStress, floatType &dFdc, + floatVector &dFdPreceedingF, double tol = 1e-9 ); + + void computeSecondOrderDruckerPragerYieldEquation( const floatVector &stressMeasure, const floatType &cohesion, + const floatVector &preceedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + floatType &yieldValue, floatVector &dFdStress, floatType &dFdc, + floatVector &dFdPreceedingF, floatMatrix &d2FdStress2, + floatMatrix &d2FdStressdPreceedingF, double tol = 1e-9 ); void computeHigherOrderDruckerPragerYieldEquation( const variableVector &stressMeasure, const variableVector &cohesion, @@ -150,7 +150,7 @@ namespace tardigradeHydra{ } - //Get the number of plastic multipliers expected for the problem + //!Get the number of plastic multipliers expected for the problem const unsigned int* getNumPlasticMultipliers( ){ return &_numPlasticMultipliers; } const unsigned int* getPlasticConfigurationIndex( ); diff --git a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index 80823780..035a4b85 100644 --- a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -1912,3 +1912,165 @@ BOOST_AUTO_TEST_CASE( testComputeDruckerPragerInternalParameters ){ } +BOOST_AUTO_TEST_CASE( testComputeSecondOrderDruckerPragerYieldEquation ){ + /*! + * Test the computation of the second order stress Drucker-Prager yield equation. + * + */ + + variableVector S = { 0.77189588, -0.84417528, 0.95929231, + -0.50465708, 0.50576944, 0.05335127, + 0.81510751, 0.76814059, -0.82146208 }; + + variableVector preceedingF = { 0.03468919, -0.31275742, -0.57541261, + -0.27865312, -0.45844965, 0.52325004, + -0.0439162 , -0.80201065, -0.44921044 }; + + variableType cohesion = 0.51; + parameterType frictionAngle = 0.46; + parameterType beta = 0.34; + + constantType answer = 3.236807543277929; + variableType result; + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, preceedingF, + frictionAngle, beta, result ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( result, answer ) ); + + //Test the Jacobian + variableType resultJ, dFdCohesion; + variableVector dFdS, dFdPreceedingF; + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, preceedingF, + frictionAngle, beta, resultJ, + dFdS, dFdCohesion, dFdPreceedingF ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( resultJ, answer ) ); + + //Test the Jacobian + variableType resultJ2, dFdcohesionJ2; + variableVector dFdSJ2, dFdPreceedingFJ2; + variableMatrix d2FdStress2J2; + variableMatrix d2FdStressdElasticRCGJ2; + variableMatrix d2FdS2J2, d2FdSdPreceedingFJ2; + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, preceedingF, + frictionAngle, beta, resultJ2, + dFdSJ2, dFdcohesionJ2, dFdPreceedingFJ2, d2FdS2J2, + d2FdSdPreceedingFJ2); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( resultJ2, answer ) ); + + //Test dFdStress + constantType eps = 1e-6; + for ( unsigned int i = 0; i < S.size(); i++ ){ + constantVector delta( S.size(), 0 ); + delta[i] = eps * fabs( S[i] ) + eps; + + floatType rp, rm; + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S + delta, cohesion, preceedingF, + frictionAngle, beta, rp ); + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S - delta, cohesion, preceedingF, + frictionAngle, beta, rm ); + + constantType gradCol = ( rp - rm ) / ( 2 * delta[i] ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradCol, dFdS[i] ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradCol, dFdSJ2[i] ) ); + + } + + //Test dFdPreceedingF + for ( unsigned int i = 0; i < preceedingF.size(); i++ ){ + constantVector delta( preceedingF.size(), 0 ); + delta[i] = eps * fabs( preceedingF[i] ) + eps; + + floatType rp, rm; + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, preceedingF + delta, + frictionAngle, beta, rp ); + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, preceedingF - delta, + frictionAngle, beta, rm ); + + constantType gradCol = ( rp - rm ) / ( 2 * delta[i] ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradCol, dFdPreceedingF[i] ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradCol, dFdPreceedingFJ2[i] ) ); + + } + + //Test dFdcohesion + for ( unsigned int i = 0; i < 1; i++ ){ + constantType delta = eps * fabs( cohesion ) + eps; + + floatType rp, rm; + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion + delta, preceedingF, + frictionAngle, beta, rp ); + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion - delta, preceedingF, + frictionAngle, beta, rm ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( ( rp - rm ) / ( 2 * delta ), dFdCohesion ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( ( rp - rm ) / ( 2 * delta ), dFdcohesionJ2 ) ); + + } + + //Test d2FdStress2 + for ( unsigned int i = 0; i < S.size(); i++ ){ + constantVector delta( S.size(), 0 ); + delta[i] = eps * fabs( S[i] ) + eps; + + floatType rp, rm; + + variableVector dFdSp, dFdSm, dFdPreceedingFp, dFdPreceedingFm; + variableType dFdcohesionp, dFdcohesionm; + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S + delta, cohesion, preceedingF, + frictionAngle, beta, rp, + dFdSp, dFdcohesionp, dFdPreceedingFp ); + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S - delta, cohesion, preceedingF, + frictionAngle, beta, rm, + dFdSm, dFdcohesionm, dFdPreceedingFm ); + + constantVector gradCol = ( dFdSp - dFdSm ) / ( 2 * delta[i] ); + + for ( unsigned int j = 0; j < gradCol.size(); j++ ){ + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradCol[j], d2FdS2J2[j][i] ) ); + } + } + + //Test d2FdStressdPreceedingF + for ( unsigned int i = 0; i < preceedingF.size(); i++ ){ + constantVector delta( preceedingF.size(), 0 ); + delta[i] = eps * fabs( preceedingF[i] ) + eps; + + floatType rp, rm; + variableVector dFdSp, dFdSm, dFdPreceedingFp, dFdPreceedingFm; + variableType dFdcohesionp, dFdcohesionm; + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, preceedingF + delta, + frictionAngle, beta, rp, + dFdSp, dFdcohesionp, dFdPreceedingFp ); + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, preceedingF - delta, + frictionAngle, beta, rm, + dFdSm, dFdcohesionm, dFdPreceedingFm ); + + constantVector gradCol = ( dFdSp - dFdSm ) / ( 2 * delta[i] ); + + for ( unsigned int j = 0; j < gradCol.size(); j++ ){ + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradCol[j], d2FdSdPreceedingFJ2[j][i] ) ); + } + } + +} + From fb2df644faaa6afb718a1415c586e8863bc83b59 Mon Sep 17 00:00:00 2001 From: Nathan Miller Date: Sun, 7 Jan 2024 22:57:01 -0700 Subject: [PATCH 06/15] FEAT: Added the calculation of the second order Drucker-Prager yield equation --- ...draMicromorphicDruckerPragerPlasticity.cpp | 222 +++++++++++++++++- ...hydraMicromorphicDruckerPragerPlasticity.h | 24 +- ...draMicromorphicDruckerPragerPlasticity.cpp | 80 +++---- 3 files changed, 271 insertions(+), 55 deletions(-) diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index 1c3bc982..37f0d780 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -49,7 +49,7 @@ namespace tardigradeHydra{ } void computeSecondOrderDruckerPragerYieldEquation( const variableVector &stressMeasure, const variableType &cohesion, - const variableVector &preceedingDeformationGradient, + const variableVector &precedingDeformationGradient, const parameterType &frictionAngle, const parameterType &beta, variableType &yieldValue ){ /*! @@ -68,7 +68,7 @@ namespace tardigradeHydra{ * * :param const variableVector &stressMeasure: The stress measure * :param const variableType &cohesion: The cohesion measure. - * :param const variableVector &preceedingDeformationGradient: The deformation gradients preceeding the configuration of the stress measure + * :param const variableVector &precedingDeformationGradient: The deformation gradients preceding the configuration of the stress measure * :param const parameterType &frictionAngle: The friction angle * :param const parameterType &beta: The beta parameter * :param variableType &yieldValue: The yield value. @@ -79,7 +79,7 @@ namespace tardigradeHydra{ //Compute the right Cauchy-Green deformation tensor floatVector rightCauchyGreen; - TARDIGRADE_ERROR_TOOLS_CATCH_NODE_POINTER( tardigradeConstitutiveTools::computeRightCauchyGreen( preceedingDeformationGradient, rightCauchyGreen ) ); + TARDIGRADE_ERROR_TOOLS_CATCH_NODE_POINTER( tardigradeConstitutiveTools::computeRightCauchyGreen( precedingDeformationGradient, rightCauchyGreen ) ); //Compute the decomposition of the stress variableType pressure; @@ -96,6 +96,222 @@ namespace tardigradeHydra{ } + void computeSecondOrderDruckerPragerYieldEquation( const variableVector &referenceStressMeasure, const variableType &cohesion, + const variableVector &precedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + variableType &yieldValue, variableVector &dFdStress, variableType &dFdc, + variableVector &dFdPrecedingF, double tol ){ + /*! + * Compute the second-order Drucker Prager Yield equation + * + * \f$F = ||dev ( stressMeasure ) || - \left( A^{\phi} \bar{c} - B^{\phi} \bar{p} \right) \leq 0 + * + * || dev ( stressMeasure ) || = \sqrt{ dev( referenceStressMeasure ) : dev( referenceStressMeasure ) } + * dev( referenceStressMeasure ) : dev( referenceStressMeasure ) = dev( referenceStressMeasure )_{IJ} dev( referenceStressMeasure )_{IJ} + * dev( referenceStressMeasure )_{IJ} = referenceStressMeasure_{IJ} - \bar{p} elasticRightCauchyGreen_{IJ}^{-1}\f$ + * + * Also compute the Jacobians + * \f$\frac{ \partial F }{ \partial stressMeasure_{IJ} } = \frac{ dev ( stressMeasure )_{AB} }{ || dev ( stressMeasure ) || } \frac{ \partial dev( stressMeasure ) \frac{ \partial dev( stressMeasure )_{AB} }{ \partial stressMeasure_{IJ} } + B^{\phi} \frac{ \partial \bar{p} }{ \partial stressMeasure_{IJ} } + * \frac{ \partial F }{ \partial \bar{c} } = -A^{\phi} + * \frac{ \partial F }{ \partial C_{IJ} } = \frac{ dev ( stressMeasure )_{AB} }{ || dev ( stressMeasure ) || } \frac{ \partial dev( stressMeasure ) \frac{ \partial dev( stressMeasure )_{AB} }{ \partial C_{IJ} } + B^{\phi} \frac{ \partial \bar{p} }{ \partial C_{IJ} } + * + * \bar{p} = \frac{1}{3} elasticRightCauchyGreen_{IJ} referenceStressMeasure_{IJ} + * A^{angle} = \beta^{angle} \cos( frictionAngle ) + * B^{angle} = \beta^{angle} \sin( frictionAngle ) + * \beta^{angle} = \frac{2 \sqrt{6} }{3 + \beta \sin( frictionAngle ) }\f$ + * + * :param const variableVector &referenceStressMeasure: The stress measure in the reference configuration + * :param const variableType &cohesion: The cohesion measure. + * :param const variableVector &precedingDeformationGradient: The preceding deformation gradient. + * :param const parameterType &frictionAngle: The friction angle + * :param const parameterType &beta: The beta parameter + * :param variableType &yieldValue: The yield value. + * :param variableVector &dFdStress: The Jacobian of the yield surface w.r.t. the stress measure. + * :param variableType &dFdc: The Jacobian of the yield surface w.r.t. the cohesion. + * :param variableVector &dFdPrecedingF: The Jacobian of the yield surface w.r.t. the preceding deformation gradient from the stress-measure's configuration + * :param double tol: The tolerance used to prevent nans in the Jacobians + */ + + parameterType AAngle, BAngle; + TARDIGRADE_ERROR_TOOLS_CATCH( computeDruckerPragerInternalParameters( frictionAngle, beta, AAngle, BAngle ) ); + + //Compute the right Cauchy-Green deformation tensor + floatVector rightCauchyGreen; + floatMatrix dRCGdPrecedingF; + + TARDIGRADE_ERROR_TOOLS_CATCH_NODE_POINTER( tardigradeConstitutiveTools::computeRightCauchyGreen( precedingDeformationGradient, rightCauchyGreen, dRCGdPrecedingF ) ); + + //Compute the decomposition of the stress + variableType pressure; + variableVector deviatoricReferenceStress; + + variableMatrix dDevStressdStress, dDevStressdRCG; + variableVector dPressuredStress, dPressuredRCG; + + TARDIGRADE_ERROR_TOOLS_CATCH_NODE_POINTER( tardigradeMicromorphicTools::computeSecondOrderReferenceStressDecomposition( referenceStressMeasure, + rightCauchyGreen, deviatoricReferenceStress, pressure, dDevStressdStress, + dDevStressdRCG, dPressuredStress, dPressuredRCG ) ); + + variableMatrix dDevStressdPrecedingF = tardigradeVectorTools::dot( dDevStressdRCG, dRCGdPrecedingF ); + variableVector dPressuredPrecedingF = tardigradeVectorTools::Tdot( dRCGdPrecedingF, dPressuredRCG ); + + //Compute the l2norm of the deviatoric stress + variableType normDevStress = tardigradeVectorTools::l2norm( deviatoricReferenceStress ); + + //Evaluate the yield equation + yieldValue = normDevStress - ( AAngle * cohesion - BAngle * pressure ); + + //Evaluate the jacobians + variableVector devStressDirection = deviatoricReferenceStress / ( normDevStress + tol ); + + dFdStress = tardigradeVectorTools::Tdot( dDevStressdStress, devStressDirection ) + + BAngle * dPressuredStress; + + dFdc = - AAngle; + + dFdPrecedingF = tardigradeVectorTools::Tdot( dDevStressdPrecedingF, devStressDirection ) + + BAngle * dPressuredPrecedingF; + } + + void computeSecondOrderDruckerPragerYieldEquation( const variableVector &stressMeasure, const variableType &cohesion, + const variableVector &precedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + variableType &yieldValue, variableVector &dFdStress, variableType &dFdc, + variableVector &dFdPrecedingF, variableMatrix &d2FdStress2, + variableMatrix &d2FdStressdPrecedingF, double tol ){ + /*! + * Compute the second-order Drucker Prager Yield equation + * + * \f$F = ||dev ( stressMeasure ) || - \left( A^{\phi} \bar{c} - B^{\phi} \bar{p} \right) \leq 0 + * + * || dev ( stressMeasure ) || = \sqrt{ dev( referenceStressMeasure ) : dev( referenceStressMeasure ) } + * dev( referenceStressMeasure ) : dev( referenceStressMeasure ) = dev( referenceStressMeasure )_{IJ} dev( referenceStressMeasure )_{IJ} + * dev( referenceStressMeasure )_{IJ} = referenceStressMeasure_{IJ} - \bar{p} elasticRightCauchyGreen_{IJ}^{-1} + * \bar{p} = \frac{1}{3} elasticRightCauchyGreen_{IJ} referenceStressMeasure_{IJ} + * A^{angle} = \beta^{angle} \cos( frictionAngle ) + * B^{angle} = \beta^{angle} \sin( frictionAngle ) + * \beta^{angle} = \frac{2 \sqrt{6} }{3 + \beta \sin( frictionAngle ) }\f$ + * + * Also compute the Jacobians + * \f$\frac{ \partial F }{ \partial stressMeasure_{IJ} } = \frac{ dev ( stressMeasure )_{AB} }{ || dev ( stressMeasure ) || } \frac{ \partial dev( stressMeasure ) \frac{ \partial dev( stressMeasure )_{AB} }{ \partial stressMeasure_{IJ} } + B^{\phi} \frac{ \partial \bar{p} }{ \partial stressMeasure_{IJ} } + * \frac{ \partial F }{ \partial \bar{c} } = -A^{\phi} + * \frac{ \partial F }{ \partial C_{IJ} } = \frac{ dev ( stressMeasure )_{AB} }{ || dev ( stressMeasure ) || } \frac{ \partial dev( stressMeasure )_{AB} }{ \partial C_{IJ} } + B^{\phi} \frac{ \partial \bar{p} }{ \partial C_{IJ} }\f$ + * + * The second deriatives of \f$\frac{ \partial F }{ \partial \stressMeasure_{IJ} }\f$ are + * \f$\frac{ \partial^2 F }{ \partial stressMeasure_{IJ} \partial stressMeasure_{KL} } = \frac{ \partial^2 || dev( stressMeasure ) || }{ \partial dev( stressMeasure )_{AB} \partial dev( stressMeasure )_{CD} } \frac{ \partial dev( stressMeasure )_{AB} } { \partial stressMeasure_{IJ} } \frac{ \partial dev( stressMeasure )_{CD} } { \partial stressMeasure_{KL} } + \frac{ dev ( stressMeasure )_{AB} }{ || dev( stressMeasure ) || } \frac{ \partial^2 dev( stressMeasure )_{AB} }{ \partial stressMeasure_{IJ} \partial stressMeasure_{KL} } + * \frac{ \partial^2 F }{ \partial stressMeasure_{IJ} \partial RCG_{KL} } = \frac{ \partial^2 || dev( stressMeasure ) || }{ \partial dev( stressMeasure )_{AB} \partial dev( stressMeasure )_{CD} } \frac{ \partial dev( stressMeasure )_{AB} } { \partial stressMeasure_{IJ} } \frac{ \partial dev( stressMeasure )_{CD} } { \partial C_{KL} } + \frac{ dev ( stressMeasure )_{AB} }{ || dev( stressMeasure ) || } \frac{ \partial^2 dev( stressMeasure )_{AB} }{ \partial stressMeasure_{IJ} \partial C_{KL} } + B^{\phi} \frac{ \partial^2 \bar{p} }{ \partial stressMeasure_{IJ} \partial C_{KL} } \f$ + * + * :param const variableVector &stressMeasure: The stress measure + * :param const variableType &cohesion: The cohesion measure. + * :param const variableVector &precedingDeformationGradient: The preceding deformation gradient + * :param const parameterType &frictionAngle: The friction angle + * :param const parameterType &beta: The beta parameter + * :param variableType &yieldValue: The yield value. + * :param variableVector &dFdStress: The Jacobian of the yield surface w.r.t. the stress measure. + * :param variableType &dFdc: The Jacobian of the yield surface w.r.t. the cohesion. + * :param variableVector &dFdElasticRCG: The Jacobian of the yield surface w.r.t. the elastic + * right Cauchy-Green deformation tensor. + * :param variableMatrix &d2FdStress2: The second derivative of the flow direction w.r.t. the stress. This + * is useful if one is using this expression as the flow potential and wants the jacobian of the flow direction \f$\frac{ \partial G }{\partial Stress_{IJ} }\f$ + * :param variableMatrix &d2FdStressdPrecedingF: The second derivative of the flow direction w.r.t. the stress and the preceding deformation gradient + * :param double tol: The tolerance used to prevent nans in the Jacobians + */ + //Assume 3D + unsigned int dim = 3; + + parameterType AAngle, BAngle; + TARDIGRADE_ERROR_TOOLS_CATCH( computeDruckerPragerInternalParameters( frictionAngle, beta, AAngle, BAngle ) ); + + //Compute the right Cauchy-Green deformation tensor + floatVector rightCauchyGreen; + floatMatrix dRCGdPrecedingF; + + TARDIGRADE_ERROR_TOOLS_CATCH_NODE_POINTER( tardigradeConstitutiveTools::computeRightCauchyGreen( precedingDeformationGradient, rightCauchyGreen, dRCGdPrecedingF ) ); + + //Compute the decomposition of the stress + variableType pressure; + variableVector deviatoricReferenceStress; + + variableMatrix dDevStressdStress, dDevStressdRCG; + variableVector dPressuredStress, dPressuredRCG; + + variableMatrix d2DevStressdStressdRCG, d2PressuredStressdRCG; + + TARDIGRADE_ERROR_TOOLS_CATCH_NODE_POINTER( tardigradeMicromorphicTools::computeSecondOrderReferenceStressDecomposition( stressMeasure, + rightCauchyGreen, deviatoricReferenceStress, pressure, dDevStressdStress, + dDevStressdRCG, dPressuredStress, dPressuredRCG, d2DevStressdStressdRCG, d2PressuredStressdRCG ) ) + + variableMatrix dDevStressdPrecedingF = tardigradeVectorTools::dot( dDevStressdRCG, dRCGdPrecedingF ); + variableVector dPressuredPrecedingF = tardigradeVectorTools::Tdot( dRCGdPrecedingF, dPressuredRCG ); + + variableMatrix d2DevStressdStressdPrecedingF( deviatoricReferenceStress.size( ), floatVector( deviatoricReferenceStress.size( ) * precedingDeformationGradient.size( ), 0 ) ); + + variableMatrix d2PressuredStressdPrecedingF( deviatoricReferenceStress.size( ), floatVector( precedingDeformationGradient.size( ), 0 ) ); + + for ( unsigned int I = 0; I < deviatoricReferenceStress.size( ); I++ ){ + + for ( unsigned int J = 0; J < deviatoricReferenceStress.size( ); J++ ){ + + for ( unsigned int K = 0; K < precedingDeformationGradient.size( ); K++ ){ + + for ( unsigned int L = 0; L < rightCauchyGreen.size( ); L++ ){ + + d2DevStressdStressdPrecedingF[ I ][ precedingDeformationGradient.size( ) * J + K ] + += d2DevStressdStressdRCG[ I ][ rightCauchyGreen.size( ) * J + L ] * dRCGdPrecedingF[ L ][ K ]; + + } + + d2PressuredStressdPrecedingF[ I ][ J ] + += d2PressuredStressdRCG[ I ][ K ] * dRCGdPrecedingF[ K ][ J ]; + + } + + } + + } + + //Compute the l2norm of the deviatoric stress + variableType normDevStress = tardigradeVectorTools::l2norm( deviatoricReferenceStress ); + + //Evaluate the yield equation + yieldValue = normDevStress - ( AAngle * cohesion - BAngle * pressure ); + + //Evaluate the jacobians + variableVector devStressDirection = deviatoricReferenceStress / ( normDevStress + tol ); + + dFdStress = tardigradeVectorTools::Tdot( dDevStressdStress, devStressDirection ) + + BAngle * dPressuredStress; + + dFdc = - AAngle; + + dFdPrecedingF = tardigradeVectorTools::Tdot( dDevStressdPrecedingF, devStressDirection ) + + BAngle * dPressuredPrecedingF; + + //Evaluate the second-order jacobians + constantMatrix EYE = tardigradeVectorTools::eye< constantType >( dim * dim ); + variableMatrix dDevStressDirectiondDevStress = ( EYE - tardigradeVectorTools::dyadic( devStressDirection, devStressDirection ) ) / ( normDevStress + tol ); + + d2FdStress2 = tardigradeVectorTools::Tdot( dDevStressdStress, tardigradeVectorTools::dot( dDevStressDirectiondDevStress, dDevStressdStress ) ); + + d2FdStressdPrecedingF = tardigradeVectorTools::Tdot( tardigradeVectorTools::dot( dDevStressDirectiondDevStress, dDevStressdStress ), dDevStressdPrecedingF ) + + BAngle * d2PressuredStressdPrecedingF; + + for ( unsigned int I = 0; I < dim; I++ ){ + for ( unsigned int J = 0; J < dim; J++ ){ + for ( unsigned int K = 0; K < dim; K++ ){ + for ( unsigned int L = 0; L < dim; L++ ){ + for ( unsigned int A = 0; A < dim; A++ ){ + for ( unsigned int B = 0; B < dim; B++ ){ + d2FdStressdPrecedingF[ dim * I + J ][ dim * K + L ] += devStressDirection[ dim * A + B ] * d2DevStressdStressdPrecedingF[ dim * A + B ][ dim * dim * dim * I + dim * dim * J + dim * K + L ]; + } + } + } + } + } + } + + } + void residual::setMacroDrivingStress( ){ /*! * Set the macro (i.e. the stress associated with the Cauchy stress) driving stress (stress in current configuration of plastic configuration) diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h index eba5c89e..d3be90f6 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h @@ -77,43 +77,43 @@ namespace tardigradeHydra{ parameterType &A, parameterType &B ); void computeSecondOrderDruckerPragerYieldEquation( const floatVector &stressMeasure, const floatType &cohesion, - const floatVector &preceedingDeformationGradient, + const floatVector &precedingDeformationGradient, const parameterType &frictionAngle, const parameterType &beta, floatType &yieldValue ); void computeSecondOrderDruckerPragerYieldEquation( const floatVector &stressMeasure, const floatType &cohesion, - const floatVector &preceedingDeformationGradient, + const floatVector &precedingDeformationGradient, const parameterType &frictionAngle, const parameterType &beta, floatType &yieldValue, floatVector &dFdStress, floatType &dFdc, - floatVector &dFdPreceedingF, double tol = 1e-9 ); + floatVector &dFdPrecedingF, double tol = 1e-9 ); void computeSecondOrderDruckerPragerYieldEquation( const floatVector &stressMeasure, const floatType &cohesion, - const floatVector &preceedingDeformationGradient, + const floatVector &precedingDeformationGradient, const parameterType &frictionAngle, const parameterType &beta, floatType &yieldValue, floatVector &dFdStress, floatType &dFdc, - floatVector &dFdPreceedingF, floatMatrix &d2FdStress2, - floatMatrix &d2FdStressdPreceedingF, double tol = 1e-9 ); + floatVector &dFdPrecedingF, floatMatrix &d2FdStress2, + floatMatrix &d2FdStressdPrecedingF, double tol = 1e-9 ); void computeHigherOrderDruckerPragerYieldEquation( const variableVector &stressMeasure, const variableVector &cohesion, - const variableVector &preceedingDeformationGradient, + const variableVector &precedingDeformationGradient, const parameterType &frictionAngle, const parameterType &beta, variableVector &yieldValue ); void computeHigherOrderDruckerPragerYieldEquation( const variableVector &stressMeasure, const variableVector &cohesion, - const variableVector &preceedingDeformationGradient, + const variableVector &precedingDeformationGradient, const parameterType &frictionAngle, const parameterType &beta, variableVector &yieldValue, variableMatrix &dFdStress, variableMatrix &dFdc, - variableMatrix &dFdPreceedingF ); + variableMatrix &dFdPrecedingF ); void computeHigherOrderDruckerPragerYieldEquation( const variableVector &stressMEasure, const variableVector &cohesion, - const variableVector &preceedingDeformationGradient, + const variableVector &precedingDeformationGradient, const parameterType &frictionAngle, const parameterType &beta, variableVector &yieldValue, variableMatrix &dFdStress, variableMatrix &dFdc, - variableMatrix &dFdPreceedingF, variableMatrix &d2FdStress2, - variableMatrix &d2FdStressdPreceedingF ); + variableMatrix &dFdPrecedingF, variableMatrix &d2FdStress2, + variableMatrix &d2FdStressdPrecedingF ); /*! * The residual for a micromorphic Drucker Prager plasticity model diff --git a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index 035a4b85..f7980f69 100644 --- a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -1922,43 +1922,43 @@ BOOST_AUTO_TEST_CASE( testComputeSecondOrderDruckerPragerYieldEquation ){ -0.50465708, 0.50576944, 0.05335127, 0.81510751, 0.76814059, -0.82146208 }; - variableVector preceedingF = { 0.03468919, -0.31275742, -0.57541261, - -0.27865312, -0.45844965, 0.52325004, - -0.0439162 , -0.80201065, -0.44921044 }; + variableVector precedingF = { 1.03482346, 0.01430697, 0.01134257, + 0.02756574, 1.03597345, 0.02115532, + 0.04903821, 0.03424149, 1.0240466 }; variableType cohesion = 0.51; parameterType frictionAngle = 0.46; parameterType beta = 0.34; - constantType answer = 3.236807543277929; + constantType answer = 1.5860256857059118; variableType result; - tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, preceedingF, + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, precedingF, frictionAngle, beta, result ); BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( result, answer ) ); //Test the Jacobian variableType resultJ, dFdCohesion; - variableVector dFdS, dFdPreceedingF; + variableVector dFdS, dFdPrecedingF; - tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, preceedingF, + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, precedingF, frictionAngle, beta, resultJ, - dFdS, dFdCohesion, dFdPreceedingF ); + dFdS, dFdCohesion, dFdPrecedingF ); BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( resultJ, answer ) ); //Test the Jacobian variableType resultJ2, dFdcohesionJ2; - variableVector dFdSJ2, dFdPreceedingFJ2; + variableVector dFdSJ2, dFdPrecedingFJ2; variableMatrix d2FdStress2J2; variableMatrix d2FdStressdElasticRCGJ2; - variableMatrix d2FdS2J2, d2FdSdPreceedingFJ2; + variableMatrix d2FdS2J2, d2FdSdPrecedingFJ2; - tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, preceedingF, + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, precedingF, frictionAngle, beta, resultJ2, - dFdSJ2, dFdcohesionJ2, dFdPreceedingFJ2, d2FdS2J2, - d2FdSdPreceedingFJ2); + dFdSJ2, dFdcohesionJ2, dFdPrecedingFJ2, d2FdS2J2, + d2FdSdPrecedingFJ2); BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( resultJ2, answer ) ); @@ -1970,10 +1970,10 @@ BOOST_AUTO_TEST_CASE( testComputeSecondOrderDruckerPragerYieldEquation ){ floatType rp, rm; - tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S + delta, cohesion, preceedingF, + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S + delta, cohesion, precedingF, frictionAngle, beta, rp ); - tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S - delta, cohesion, preceedingF, + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S - delta, cohesion, precedingF, frictionAngle, beta, rm ); constantType gradCol = ( rp - rm ) / ( 2 * delta[i] ); @@ -1984,24 +1984,24 @@ BOOST_AUTO_TEST_CASE( testComputeSecondOrderDruckerPragerYieldEquation ){ } - //Test dFdPreceedingF - for ( unsigned int i = 0; i < preceedingF.size(); i++ ){ - constantVector delta( preceedingF.size(), 0 ); - delta[i] = eps * fabs( preceedingF[i] ) + eps; + //Test dFdPrecedingF + for ( unsigned int i = 0; i < precedingF.size(); i++ ){ + constantVector delta( precedingF.size(), 0 ); + delta[i] = eps * fabs( precedingF[i] ) + eps; floatType rp, rm; - tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, preceedingF + delta, + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, precedingF + delta, frictionAngle, beta, rp ); - tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, preceedingF - delta, + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, precedingF - delta, frictionAngle, beta, rm ); constantType gradCol = ( rp - rm ) / ( 2 * delta[i] ); - BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradCol, dFdPreceedingF[i] ) ); + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradCol, dFdPrecedingF[i] ) ); - BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradCol, dFdPreceedingFJ2[i] ) ); + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradCol, dFdPrecedingFJ2[i] ) ); } @@ -2011,10 +2011,10 @@ BOOST_AUTO_TEST_CASE( testComputeSecondOrderDruckerPragerYieldEquation ){ floatType rp, rm; - tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion + delta, preceedingF, + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion + delta, precedingF, frictionAngle, beta, rp ); - tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion - delta, preceedingF, + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion - delta, precedingF, frictionAngle, beta, rm ); BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( ( rp - rm ) / ( 2 * delta ), dFdCohesion ) ); @@ -2030,16 +2030,16 @@ BOOST_AUTO_TEST_CASE( testComputeSecondOrderDruckerPragerYieldEquation ){ floatType rp, rm; - variableVector dFdSp, dFdSm, dFdPreceedingFp, dFdPreceedingFm; + variableVector dFdSp, dFdSm, dFdPrecedingFp, dFdPrecedingFm; variableType dFdcohesionp, dFdcohesionm; - tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S + delta, cohesion, preceedingF, + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S + delta, cohesion, precedingF, frictionAngle, beta, rp, - dFdSp, dFdcohesionp, dFdPreceedingFp ); + dFdSp, dFdcohesionp, dFdPrecedingFp ); - tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S - delta, cohesion, preceedingF, + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S - delta, cohesion, precedingF, frictionAngle, beta, rm, - dFdSm, dFdcohesionm, dFdPreceedingFm ); + dFdSm, dFdcohesionm, dFdPrecedingFm ); constantVector gradCol = ( dFdSp - dFdSm ) / ( 2 * delta[i] ); @@ -2048,27 +2048,27 @@ BOOST_AUTO_TEST_CASE( testComputeSecondOrderDruckerPragerYieldEquation ){ } } - //Test d2FdStressdPreceedingF - for ( unsigned int i = 0; i < preceedingF.size(); i++ ){ - constantVector delta( preceedingF.size(), 0 ); - delta[i] = eps * fabs( preceedingF[i] ) + eps; + //Test d2FdStressdPrecedingF + for ( unsigned int i = 0; i < precedingF.size(); i++ ){ + constantVector delta( precedingF.size(), 0 ); + delta[i] = eps * fabs( precedingF[i] ) + eps; floatType rp, rm; - variableVector dFdSp, dFdSm, dFdPreceedingFp, dFdPreceedingFm; + variableVector dFdSp, dFdSm, dFdPrecedingFp, dFdPrecedingFm; variableType dFdcohesionp, dFdcohesionm; - tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, preceedingF + delta, + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, precedingF + delta, frictionAngle, beta, rp, - dFdSp, dFdcohesionp, dFdPreceedingFp ); + dFdSp, dFdcohesionp, dFdPrecedingFp ); - tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, preceedingF - delta, + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( S, cohesion, precedingF - delta, frictionAngle, beta, rm, - dFdSm, dFdcohesionm, dFdPreceedingFm ); + dFdSm, dFdcohesionm, dFdPrecedingFm ); constantVector gradCol = ( dFdSp - dFdSm ) / ( 2 * delta[i] ); for ( unsigned int j = 0; j < gradCol.size(); j++ ){ - BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradCol[j], d2FdSdPreceedingFJ2[j][i] ) ); + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradCol[j], d2FdSdPrecedingFJ2[j][i] ) ); } } From 4232b0e3ae399c20e04e2e14fb28bc959b02f66b Mon Sep 17 00:00:00 2001 From: Nathan Miller Date: Sun, 7 Jan 2024 22:59:26 -0700 Subject: [PATCH 07/15] MAINT: Updated the formatting --- ...tardigrade_hydraMicromorphicDruckerPragerPlasticity.h | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h index d3be90f6..5cee432b 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h @@ -76,18 +76,21 @@ namespace tardigradeHydra{ void computeDruckerPragerInternalParameters( const parameterType &frictionAngle, const parameterType &beta, parameterType &A, parameterType &B ); - void computeSecondOrderDruckerPragerYieldEquation( const floatVector &stressMeasure, const floatType &cohesion, + void computeSecondOrderDruckerPragerYieldEquation( const floatVector &stressMeasure, + const floatType &cohesion, const floatVector &precedingDeformationGradient, const parameterType &frictionAngle, const parameterType &beta, floatType &yieldValue ); - void computeSecondOrderDruckerPragerYieldEquation( const floatVector &stressMeasure, const floatType &cohesion, + void computeSecondOrderDruckerPragerYieldEquation( const floatVector &stressMeasure, + const floatType &cohesion, const floatVector &precedingDeformationGradient, const parameterType &frictionAngle, const parameterType &beta, floatType &yieldValue, floatVector &dFdStress, floatType &dFdc, floatVector &dFdPrecedingF, double tol = 1e-9 ); - void computeSecondOrderDruckerPragerYieldEquation( const floatVector &stressMeasure, const floatType &cohesion, + void computeSecondOrderDruckerPragerYieldEquation( const floatVector &stressMeasure, + const floatType &cohesion, const floatVector &precedingDeformationGradient, const parameterType &frictionAngle, const parameterType &beta, floatType &yieldValue, floatVector &dFdStress, floatType &dFdc, From 6c3683f07c5a6ca80da06a091ec0f36fddd75f68 Mon Sep 17 00:00:00 2001 From: Nathan Miller Date: Mon, 8 Jan 2024 00:24:37 -0700 Subject: [PATCH 08/15] FEAT: Added the Drucker Prager yield surface for a higher order stress measure --- ...draMicromorphicDruckerPragerPlasticity.cpp | 291 ++++++++++++++++++ ...draMicromorphicDruckerPragerPlasticity.cpp | 182 +++++++++++ 2 files changed, 473 insertions(+) diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index 37f0d780..82ee5a5a 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -312,6 +312,297 @@ namespace tardigradeHydra{ } + void computeHigherOrderDruckerPragerYieldEquation( const variableVector &stressMeasure, + const variableVector &cohesion, + const variableVector &precedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + variableVector &yieldValue ){ + /*! + * Compute the higher-order Drucker Prager Yield equation + * + * \f$F_K = ||dev ( M ) ||_K - \left( A^{\phi} \bar{c}_K - B^{\phi} \bar{p}_K \right) \leq 0 + * + * || dev ( stressMeasure ) ||_K = \sqrt{ dev( M )_{IJK} : dev( M )_{IJK} }\f$ + * where the K's aren't summed. + * \f$dev( M )_{IJK} = M_{IJK} - \bar{p}_K elasticRightCauchyGreen_{IJ}^{-1} + * \bar{p} = \frac{1}{3} elasticRightCauchyGreen_{IJ} M_{IJK} + * A^{angle} = \beta^{angle} \cos( frictionAngle ) + * B^{angle} = \beta^{angle} \sin( frictionAngle ) + * \beta^{angle} = \frac{2 \sqrt{6} }{3 + \beta \sin( frictionAngle ) }\f$ + * + * :param const variableVector &stressMeasure: The stress measure + * :param const variableVector &cohesion: The cohesion measure. + * :param const variableVector &precedingDeformationGradient: The preceding deformation gradient + * :param const parameterType &frictionAngle: The friction angle + * :param const parameterType &beta: The beta parameter + * :param variableVector &yieldValue: The yield value. + */ + + parameterType AAngle, BAngle; + TARDIGRADE_ERROR_TOOLS_CATCH( computeDruckerPragerInternalParameters( frictionAngle, beta, AAngle, BAngle ) ) + + //Compute the right Cauchy-Green deformation tensor + floatVector rightCauchyGreen; + TARDIGRADE_ERROR_TOOLS_CATCH_NODE_POINTER( tardigradeConstitutiveTools::computeRightCauchyGreen( precedingDeformationGradient, rightCauchyGreen ) ); + + //Compute the decomposition of the stress + variableVector pressure; + variableVector deviatoricReferenceStress; + + TARDIGRADE_ERROR_TOOLS_CATCH_NODE_POINTER( tardigradeMicromorphicTools::computeHigherOrderReferenceStressDecomposition( stressMeasure, + rightCauchyGreen, deviatoricReferenceStress, pressure ) ); + + //Compute the l2norm of the deviatoric stress + variableVector normDevStress; + TARDIGRADE_ERROR_TOOLS_CATCH_NODE_POINTER( tardigradeMicromorphicTools::computeHigherOrderStressNorm( deviatoricReferenceStress, normDevStress ) ); + + //Evaluate the yield equation + yieldValue = normDevStress - ( AAngle * cohesion - BAngle * pressure ); + + } + + void computeHigherOrderDruckerPragerYieldEquation( const variableVector &stressMeasure, + const variableVector &cohesion, + const variableVector &precedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + variableVector &yieldValue, variableMatrix &dFdStress, variableMatrix &dFdc, + variableMatrix &dFdPrecedingF ){ + /*! + * Compute the higher-order Drucker Prager Yield equation + * + * \f$F_K = ||dev ( M ) ||_K - \left( A^{\phi} \bar{c}_K - B^{\phi} \bar{p}_K \right) \leq 0 + * + * || dev ( stressMeasure ) ||_K = \sqrt{ dev( M )_{IJK} : dev( M )_{IJK} }\f$ + * where the K's aren't summed. + * \f$dev( M )_{IJK} = M_{IJK} - \bar{p}_K elasticRightCauchyGreen_{IJ}^{-1} + * \bar{p} = \frac{1}{3} elasticRightCauchyGreen_{IJ} M_{IJK} + * A^{angle} = \beta^{angle} \cos( frictionAngle ) + * B^{angle} = \beta^{angle} \sin( frictionAngle ) + * \beta^{angle} = \frac{2 \sqrt{6} }{3 + \beta \sin( frictionAngle ) }\f$ + * + * Also computes the Jacobians + * + * :param const variableVector &stressMeasure: The stress measure + * :param const variableVector &cohesion: The cohesion measure. + * :param const variableVector &precedingDeformationGradient: The preceding deformation gradient + * :param const parameterType &frictionAngle: The friction angle + * :param const parameterType &beta: The beta parameter + * :param variableVector &yieldValue: The yield value. + * :param variableMatrix &dFdStress: The Jacobian of the yield function w.r.t. the reference higher order stress. + * :param variableMatrix &dFdc: The Jacobian of the yield function w.r.t. the cohesion. + * :param variableMatrix &dFdElasticRCG: The Jacobian of the yield function w.r.t. the preceding deformation gradient + * deformation tensor. + */ + + parameterType AAngle, BAngle; + TARDIGRADE_ERROR_TOOLS_CATCH( computeDruckerPragerInternalParameters( frictionAngle, beta, AAngle, BAngle ) ); + + //Compute the right Cauchy-Green deformation tensor + floatVector rightCauchyGreen; + floatMatrix dRCGdPrecedingF; + TARDIGRADE_ERROR_TOOLS_CATCH_NODE_POINTER( tardigradeConstitutiveTools::computeRightCauchyGreen( precedingDeformationGradient, rightCauchyGreen, dRCGdPrecedingF ) ); + + //Compute the decomposition of the stress + variableVector pressure; + variableVector deviatoricReferenceStress; + + variableMatrix dDevStressdStress, dDevStressdRCG; + variableMatrix dPressuredStress, dPressuredRCG; + + TARDIGRADE_ERROR_TOOLS_CATCH_NODE_POINTER( tardigradeMicromorphicTools::computeHigherOrderReferenceStressDecomposition( stressMeasure, + rightCauchyGreen, deviatoricReferenceStress, pressure, dDevStressdStress, + dDevStressdRCG, dPressuredStress, dPressuredRCG ) ); + + floatMatrix dDevStressdPrecedingF = tardigradeVectorTools::dot( dDevStressdRCG, dRCGdPrecedingF ); + floatMatrix dPressuredPrecedingF = tardigradeVectorTools::dot( dPressuredRCG, dRCGdPrecedingF ); + + //Compute the l2norm of the deviatoric stress + variableVector normDevStress; + variableMatrix dNormDevStressdDevStress; + TARDIGRADE_ERROR_TOOLS_CATCH_NODE_POINTER( tardigradeMicromorphicTools::computeHigherOrderStressNorm( deviatoricReferenceStress, normDevStress, dNormDevStressdDevStress ) ); + + //Evaluate the yield equation + yieldValue = normDevStress - ( AAngle * cohesion - BAngle * pressure ); + + //Construct the Jacobians + dFdStress = tardigradeVectorTools::dot( dNormDevStressdDevStress, dDevStressdStress ) + + BAngle * dPressuredStress; + + dFdc = -AAngle * tardigradeVectorTools::eye< constantType >( cohesion.size() ); + + dFdPrecedingF = tardigradeVectorTools::dot( dNormDevStressdDevStress, dDevStressdPrecedingF ) + + BAngle * dPressuredPrecedingF; + + } + + void computeHigherOrderDruckerPragerYieldEquation( const variableVector &stressMeasure, + const variableVector &cohesion, + const variableVector &precedingDeformationGradient, + const parameterType &frictionAngle, const parameterType &beta, + variableVector &yieldValue, variableMatrix &dFdStress, variableMatrix &dFdc, + variableMatrix &dFdPrecedingF, variableMatrix &d2FdStress2, + variableMatrix &d2FdStressdPrecedingF ){ + /*! + * Compute the higher-order Drucker Prager Yield equation + * + * \f$F_K = ||dev ( M ) ||_K - \left( A^{\phi} \bar{c}_K - B^{\phi} \bar{p}_K \right) \leq 0 + * + * || dev ( stressMeasure ) ||_K = \sqrt{ dev( M )_{IJK} : dev( M )_{IJK} }\f$ + * where the K's aren't summed. + * \f$dev( M )_{IJK} = M_{IJK} - \bar{p}_K elasticRightCauchyGreen_{IJ}^{-1} + * \bar{p} = \frac{1}{3} elasticRightCauchyGreen_{IJ} M_{IJK} + * A^{angle} = \beta^{angle} \cos( frictionAngle ) + * B^{angle} = \beta^{angle} \sin( frictionAngle ) + * \beta^{angle} = \frac{2 \sqrt{6} }{3 + \beta \sin( frictionAngle ) }\f$ + * + * Also computes the Jacobians + * + * :param const variableVector &stressMeasure: The stress measure + * :param const variableVector &cohesion: The cohesion measure. + * :param const variableVector &precedingDeformationGradient: The preceding deformation gradient + * :param const parameterType &frictionAngle: The friction angle + * :param const parameterType &beta: The beta parameter + * :param variableVector &yieldValue: The yield value. + * :param variableMatrix &dFdStress: The Jacobian of the yield function w.r.t. the reference higher order stress. + * :param variableMatrix &dFdc: The Jacobian of the yield function w.r.t. the cohesion. + * :param variableMatrix &dFdPrecedingF: The Jacobian of the yield function w.r.t. the preceding deformation gradient + * :param variableMatrix &d2FdStress2: The second order Jacobian of the yield function w.r.t. the reference + * higher order stress. + * :param variableMatrix &d2FdStressdPrecedingF: The second order Jacobian of the yield function w.r.t. the + * reference higher order stress and the preceding deformation gradient + */ + + //Assume 3D + unsigned int dim = 3; + + parameterType AAngle, BAngle; + TARDIGRADE_ERROR_TOOLS_CATCH( computeDruckerPragerInternalParameters( frictionAngle, beta, AAngle, BAngle ) ); + + //Compute the right Cauchy-Green deformation tensor + floatVector rightCauchyGreen; + floatMatrix dRCGdPrecedingF; + TARDIGRADE_ERROR_TOOLS_CATCH_NODE_POINTER( tardigradeConstitutiveTools::computeRightCauchyGreen( precedingDeformationGradient, rightCauchyGreen, dRCGdPrecedingF ) ); + + //Compute the decomposition of the stress + variableVector pressure; + variableVector deviatoricReferenceStress; + + variableMatrix dDevStressdStress, dDevStressdRCG; + variableMatrix dPressuredStress, dPressuredRCG; + + variableMatrix d2DevStressdStressdRCG, d2PressuredStressdRCG; + + TARDIGRADE_ERROR_TOOLS_CATCH_NODE_POINTER( tardigradeMicromorphicTools::computeHigherOrderReferenceStressDecomposition( stressMeasure, + rightCauchyGreen, deviatoricReferenceStress, pressure, dDevStressdStress, + dDevStressdRCG, dPressuredStress, dPressuredRCG, d2DevStressdStressdRCG, d2PressuredStressdRCG ) ) + + floatMatrix dDevStressdPrecedingF = tardigradeVectorTools::dot( dDevStressdRCG, dRCGdPrecedingF ); + floatMatrix dPressuredPrecedingF = tardigradeVectorTools::dot( dPressuredRCG, dRCGdPrecedingF ); + + variableMatrix d2DevStressdStressdPrecedingF( deviatoricReferenceStress.size( ), floatVector( deviatoricReferenceStress.size( ) * precedingDeformationGradient.size( ), 0 ) ); + + variableMatrix d2PressuredStressdPrecedingF( pressure.size( ), floatVector( deviatoricReferenceStress.size( ) * precedingDeformationGradient.size( ), 0 ) ); + + for ( unsigned int I = 0; I < deviatoricReferenceStress.size( ); I++ ){ + + for ( unsigned int J = 0; J < stressMeasure.size( ); J++ ){ + + for ( unsigned int K = 0; K < precedingDeformationGradient.size( ); K++ ){ + + for ( unsigned int L = 0; L < rightCauchyGreen.size( ); L++ ){ + + d2DevStressdStressdPrecedingF[ I ][ precedingDeformationGradient.size( ) * J + K ] + += d2DevStressdStressdRCG[ I ][ rightCauchyGreen.size( ) * J + L ] * dRCGdPrecedingF[ L ][ K ]; + + } + + } + + } + + } + + for ( unsigned int I = 0; I < pressure.size( ); I++ ){ + + for ( unsigned int J = 0; J < stressMeasure.size( ); J++ ){ + + for ( unsigned int K = 0; K < precedingDeformationGradient.size( ); K++ ){ + + for ( unsigned int L = 0; L < rightCauchyGreen.size( ); L++ ){ + + d2PressuredStressdPrecedingF[ I ][ precedingDeformationGradient.size( ) * J + K ] + += d2PressuredStressdRCG[ I ][ rightCauchyGreen.size( ) * J + L ] * dRCGdPrecedingF[ L ][ K ]; + + } + + } + + } + + } + + //Compute the l2norm of the deviatoric stress + variableVector normDevStress; + variableMatrix dNormDevStressdDevStress; + variableMatrix d2NormDevStressdDevStress2; + TARDIGRADE_ERROR_TOOLS_CATCH_NODE_POINTER( tardigradeMicromorphicTools::computeHigherOrderStressNorm( deviatoricReferenceStress, normDevStress, + dNormDevStressdDevStress, + d2NormDevStressdDevStress2 ) ) + + //Evaluate the yield equation + yieldValue = normDevStress - ( AAngle * cohesion - BAngle * pressure ); + + //Construct the Jacobians + dFdStress = tardigradeVectorTools::dot( dNormDevStressdDevStress, dDevStressdStress ) + + BAngle * dPressuredStress; + + dFdc = -AAngle * tardigradeVectorTools::eye< constantType >( cohesion.size() ); + + dFdPrecedingF = tardigradeVectorTools::dot( dNormDevStressdDevStress, dDevStressdPrecedingF ) + + BAngle * dPressuredPrecedingF; + + //Construct the second-order jacobians + d2FdStress2 = variableMatrix( dim, variableVector( dim * dim * dim * dim * dim * dim, 0 ) ); + d2FdStressdPrecedingF = tardigradeVectorTools::dot( dNormDevStressdDevStress, d2DevStressdStressdPrecedingF ) + + BAngle * d2PressuredStressdPrecedingF; + + for ( unsigned int K = 0; K < 3; K++ ){ + for ( unsigned int L = 0; L < 3; L++ ){ + for ( unsigned int M = 0; M < 3; M++ ){ + for ( unsigned int N = 0; N < 3; N++ ){ + for ( unsigned int O = 0; O < 3; O++ ){ + for ( unsigned int P = 0; P < 3; P++ ){ + for ( unsigned int Q = 0; Q < 3; Q++ ){ + for ( unsigned int A = 0; A < 3; A++ ){ + for ( unsigned int B = 0; B < 3; B++ ){ + for ( unsigned int C = 0; C < 3; C++ ){ + for ( unsigned int D = 0; D < 3; D++ ){ + for ( unsigned int E = 0; E < 3; E++ ){ + d2FdStressdPrecedingF[ K ][ dim * dim * dim * dim * L + dim * dim * dim * M + dim * dim * N + dim * O + P ] + += d2NormDevStressdDevStress2[ K ][ dim * dim * dim * dim * dim * Q + dim * dim * dim * dim * A + dim * dim * dim * B + dim * dim * C + dim * D + E ] + * dDevStressdStress[ dim * dim * Q + dim * A + B ][ dim * dim * L + dim * M + N ] + * dDevStressdPrecedingF[ dim * dim * C + dim * D + E ][ dim * O + P ]; + for ( unsigned int F = 0; F < 3; F++ ){ + d2FdStress2[ K ][ dim * dim * dim * dim * dim * L + dim * dim * dim * dim * M + dim * dim * dim * N + dim * dim * O + dim * P + Q ] + += d2NormDevStressdDevStress2[ K ][ dim * dim * dim * dim * dim * A + dim * dim * dim * dim * B + dim * dim * dim * C + dim * dim * D + dim * E + F ] + * dDevStressdStress[ dim * dim * A + dim * B + C ][ dim * dim * L + dim * M + N ] + * dDevStressdStress[ dim * dim * D + dim * E + F ][ dim * dim * O + dim * P + Q ]; + } + } + } + } + } + } + } + } + } + } + } + } + } + + } + void residual::setMacroDrivingStress( ){ /*! * Set the macro (i.e. the stress associated with the Cauchy stress) driving stress (stress in current configuration of plastic configuration) diff --git a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index f7980f69..6bc55aef 100644 --- a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -2074,3 +2074,185 @@ BOOST_AUTO_TEST_CASE( testComputeSecondOrderDruckerPragerYieldEquation ){ } +BOOST_AUTO_TEST_CASE( testComputeHigherOrderDruckerPragerYieldEquation ){ + /*! + * Test the computation of the higher order stress Drucker-Prager yield equation. + * + */ + + variableVector M = { 0.80732114, 0.79202055, 0.17990022, 0.97454675, 0.703207 , + -0.58236697, 0.53324571, -0.93438873, -0.40650796, 0.14071918, + 0.66933708, -0.67854069, -0.30317772, -0.93821882, 0.97270622, + 0.00295302, -0.12441126, 0.30539971, -0.0580227 , 0.89696105, + 0.17567709, -0.9592962 , 0.63535407, 0.95437804, -0.64531877, + 0.69978907, 0.81327586 }; + + variableVector precedingF = { 1.03482346, 0.01430697, 0.01134257, + 0.02756574, 1.03597345, 0.02115532, + 0.04903821, 0.03424149, 1.0240466 }; + + variableVector cohesion = { 0.51, 0.71, .14 }; + parameterType frictionAngle = 0.46; + parameterType beta = 0.34; + + constantVector answer = { 1.08862863, 1.39180044, 1.82194477 }; + variableVector result; + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeHigherOrderDruckerPragerYieldEquation( M, cohesion, precedingF, + frictionAngle, beta, result ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( result, answer ) ); + + //Test the Jacobians + + variableVector resultJ; + variableMatrix dFdStress, dFdc, dFdPrecedingF; + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeHigherOrderDruckerPragerYieldEquation( M, cohesion, precedingF, + frictionAngle, beta, resultJ, + dFdStress, dFdc, dFdPrecedingF ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( resultJ, answer ) ); + + variableVector resultJ2; + variableMatrix dFdStressJ2, dFdcJ2, dFdPrecedingFJ2; + variableMatrix d2FdStress2J2, d2FdStressdPrecedingFJ2; + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeHigherOrderDruckerPragerYieldEquation( M, cohesion, precedingF, + frictionAngle, beta, resultJ2, + dFdStressJ2, dFdcJ2, dFdPrecedingFJ2, + d2FdStress2J2, d2FdStressdPrecedingFJ2 ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( resultJ2, answer ) ); + + //Test derivatives w.r.t stress + constantType eps = 1e-6; + for ( unsigned int i = 0; i < M.size(); i++ ){ + constantVector delta( M.size(), 0 ); + delta[i] = eps * fabs( M[i] ) + eps; + + variableVector resultP, resultM; + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeHigherOrderDruckerPragerYieldEquation( M + delta, cohesion, precedingF, + frictionAngle, beta, resultP ); + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeHigherOrderDruckerPragerYieldEquation( M - delta, cohesion, precedingF, + frictionAngle, beta, resultM ); + + constantVector gradCol = ( resultP - resultM ) / ( 2 * delta[i] ); + + for ( unsigned int j = 0; j < gradCol.size(); j++ ){ + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradCol[j], dFdStress[j][i] ) ); + } + + for ( unsigned int j = 0; j < gradCol.size(); j++ ){ + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradCol[j], dFdStressJ2[j][i] ) ); + } + + variableMatrix dFdStressP, dFdcP, dFdPrecedingFP; + variableMatrix dFdStressM, dFdcM, dFdPrecedingFM; + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeHigherOrderDruckerPragerYieldEquation( M + delta, cohesion, precedingF, + frictionAngle, beta, resultP, + dFdStressP, dFdcP, dFdPrecedingFP ); + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeHigherOrderDruckerPragerYieldEquation( M - delta, cohesion, precedingF, + frictionAngle, beta, resultM, + dFdStressM, dFdcM, dFdPrecedingFM ); + + constantMatrix gradMat = ( dFdStressP - dFdStressM ) / ( 2 * delta[i] ); + + unsigned int n, o, p; + + for ( unsigned int j = 0; j < 3; j++ ){ + for ( unsigned int k = 0; k < 3; k++ ){ + for ( unsigned int l = 0; l < 3; l++ ){ + for ( unsigned int m = 0; m < 3; m++ ){ + n = ( int )( i / 9 ); + o = ( int )( (i - 9 * n ) / 3 ); + p = ( i - 9 * n - 3 * o ) % 3; + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradMat[ j ][ 9 * k + 3 * l + m ], + d2FdStress2J2[ j ][ 243 * k + 81 * l + 27 * m + 9 * n + 3 * o + p ] ) ); + } + } + } + } + } + + //Test derivatives w.r.t. the cohesion + for ( unsigned int i = 0; i < cohesion.size(); i++ ){ + constantVector delta( cohesion.size(), 0 ); + delta[i] = eps * fabs( cohesion[i] ) + eps; + + variableVector resultP, resultM; + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeHigherOrderDruckerPragerYieldEquation( M, cohesion + delta, precedingF, + frictionAngle, beta, resultP ); + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeHigherOrderDruckerPragerYieldEquation( M, cohesion - delta, precedingF, + frictionAngle, beta, resultM ); + + constantVector gradCol = ( resultP - resultM ) / ( 2 * delta[i] ); + + for ( unsigned int j = 0; j < gradCol.size(); j++ ){ + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradCol[j], dFdc[j][i] ) ); + } + + for ( unsigned int j = 0; j < gradCol.size(); j++ ){ + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradCol[j], dFdcJ2[j][i] ) ); + } + } + + //Test derivatives w.r.t. the preceding deformation gradient + for ( unsigned int i = 0; i < precedingF.size(); i++ ){ + constantVector delta( precedingF.size(), 0 ); + delta[i] = eps * fabs( precedingF[i] ) + eps; + + variableVector resultP, resultM; + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeHigherOrderDruckerPragerYieldEquation( M, cohesion, precedingF + delta, + frictionAngle, beta, resultP ); + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeHigherOrderDruckerPragerYieldEquation( M, cohesion, precedingF - delta, + frictionAngle, beta, resultM ); + + constantVector gradCol = ( resultP - resultM ) / ( 2 * delta[i] ); + + for ( unsigned int j = 0; j < gradCol.size(); j++ ){ + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradCol[j], dFdPrecedingF[j][i] ) ); + } + + for ( unsigned int j = 0; j < gradCol.size(); j++ ){ + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradCol[j], dFdPrecedingF[j][i] ) ); + } + + variableMatrix dFdStressP, dFdcP, dFdPrecedingFP; + variableMatrix dFdStressM, dFdcM, dFdPrecedingFM; + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeHigherOrderDruckerPragerYieldEquation( M, cohesion, precedingF + delta, + frictionAngle, beta, resultP, + dFdStressP, dFdcP, dFdPrecedingFP ); + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeHigherOrderDruckerPragerYieldEquation( M, cohesion, precedingF - delta, + frictionAngle, beta, resultM, + dFdStressM, dFdcM, dFdPrecedingFM ); + + constantMatrix gradMat = ( dFdStressP - dFdStressM ) / ( 2 * delta[i] ); + + unsigned int n, o; + + for ( unsigned int j = 0; j < 3; j++ ){ + for ( unsigned int k = 0; k < 3; k++ ){ + for ( unsigned int l = 0; l < 3; l++ ){ + for ( unsigned int m = 0; m < 3; m++ ){ + n = ( int )( i / 3 ); + o = ( i % 3 ); + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( gradMat[ j ][ 9 * k + 3 * l + m ], + d2FdStressdPrecedingFJ2[ j ][ 81 * k + 27 * l + 9 * m + 3 * n + o ] ) ); + } + } + } + } + } + +} From 8e703e0db55e0b3365d9e9f245ecade05a3f2a45 Mon Sep 17 00:00:00 2001 From: Nathan Miller Date: Mon, 8 Jan 2024 08:50:49 -0700 Subject: [PATCH 09/15] FEAT: Added the calculation of the flow potential gradients --- ...draMicromorphicDruckerPragerPlasticity.cpp | 226 +++++++++++ ...hydraMicromorphicDruckerPragerPlasticity.h | 52 ++- ...draMicromorphicDruckerPragerPlasticity.cpp | 355 +++++++++++++++++- 3 files changed, 625 insertions(+), 8 deletions(-) diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index 82ee5a5a..0a946f34 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -1429,6 +1429,232 @@ namespace tardigradeHydra{ } + void residual::setdMacroFlowdc( ){ + /*! + * Set the derivative of the macro flow potential w.r.t. the cohesion + */ + + setFlowPotentialGradients( false ); + + } + + void residual::setdMicroFlowdc( ){ + /*! + * Set the derivative of the micro flow potential w.r.t. the cohesion + */ + + setFlowPotentialGradients( false ); + + } + + void residual::setdMicroGradientFlowdc( ){ + /*! + * Set the derivative of the micro gradient flow potential w.r.t. the cohesion + */ + + setFlowPotentialGradients( false ); + + } + + void residual::setdMacroFlowdDrivingStress( ){ + /*! + * Set the derivative of the macro flow potential w.r.t. the macro driving stress + */ + + setFlowPotentialGradients( false ); + + } + + void residual::setdMicroFlowdDrivingStress( ){ + /*! + * Set the derivative of the micro flow potential w.r.t. the micro driving stress + */ + + setFlowPotentialGradients( false ); + + } + + void residual::setdMicroGradientFlowdDrivingStress( ){ + /*! + * Set the derivative of the micro gradient flow potential w.r.t. the micro-gradient driving stress + */ + + setFlowPotentialGradients( false ); + + } + + void residual::setPreviousdMacroFlowdc( ){ + /*! + * Set the previous derivative of the macro flow potential w.r.t. the cohesion + */ + + setFlowPotentialGradients( true ); + + } + + void residual::setPreviousdMicroFlowdc( ){ + /*! + * Set the previous derivative of the micro flow potential w.r.t. the cohesion + */ + + setFlowPotentialGradients( true ); + + } + + void residual::setPreviousdMicroGradientFlowdc( ){ + /*! + * Set the previous derivative of the micro gradient flow potential w.r.t. the cohesion + */ + + setFlowPotentialGradients( true ); + + } + + void residual::setPreviousdMacroFlowdDrivingStress( ){ + /*! + * Set the previous derivative of the macro flow potential w.r.t. the macro driving stress + */ + + setFlowPotentialGradients( true ); + + } + + void residual::setPreviousdMicroFlowdDrivingStress( ){ + /*! + * Set the previous derivative of the micro flow potential w.r.t. the micro driving stress + */ + + setFlowPotentialGradients( true ); + + } + + void residual::setPreviousdMicroGradientFlowdDrivingStress( ){ + /*! + * Set the previous derivative of the micro gradient flow potential w.r.t. the micro-gradient driving stress + */ + + setFlowPotentialGradients( true ); + + } + + void residual::setFlowPotentialGradients( const bool isPrevious ){ + /*! + * Set the gradients of the flow potential + * + * \param &isPrevious: Flag for whether to set the current (false) or previous (true) gradients + */ + + const floatType *macroCohesion; + + const floatType *microCohesion; + + const floatVector *microGradientCohesion; + + const floatVector *macroDrivingStress; + + const floatVector *microDrivingStress; + + const floatVector *microGradientDrivingStress; + + floatVector precedingDeformationGradient; + + const floatVector *macroFlowParameters = get_macroFlowParameters( ); + + const floatVector *microFlowParameters = get_microFlowParameters( ); + + const floatVector *microGradientFlowParameters = get_microGradientFlowParameters( ); + + if ( isPrevious ){ + + precedingDeformationGradient = hydra->getPreviousPrecedingConfiguration( *getPlasticConfigurationIndex( ) ); + + macroCohesion = get_previousMacroCohesion( ); + + microCohesion = get_previousMicroCohesion( ); + + microGradientCohesion = get_previousMicroGradientCohesion( ); + + macroDrivingStress = get_previousMacroDrivingStress( ); + + microDrivingStress = get_previousSymmetricMicroDrivingStress( ); + + microGradientDrivingStress = get_previousHigherOrderDrivingStress( ); + + } + else{ + + precedingDeformationGradient = hydra->getPrecedingConfiguration( *getPlasticConfigurationIndex( ) ); + + macroCohesion = get_macroCohesion( ); + + microCohesion = get_microCohesion( ); + + microGradientCohesion = get_microGradientCohesion( ); + + macroDrivingStress = get_macroDrivingStress( ); + + microDrivingStress = get_symmetricMicroDrivingStress( ); + + microGradientDrivingStress = get_higherOrderDrivingStress( ); + + } + + floatType tempYield; + + floatVector tempVectorYield; + + floatType dMacroFlowdCohesion, dMicroFlowdCohesion; + + floatVector dMacroFlowdDrivingStress, dMicroFlowdDrivingStress, + dMacroFlowdPrecedingF, dMicroFlowdPrecedingF; + + floatMatrix dMicroGradientFlowdDrivingStress, dMicroGradientFlowdCohesion, dMicroGradientFlowdPrecedingF; + + TARDIGRADE_ERROR_TOOLS_CATCH( computeSecondOrderDruckerPragerYieldEquation( *macroDrivingStress, *macroCohesion, precedingDeformationGradient, + ( *macroFlowParameters )[ 0 ], ( *macroFlowParameters )[ 1 ], + tempYield, dMacroFlowdDrivingStress, dMacroFlowdCohesion, dMacroFlowdPrecedingF ) ); + + TARDIGRADE_ERROR_TOOLS_CATCH( computeSecondOrderDruckerPragerYieldEquation( *microDrivingStress, *microCohesion, precedingDeformationGradient, + ( *microFlowParameters )[ 0 ], ( *microFlowParameters )[ 1 ], + tempYield, dMicroFlowdDrivingStress, dMicroFlowdCohesion, dMicroFlowdPrecedingF ) ); + + TARDIGRADE_ERROR_TOOLS_CATCH( computeHigherOrderDruckerPragerYieldEquation( *microGradientDrivingStress, *microGradientCohesion, precedingDeformationGradient, + ( *microGradientFlowParameters )[ 0 ], ( *microGradientFlowParameters )[ 1 ], + tempVectorYield, dMicroGradientFlowdDrivingStress, dMicroGradientFlowdCohesion, dMicroGradientFlowdPrecedingF ) ); + + if ( isPrevious ){ + + set_previousdMacroFlowdc( dMacroFlowdCohesion ); + + set_previousdMicroFlowdc( dMicroFlowdCohesion ); + + set_previousdMicroGradientFlowdc( dMicroGradientFlowdCohesion ); + + set_previousdMacroFlowdDrivingStress( dMacroFlowdDrivingStress ); + + set_previousdMicroFlowdDrivingStress( dMicroFlowdDrivingStress ); + + set_previousdMicroGradientFlowdDrivingStress( dMicroGradientFlowdDrivingStress ); + + } + else{ + + set_dMacroFlowdc( dMacroFlowdCohesion ); + + set_dMicroFlowdc( dMicroFlowdCohesion ); + + set_dMicroGradientFlowdc( dMicroGradientFlowdCohesion ); + + set_dMacroFlowdDrivingStress( dMacroFlowdDrivingStress ); + + set_dMicroFlowdDrivingStress( dMicroFlowdDrivingStress ); + + set_dMicroGradientFlowdDrivingStress( dMicroGradientFlowdDrivingStress ); + + } + + } + void residual::setMacroCohesion( ){ /*! * Set the macro cohesion diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h index 5cee432b..35d939ef 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h @@ -110,7 +110,7 @@ namespace tardigradeHydra{ variableVector &yieldValue, variableMatrix &dFdStress, variableMatrix &dFdc, variableMatrix &dFdPrecedingF ); - void computeHigherOrderDruckerPragerYieldEquation( const variableVector &stressMEasure, + void computeHigherOrderDruckerPragerYieldEquation( const variableVector &stressMeasure, const variableVector &cohesion, const variableVector &precedingDeformationGradient, const parameterType &frictionAngle, const parameterType &beta, @@ -274,6 +274,32 @@ namespace tardigradeHydra{ virtual void setCohesionsJacobians( const bool isPrevious ); + virtual void setdMacroFlowdc( ); + + virtual void setdMicroFlowdc( ); + + virtual void setdMicroGradientFlowdc( ); + + virtual void setdMacroFlowdDrivingStress( ); + + virtual void setdMicroFlowdDrivingStress( ); + + virtual void setdMicroGradientFlowdDrivingStress( ); + + virtual void setPreviousdMacroFlowdc( ); + + virtual void setPreviousdMicroFlowdc( ); + + virtual void setPreviousdMicroGradientFlowdc( ); + + virtual void setPreviousdMacroFlowdDrivingStress( ); + + virtual void setPreviousdMicroFlowdDrivingStress( ); + + virtual void setPreviousdMicroGradientFlowdDrivingStress( ); + + virtual void setFlowPotentialGradients( const bool isPrevious ); + private: unsigned int _plasticConfigurationIndex; //! The index of the plastic configuration @@ -392,6 +418,30 @@ namespace tardigradeHydra{ TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousdMicroGradientCohesiondStateVariables, floatMatrix, setPreviousdMicroGradientCohesiondStateVariables ) + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, dMacroFlowdc, floatType, setdMacroFlowdc ) + + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, dMacroFlowdDrivingStress, floatVector, setdMacroFlowdDrivingStress ) + + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, dMicroFlowdc, floatType, setdMicroFlowdc ) + + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, dMicroFlowdDrivingStress, floatVector, setdMicroFlowdDrivingStress ) + + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, dMicroGradientFlowdc, floatMatrix, setdMicroGradientFlowdc ) + + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, dMicroGradientFlowdDrivingStress, floatMatrix, setdMicroGradientFlowdDrivingStress ) + + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousdMacroFlowdc, floatType, setPreviousdMacroFlowdc ) + + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousdMacroFlowdDrivingStress, floatVector, setPreviousdMacroFlowdDrivingStress ) + + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousdMicroFlowdc, floatType, setPreviousdMicroFlowdc ) + + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousdMicroFlowdDrivingStress, floatVector, setPreviousdMicroFlowdDrivingStress ) + + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousdMicroGradientFlowdc, floatMatrix, setPreviousdMicroGradientFlowdc ) + + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousdMicroGradientFlowdDrivingStress, floatMatrix, setPreviousdMicroGradientFlowdDrivingStress ) + }; } diff --git a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index 6bc55aef..097532ea 100644 --- a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -1485,7 +1485,7 @@ BOOST_AUTO_TEST_CASE( test_setDrivingStresses ){ BOOST_AUTO_TEST_CASE( test_setCohesion ){ /*! - * Test of the extraction of the parameters + * Test setting the cohesion values */ // Form a hydra class and residual class @@ -1625,12 +1625,12 @@ BOOST_AUTO_TEST_CASE( test_setCohesion ){ floatVector plasticParameters = { 2, 0.53895133, 0.37172145, 2, 0.37773052, 0.92739145, 2, 0.53186824, 0.75454313, - 3, 0.95338442, 0.74042148, 0.09916127, - 3, 0.38093104, 0.49241325, 0.46187452, - 5, 0.82121039, 0.90566759, 0.50466975, 0.04830311, 0.85951495, - 3, 0.01166325, 0.05331896, 0.28081774, - 3, 0.32982199, 0.60161431, 0.33157768, - 5, 0.58881096, 0.11473813, 0.58001078, 0.83382529, 0.3260217 }; + 2, 0.95338442, 0.74042148, + 2, 0.38093104, 0.49241325, + 2, 0.82121039, 0.90566759, + 2, 0.01166325, 0.05331896, + 2, 0.32982199, 0.60161431, + 2, 0.58881096, 0.11473813 }; private: @@ -2256,3 +2256,344 @@ BOOST_AUTO_TEST_CASE( testComputeHigherOrderDruckerPragerYieldEquation ){ } } + +BOOST_AUTO_TEST_CASE( test_setFlowDerivatives ){ + /*! + * Test setting the cohesion values + */ + + // Form a hydra class and residual class + floatType time = 1.23; + + floatType deltaTime = 2.34; + + floatType temperature = 3.45; + + floatType previousTemperature = 4.56; + + variableVector deformationGradient = { 1.00157757, 0.00159138, 0.00672005, + 0.01747159, 1.01122277, 0.00555118, + 0.01112217, -0.00885205, 0.99308943 }; + + variableVector microDeformation = { 0.99460588, -0.0078411 , 0.01145249, + -0.00307139, 0.97798389, -0.00509779, + 0.01189977, -0.01587541, 0.98377259 }; + + variableVector gradientMicroDeformation = { -1.35868385e-02, -1.03142977e-02, 6.54880619e-03, -2.03947530e-02, + -3.31494137e-03, -3.45686183e-03, -3.15745117e-04, -3.70848549e-03, + -9.38693885e-03, -3.68243465e-03, 1.96694582e-02, 2.22080009e-02, + 9.18337942e-05, 6.19764759e-03, -1.92190802e-02, -9.13572591e-03, + -4.25868940e-03, 1.83154579e-02, -1.24772317e-02, -8.48286787e-04, + 2.42779893e-02, 9.74255963e-04, 5.64472629e-03, -1.89685667e-02, + 1.63170400e-02, 5.15300642e-03, 2.25340032e-03 }; + + floatVector previousDeformationGradient = { 9.94656270e-01, 4.82152400e-02, 3.31984800e-02, 2.81918700e-02, + 1.02086536e+00, -1.77592100e-02, -2.24798000e-03, -1.28410000e-04, + 9.77165250e-01 }; + + floatVector previousMicroDeformation = { 0.96917405, -0.01777599, 0.00870406, -0.02163002, 0.9998683 , + -0.01669352, 0.03355217, 0.04427456, 1.01778466 }; + + floatVector previousGradientMicroDeformation = { 0.05043761, 0.02160516, -0.0565408 , 0.01218304, -0.05851034, + 0.00485749, -0.00962607, -0.03455912, 0.04490067, 0.01552915, + -0.02878364, 0.00595866, 0.04750406, -0.02377005, -0.05041534, + -0.02922214, 0.06280788, 0.02850865, -0.00226005, 0.0146049 , + 0.01560184, 0.03224767, 0.05822091, -0.05294424, -0.03518206, + 0.01831308, 0.03774438 }; + + floatVector previousStateVariables = {-0.02495446, -0.00169657, 0.04855598, 0.00194851, 0.01128945, + -0.03793713, 0.03263408, 0.01030601, 0.0045068 , -0.01572362, + -0.01958792, -0.00829778, 0.01813008, 0.03754568, 0.00104223, + 0.01693138, 0.00859366, 0.01249035, 0.01746891, 0.03423424, + -0.0416805 , 0.02636828, -0.02563336, -0.0305777 , 0.0072457 , + -0.04042875, 0.03853268, 0.0127249 , 0.02234164, -0.04838708, + 0.00944319, 0.00567852, -0.03410404, -0.03469295, 0.01955295, + -0.01812336, 0.01919703, 0.00543832, -0.01110494, 0.04251325, + 0.034167 , -0.01426024, -0.04564085, -0.01952319, -0.01018143, + 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0 }; + + floatVector parameters = { 2, 0.53895133, 0.37172145, + 2, 0.37773052, 0.92739145, + 2, 0.53186824, 0.75454313, + 2, 0.95338442, 0.74042148, + 2, 0.38093104, 0.49241325, + 2, 0.82121039, 0.90566759, + 2, 0.01166325, 0.05331896, + 2, 0.32982199, 0.60161431, + 2, 0.58881096, 0.11473813 }; + + floatVector unknownVector = { 4.31202550e-01, -1.78960429e-01, -6.17986089e-01, 9.34988614e-01, + 3.01500733e-01, 7.30919703e-01, -9.49515284e-01, -4.66188370e-01, + 4.14220065e-03, -8.65102730e-01, 9.86066522e-01, -5.27075208e-01, + -2.51415635e-01, -5.71976170e-01, -7.89108268e-01, -5.35040429e-01, + -3.98779729e-01, 2.68884536e-01, -4.37530437e-01, -2.75446478e-01, + -9.88114313e-01, -2.68561748e-01, 6.77719634e-02, -6.75968326e-01, + 1.94866217e-01, -4.13695063e-01, 2.64100990e-01, -9.47606789e-01, + 7.75186921e-01, -9.67762739e-01, -7.46083938e-01, 5.54324923e-01, + -9.08209536e-01, 4.21997387e-01, 9.42092281e-01, 7.43365866e-01, + 4.20323303e-01, 9.17019486e-01, -1.40373324e-01, 7.45757829e-01, + -2.88084664e-01, 8.59527306e-01, -7.02444688e-01, 8.80058030e-01, + 6.65432395e-01, 1.01730274e+00, -1.88038495e-02, 4.82434492e-03, + -2.41803760e-02, 1.01105922e+00, -2.46131243e-02, -2.07588861e-02, + -1.37250795e-02, 1.01875623e+00, 9.93178816e-01, 1.99799676e-03, + 3.40516069e-03, -1.37268320e-02, 1.00360734e+00, 8.04758975e-03, + -1.00877303e-02, -4.06865705e-03, 9.97654446e-01, 2.16175331e-02, + 4.37468737e-03, 2.24126186e-02, 2.80173769e-03, 2.80710425e-05, + -2.48233895e-02, -9.55547806e-04, 2.13727499e-02, -1.50817155e-02, + -2.23954433e-02, -4.66105533e-03, -6.38017597e-03, 1.78576529e-02, + -2.36694442e-02, 2.10074615e-02, 9.04514995e-03, 2.02112997e-02, + 5.37645354e-03, 1.55976656e-02, -8.22280632e-03, -7.52168860e-03, + -5.50628848e-03, 1.27398541e-02, -6.53544128e-03, -1.28890097e-02, + 2.18834178e-02, 2.04005542e-02, + 0.01, 0.02, 0.03, 0.04, 0.05, + 0.06, 0.07, 0.08, 0.09, 0.10 }; + + unsigned int numConfigurations = 2; + + unsigned int numNonLinearSolveStateVariables = 10; + + std::vector< unsigned int > stateVariableIndices = { 0, 1, 2, 3, 4, 5, 6, 7, 8, 9 }; + + unsigned int dimension = 3; + + unsigned int configuration_unknown_count = 45; + + floatType tolr = 1e-2; + + floatType tola = 1e-3; + + unsigned int maxIterations = 24; + + unsigned int maxLSIterations = 45; + + floatType lsAlpha = 2.3; + + class stressMock : public tardigradeHydra::residualBaseMicromorphic{ + + public: + + using tardigradeHydra::residualBaseMicromorphic::residualBaseMicromorphic; + + }; + + class residualMock : public tardigradeHydra::micromorphicDruckerPragerPlasticity::residual{ + + public: + + using tardigradeHydra::micromorphicDruckerPragerPlasticity::residual::residual; + + floatType macroCohesion = 1.23; + + floatType microCohesion = 2.34; + + floatVector microGradientCohesion = { 3.4, 4.5, 5.6 }; + + floatVector macroDrivingStress = { 1, 2, 3, 4, 5, 6, 7, 8, 9 }; + + floatVector microDrivingStress = { 10, 11, 12, 13, 14, 15, 16, 17, 18 }; + + floatVector microGradientDrivingStress = { 19, 20, 21, 22, 23, 24, 25, 26, 27, + 28, 29, 30, 31, 32, 33, 34, 35, 36, + 37, 38, 39, 40, 41, 42, 43, 44, 45 }; + + floatType previousMacroCohesion = -1.23; + + floatType previousMicroCohesion = -2.34; + + floatVector previousMicroGradientCohesion = { -3.4, -4.5, -5.6 }; + + floatVector previousMacroDrivingStress = { -1, -2, -3, -4, -5, -6, -7, -8, -9 }; + + floatVector previousMicroDrivingStress = { -10, -11, -12, -13, -14, -15, -16, -17, -18 }; + + floatVector previousMicroGradientDrivingStress = { -19, -20, -21, -22, -23, -24, -25, -26, -27, + -28, -29, -30, -31, -32, -33, -34, -35, -36, + -37, -38, -39, -40, -41, -42, -43, -44, -45 }; + + floatVector plasticParameters = { 2, 0.53895133, 0.37172145, + 2, 0.37773052, 0.92739145, + 2, 0.53186824, 0.75454313, + 2, 0.95338442, 0.74042148, + 2, 0.38093104, 0.49241325, + 2, 0.82121039, 0.90566759, + 2, 0.01166325, 0.05331896, + 2, 0.32982199, 0.60161431, + 2, 0.58881096, 0.11473813 }; + + virtual void setCohesions( const bool isPrevious ) override{ + + if ( isPrevious ){ + + set_previousMacroCohesion( previousMacroCohesion ); + + set_previousMicroCohesion( previousMicroCohesion ); + + set_previousMicroGradientCohesion( previousMicroGradientCohesion ); + + } + else{ + + set_macroCohesion( macroCohesion ); + + set_microCohesion( microCohesion ); + + set_microGradientCohesion( microGradientCohesion ); + + } + + } + + virtual void setDrivingStresses( const bool isPrevious ) override{ + + if ( isPrevious ){ + + set_previousMacroDrivingStress( previousMacroDrivingStress ); + + set_previousSymmetricMicroDrivingStress( previousMicroDrivingStress ); + + set_previousHigherOrderDrivingStress( previousMicroGradientDrivingStress ); + + } + else{ + + set_macroDrivingStress( macroDrivingStress ); + + set_symmetricMicroDrivingStress( microDrivingStress ); + + set_higherOrderDrivingStress( microGradientDrivingStress ); + + } + + } + + }; + + class hydraBaseMicromorphicMock : public tardigradeHydra::hydraBaseMicromorphic{ + + public: + + using tardigradeHydra::hydraBaseMicromorphic::hydraBaseMicromorphic; + + stressMock elasticity; + + residualMock plasticity; + + std::vector< unsigned int > stateVariableIndices = { 0, 1, 2, 3, 4, 5, 6, 7, 8, 9 }; + + floatVector plasticParameters = { 2, 0.53895133, 0.37172145, + 2, 0.37773052, 0.92739145, + 2, 0.53186824, 0.75454313, + 2, 0.95338442, 0.74042148, + 2, 0.38093104, 0.49241325, + 2, 0.82121039, 0.90566759, + 2, 0.01166325, 0.05331896, + 2, 0.32982199, 0.60161431, + 2, 0.58881096, 0.11473813 }; + + private: + + using tardigradeHydra::hydraBaseMicromorphic::setResidualClasses; + + virtual void setResidualClasses( ) override{ + + std::vector< tardigradeHydra::residualBase* > residuals( 2 ); + + elasticity = stressMock( this, 45 ); + + plasticity = residualMock( this, 55, 1, stateVariableIndices, plasticParameters ); + + residuals[ 0 ] = &elasticity; + + residuals[ 1 ] = &plasticity; + + setResidualClasses( residuals ); + + } + + }; + + hydraBaseMicromorphicMock hydra( time, deltaTime, temperature, previousTemperature, deformationGradient, previousDeformationGradient, + microDeformation, previousMicroDeformation, gradientMicroDeformation, previousGradientMicroDeformation, + previousStateVariables, parameters, + numConfigurations, numNonLinearSolveStateVariables, + dimension, configuration_unknown_count, + tolr, tola, maxIterations, maxLSIterations, lsAlpha ); + + tardigradeHydra::unit_test::hydraBaseTester::updateUnknownVector( hydra, unknownVector ); + + residualMock R( &hydra, 55, 1, stateVariableIndices, parameters ); + + floatType tempYield, dMacroFlowdCohesion, dMicroFlowdCohesion; + + floatVector tempVectorYield, dMacroFlowdDrivingStress, dMicroFlowdDrivingStress, dMacroFlowdPrecedingF, dMicroFlowdPrecedingF; + + floatMatrix dMicroGradientFlowdDrivingStress, dMicroGradientFlowdCohesion, dMicroGradientFlowdPrecedingF; + + floatType previousdMacroFlowdCohesion, previousdMicroFlowdCohesion; + + floatVector previousdMacroFlowdDrivingStress, previousdMicroFlowdDrivingStress, previousdMacroFlowdPrecedingF, previousdMicroFlowdPrecedingF; + + floatMatrix previousdMicroGradientFlowdDrivingStress, previousdMicroGradientFlowdCohesion, previousdMicroGradientFlowdPrecedingF; + + floatVector Fp = floatVector( unknownVector.begin( ) + configuration_unknown_count, + unknownVector.begin( ) + configuration_unknown_count + 9 ); + + floatVector eye( 9 ); + tardigradeVectorTools::eye( eye ); + + floatVector previousFp = floatVector( previousStateVariables.begin( ), + previousStateVariables.begin( ) + 9 ) + eye; + + floatVector precedingDeformationGradient = tardigradeVectorTools::matrixMultiply( deformationGradient, tardigradeVectorTools::inverse( Fp, 3, 3 ), 3, 3, 3, 3 ); + + floatVector previousPrecedingDeformationGradient = tardigradeVectorTools::matrixMultiply( previousDeformationGradient, tardigradeVectorTools::inverse( previousFp, 3, 3 ), 3, 3, 3, 3 ); + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( R.macroDrivingStress, R.macroCohesion, precedingDeformationGradient, + R.plasticParameters[ 10 ], R.plasticParameters[ 11 ], + tempYield, dMacroFlowdDrivingStress, dMacroFlowdCohesion, dMacroFlowdPrecedingF ); + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( R.microDrivingStress, R.microCohesion, precedingDeformationGradient, + R.plasticParameters[ 13 ], R.plasticParameters[ 14 ], + tempYield, dMicroFlowdDrivingStress, dMicroFlowdCohesion, dMicroFlowdPrecedingF ); + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeHigherOrderDruckerPragerYieldEquation( R.microGradientDrivingStress, R.microGradientCohesion, precedingDeformationGradient, + R.plasticParameters[ 16 ], R.plasticParameters[ 17 ], + tempVectorYield, dMicroGradientFlowdDrivingStress, dMicroGradientFlowdCohesion, dMicroGradientFlowdPrecedingF ); + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( R.previousMacroDrivingStress, R.previousMacroCohesion, previousPrecedingDeformationGradient, + R.plasticParameters[ 10 ], R.plasticParameters[ 11 ], + tempYield, previousdMacroFlowdDrivingStress, previousdMacroFlowdCohesion, previousdMacroFlowdPrecedingF ); + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeSecondOrderDruckerPragerYieldEquation( R.previousMicroDrivingStress, R.previousMicroCohesion, previousPrecedingDeformationGradient, + R.plasticParameters[ 13 ], R.plasticParameters[ 14 ], + tempYield, previousdMicroFlowdDrivingStress, previousdMicroFlowdCohesion, previousdMicroFlowdPrecedingF ); + + tardigradeHydra::micromorphicDruckerPragerPlasticity::computeHigherOrderDruckerPragerYieldEquation( R.previousMicroGradientDrivingStress, R.previousMicroGradientCohesion, previousPrecedingDeformationGradient, + R.plasticParameters[ 16 ], R.plasticParameters[ 17 ], + tempVectorYield, previousdMicroGradientFlowdDrivingStress, previousdMicroGradientFlowdCohesion, previousdMicroGradientFlowdPrecedingF ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMacroFlowdCohesion, *R.get_dMacroFlowdc( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMicroFlowdCohesion, *R.get_dMicroFlowdc( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMicroGradientFlowdCohesion, *R.get_dMicroGradientFlowdc( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMacroFlowdDrivingStress, *R.get_dMacroFlowdDrivingStress( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMicroFlowdDrivingStress, *R.get_dMicroFlowdDrivingStress( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMicroGradientFlowdDrivingStress, *R.get_dMicroGradientFlowdDrivingStress( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousdMacroFlowdCohesion, *R.get_previousdMacroFlowdc( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousdMicroFlowdCohesion, *R.get_previousdMicroFlowdc( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousdMicroGradientFlowdCohesion, *R.get_previousdMicroGradientFlowdc( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousdMacroFlowdDrivingStress, *R.get_previousdMacroFlowdDrivingStress( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousdMicroFlowdDrivingStress, *R.get_previousdMicroFlowdDrivingStress( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousdMicroGradientFlowdDrivingStress, *R.get_previousdMicroGradientFlowdDrivingStress( ) ) ); + +} From 938c29a6961c5e7aeec93884ecc46e72d7146816 Mon Sep 17 00:00:00 2001 From: Nathan Miller Date: Mon, 8 Jan 2024 09:28:41 -0700 Subject: [PATCH 10/15] WIP: Initial commit of the jacobians of the flow directions --- ...draMicromorphicDruckerPragerPlasticity.cpp | 280 ++++++++++++++++++ ...hydraMicromorphicDruckerPragerPlasticity.h | 74 +++++ ...draMicromorphicDruckerPragerPlasticity.cpp | 39 ++- 3 files changed, 387 insertions(+), 6 deletions(-) diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index 0a946f34..73aa8cbc 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -1655,6 +1655,286 @@ namespace tardigradeHydra{ } + void residual::setd2MacroFlowdDrivingStressdStress( ){ + /*! + * Set the Jacobian of the derivative of the macro flow potential w.r.t. the macro driving stress w.r.t. the macro stress + */ + + setFlowPotentialGradientsJacobians( false ); + + } + + void residual::setd2MacroFlowdDrivingStressdF( ){ + /*! + * Set the Jacobian of the derivative of the macro flow potential w.r.t. the macro driving stress w.r.t. the deformation gradient + */ + + setFlowPotentialGradientsJacobians( false ); + + } + + void residual::setd2MacroFlowdDrivingStressdFn( ){ + /*! + * Set the Jacobian of the derivative of the macro flow potential w.r.t. the macro driving stress w.r.t. the sub deformation gradients + */ + + setFlowPotentialGradientsJacobians( false ); + + } + + void residual::setd2MicroFlowdDrivingStressdStress( ){ + /*! + * Set the Jacobian of the derivative of the micro flow potential w.r.t. the micro driving stress w.r.t. the micro stress + */ + + setFlowPotentialGradientsJacobians( false ); + + } + + void residual::setd2MicroFlowdDrivingStressdF( ){ + /*! + * Set the Jacobian of the derivative of the micro flow potential w.r.t. the micro driving stress w.r.t. the deformation gradient + */ + + setFlowPotentialGradientsJacobians( false ); + + } + + void residual::setd2MicroFlowdDrivingStressdFn( ){ + /*! + * Set the Jacobian of the derivative of the micro flow potential w.r.t. the micro driving stress w.r.t. the sub deformation gradients + */ + + setFlowPotentialGradientsJacobians( false ); + + } + + void residual::setd2MicroGradientFlowdDrivingStressdStress( ){ + /*! + * Set the Jacobian of the derivative of the micro gradient flow potential w.r.t. the micro-gradient driving stress w.r.t. the micro gradient stress + */ + + setFlowPotentialGradientsJacobians( false ); + + } + + void residual::setd2MicroGradientFlowdDrivingStressdF( ){ + /*! + * Set the Jacobian of the derivative of the micro gradient flow potential w.r.t. the micro-gradient driving stress w.r.t. the deformation gradient + */ + + setFlowPotentialGradientsJacobians( false ); + + } + + void residual::setd2MicroGradientFlowdDrivingStressdFn( ){ + /*! + * Set the Jacobian of the derivative of the micro gradient flow potential w.r.t. the micro-gradient driving stress w.r.t. the sub deformation gradients + */ + + setFlowPotentialGradientsJacobians( false ); + + } + + void residual::setPreviousd2MacroFlowdDrivingStressdStress( ){ + /*! + * Set the Jacobian of the previous derivative of the macro flow potential w.r.t. the macro driving stress w.r.t. the macro stress + */ + + setFlowPotentialGradientsJacobians( true ); + + } + + void residual::setPreviousd2MacroFlowdDrivingStressdF( ){ + /*! + * Set the Jacobian of the previous derivative of the macro flow potential w.r.t. the macro driving stress w.r.t. the deformation gradient + */ + + setFlowPotentialGradientsJacobians( true ); + + } + + void residual::setPreviousd2MacroFlowdDrivingStressdFn( ){ + /*! + * Set the Jacobian of the previous derivative of the macro flow potential w.r.t. the macro driving stress w.r.t. the sub deformation gradients + */ + + setFlowPotentialGradientsJacobians( true ); + + } + + void residual::setPreviousd2MicroFlowdDrivingStressdStress( ){ + /*! + * Set the Jacobian of the previous derivative of the micro flow potential w.r.t. the micro driving stress w.r.t. the micro stress + */ + + setFlowPotentialGradientsJacobians( true ); + + } + + void residual::setPreviousd2MicroFlowdDrivingStressdF( ){ + /*! + * Set the Jacobian of the previous derivative of the micro flow potential w.r.t. the micro driving stress w.r.t. the deformation gradient + */ + + setFlowPotentialGradientsJacobians( true ); + + } + + void residual::setPreviousd2MicroFlowdDrivingStressdFn( ){ + /*! + * Set the Jacobian of the previous derivative of the micro flow potential w.r.t. the micro driving stress w.r.t. the sub deformation gradients + */ + + setFlowPotentialGradientsJacobians( true ); + + } + + void residual::setPreviousd2MicroGradientFlowdDrivingStressdStress( ){ + /*! + * Set the Jacobian of the previous derivative of the micro gradient flow potential w.r.t. the micro-gradient driving stress w.r.t. the micro gradient stress + */ + + setFlowPotentialGradientsJacobians( true ); + + } + + void residual::setPreviousd2MicroGradientFlowdDrivingStressdF( ){ + /*! + * Set the Jacobian of the previous derivative of the micro gradient flow potential w.r.t. the micro-gradient driving stress w.r.t. the deformation gradient + */ + + setFlowPotentialGradientsJacobians( true ); + + } + + void residual::setPreviousd2MicroGradientFlowdDrivingStressdFn( ){ + /*! + * Set the Jacobian of the previous derivative of the micro gradient flow potential w.r.t. the micro-gradient driving stress w.r.t. the sub deformation gradients + */ + + setFlowPotentialGradientsJacobians( true ); + + } + + void residual::setFlowPotentialGradientsJacobians( const bool isPrevious ){ + /*! + * Set the Jacobians of the flow potential gradients + * + * \param isPrevious: A flag for whether to set the current (false) or previous (true) derivatives + */ + + const floatType *macroCohesion; + + const floatType *microCohesion; + + const floatVector *microGradientCohesion; + + const floatVector *macroDrivingStress; + + const floatVector *microDrivingStress; + + const floatVector *microGradientDrivingStress; + + floatVector precedingDeformationGradient; + + const floatVector *macroFlowParameters = get_macroFlowParameters( ); + + const floatVector *microFlowParameters = get_microFlowParameters( ); + + const floatVector *microGradientFlowParameters = get_microGradientFlowParameters( ); + + if ( isPrevious ){ + + precedingDeformationGradient = hydra->getPreviousPrecedingConfiguration( *getPlasticConfigurationIndex( ) ); + + macroCohesion = get_previousMacroCohesion( ); + + microCohesion = get_previousMicroCohesion( ); + + microGradientCohesion = get_previousMicroGradientCohesion( ); + + macroDrivingStress = get_previousMacroDrivingStress( ); + + microDrivingStress = get_previousSymmetricMicroDrivingStress( ); + + microGradientDrivingStress = get_previousHigherOrderDrivingStress( ); + + } + else{ + + precedingDeformationGradient = hydra->getPrecedingConfiguration( *getPlasticConfigurationIndex( ) ); + + macroCohesion = get_macroCohesion( ); + + microCohesion = get_microCohesion( ); + + microGradientCohesion = get_microGradientCohesion( ); + + macroDrivingStress = get_macroDrivingStress( ); + + microDrivingStress = get_symmetricMicroDrivingStress( ); + + microGradientDrivingStress = get_higherOrderDrivingStress( ); + + } + + floatType tempYield; + + floatVector tempVectorYield; + + floatType dMacroFlowdCohesion, dMicroFlowdCohesion; + + floatVector dMacroFlowdDrivingStress, dMicroFlowdDrivingStress, + dMacroFlowdPrecedingF, dMicroFlowdPrecedingF; + + floatMatrix dMicroGradientFlowdDrivingStress, dMicroGradientFlowdCohesion, dMicroGradientFlowdPrecedingF; + + TARDIGRADE_ERROR_TOOLS_CATCH( computeSecondOrderDruckerPragerYieldEquation( *macroDrivingStress, *macroCohesion, precedingDeformationGradient, + ( *macroFlowParameters )[ 0 ], ( *macroFlowParameters )[ 1 ], + tempYield, dMacroFlowdDrivingStress, dMacroFlowdCohesion, dMacroFlowdPrecedingF ) ); + + TARDIGRADE_ERROR_TOOLS_CATCH( computeSecondOrderDruckerPragerYieldEquation( *microDrivingStress, *microCohesion, precedingDeformationGradient, + ( *microFlowParameters )[ 0 ], ( *microFlowParameters )[ 1 ], + tempYield, dMicroFlowdDrivingStress, dMicroFlowdCohesion, dMicroFlowdPrecedingF ) ); + + TARDIGRADE_ERROR_TOOLS_CATCH( computeHigherOrderDruckerPragerYieldEquation( *microGradientDrivingStress, *microGradientCohesion, precedingDeformationGradient, + ( *microGradientFlowParameters )[ 0 ], ( *microGradientFlowParameters )[ 1 ], + tempVectorYield, dMicroGradientFlowdDrivingStress, dMicroGradientFlowdCohesion, dMicroGradientFlowdPrecedingF ) ); + + if ( isPrevious ){ + + set_previousdMacroFlowdc( dMacroFlowdCohesion ); + + set_previousdMicroFlowdc( dMicroFlowdCohesion ); + + set_previousdMicroGradientFlowdc( dMicroGradientFlowdCohesion ); + + set_previousdMacroFlowdDrivingStress( dMacroFlowdDrivingStress ); + + set_previousdMicroFlowdDrivingStress( dMicroFlowdDrivingStress ); + + set_previousdMicroGradientFlowdDrivingStress( dMicroGradientFlowdDrivingStress ); + + } + else{ + + set_dMacroFlowdc( dMacroFlowdCohesion ); + + set_dMicroFlowdc( dMicroFlowdCohesion ); + + set_dMicroGradientFlowdc( dMicroGradientFlowdCohesion ); + + set_dMacroFlowdDrivingStress( dMacroFlowdDrivingStress ); + + set_dMicroFlowdDrivingStress( dMicroFlowdDrivingStress ); + + set_dMicroGradientFlowdDrivingStress( dMicroGradientFlowdDrivingStress ); + + } + + } + void residual::setMacroCohesion( ){ /*! * Set the macro cohesion diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h index 35d939ef..cc35289f 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h @@ -300,6 +300,44 @@ namespace tardigradeHydra{ virtual void setFlowPotentialGradients( const bool isPrevious ); + virtual void setd2MacroFlowdDrivingStressdStress( ); + + virtual void setd2MacroFlowdDrivingStressdF( ); + + virtual void setd2MacroFlowdDrivingStressdFn( ); + + virtual void setd2MicroFlowdDrivingStressdStress( ); + + virtual void setd2MicroFlowdDrivingStressdF( ); + + virtual void setd2MicroFlowdDrivingStressdFn( ); + + virtual void setd2MicroGradientFlowdDrivingStressdStress( ); + + virtual void setd2MicroGradientFlowdDrivingStressdF( ); + + virtual void setd2MicroGradientFlowdDrivingStressdFn( ); + + virtual void setPreviousd2MacroFlowdDrivingStressdStress( ); + + virtual void setPreviousd2MacroFlowdDrivingStressdF( ); + + virtual void setPreviousd2MacroFlowdDrivingStressdFn( ); + + virtual void setPreviousd2MicroFlowdDrivingStressdStress( ); + + virtual void setPreviousd2MicroFlowdDrivingStressdF( ); + + virtual void setPreviousd2MicroFlowdDrivingStressdFn( ); + + virtual void setPreviousd2MicroGradientFlowdDrivingStressdStress( ); + + virtual void setPreviousd2MicroGradientFlowdDrivingStressdF( ); + + virtual void setPreviousd2MicroGradientFlowdDrivingStressdFn( ); + + virtual void setFlowPotentialGradientsJacobians( const bool isPrevious ); + private: unsigned int _plasticConfigurationIndex; //! The index of the plastic configuration @@ -422,26 +460,62 @@ namespace tardigradeHydra{ TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, dMacroFlowdDrivingStress, floatVector, setdMacroFlowdDrivingStress ) + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, d2MacroFlowdDrivingStressdStress, floatMatrix, setd2MacroFlowdDrivingStressdStress ) + + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, d2MacroFlowdDrivingStressdF, floatMatrix, setd2MacroFlowdDrivingStressdF ) + + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, d2MacroFlowdDrivingStressdFn, floatMatrix, setd2MacroFlowdDrivingStressdFn ) + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, dMicroFlowdc, floatType, setdMicroFlowdc ) TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, dMicroFlowdDrivingStress, floatVector, setdMicroFlowdDrivingStress ) + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, d2MicroFlowdDrivingStressdStress, floatMatrix, setd2MicroFlowdDrivingStressdStress ) + + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, d2MicroFlowdDrivingStressdF, floatMatrix, setd2MicroFlowdDrivingStressdF ) + + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, d2MicroFlowdDrivingStressdFn, floatMatrix, setd2MicroFlowdDrivingStressdFn ) + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, dMicroGradientFlowdc, floatMatrix, setdMicroGradientFlowdc ) TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, dMicroGradientFlowdDrivingStress, floatMatrix, setdMicroGradientFlowdDrivingStress ) + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, d2MicroGradientFlowdDrivingStressdStress, floatMatrix, setd2MicroGradientFlowdDrivingStressdStress ) + + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, d2MicroGradientFlowdDrivingStressdF, floatMatrix, setd2MicroGradientFlowdDrivingStressdF ) + + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, d2MicroGradientFlowdDrivingStressdFn, floatMatrix, setd2MicroGradientFlowdDrivingStressdFn ) + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousdMacroFlowdc, floatType, setPreviousdMacroFlowdc ) TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousdMacroFlowdDrivingStress, floatVector, setPreviousdMacroFlowdDrivingStress ) + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousd2MacroFlowdDrivingStressdStress, floatMatrix, setPreviousd2MacroFlowdDrivingStressdStress ) + + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousd2MacroFlowdDrivingStressdF, floatMatrix, setPreviousd2MacroFlowdDrivingStressdF ) + + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousd2MacroFlowdDrivingStressdFn, floatMatrix, setPreviousd2MacroFlowdDrivingStressdFn ) + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousdMicroFlowdc, floatType, setPreviousdMicroFlowdc ) TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousdMicroFlowdDrivingStress, floatVector, setPreviousdMicroFlowdDrivingStress ) + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousd2MicroFlowdDrivingStressdStress, floatMatrix, setPreviousd2MicroFlowdDrivingStressdStress ) + + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousd2MicroFlowdDrivingStressdF, floatMatrix, setPreviousd2MicroFlowdDrivingStressdF ) + + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousd2MicroFlowdDrivingStressdFn, floatMatrix, setPreviousd2MicroFlowdDrivingStressdFn ) + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousdMicroGradientFlowdc, floatMatrix, setPreviousdMicroGradientFlowdc ) TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousdMicroGradientFlowdDrivingStress, floatMatrix, setPreviousdMicroGradientFlowdDrivingStress ) + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousd2MicroGradientFlowdDrivingStressdStress, floatMatrix, setPreviousd2MicroGradientFlowdDrivingStressdStress ) + + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousd2MicroGradientFlowdDrivingStressdF, floatMatrix, setPreviousd2MicroGradientFlowdDrivingStressdF ) + + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousd2MicroGradientFlowdDrivingStressdFn, floatMatrix, setPreviousd2MicroGradientFlowdDrivingStressdFn ) + }; } diff --git a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index 097532ea..d6bfa18e 100644 --- a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -2523,6 +2523,8 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives ){ residualMock R( &hydra, 55, 1, stateVariableIndices, parameters ); + residualMock RJ( &hydra, 55, 1, stateVariableIndices, parameters ); + floatType tempYield, dMacroFlowdCohesion, dMicroFlowdCohesion; floatVector tempVectorYield, dMacroFlowdDrivingStress, dMicroFlowdDrivingStress, dMacroFlowdPrecedingF, dMicroFlowdPrecedingF; @@ -2572,17 +2574,19 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives ){ R.plasticParameters[ 16 ], R.plasticParameters[ 17 ], tempVectorYield, previousdMicroGradientFlowdDrivingStress, previousdMicroGradientFlowdCohesion, previousdMicroGradientFlowdPrecedingF ); - BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMacroFlowdCohesion, *R.get_dMacroFlowdc( ) ) ); + RJ.get_d2MacroFlowdDrivingStressdStress( ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMacroFlowdCohesion, *R.get_dMacroFlowdc( ) ) ); - BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMicroFlowdCohesion, *R.get_dMicroFlowdc( ) ) ); + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMicroFlowdCohesion, *R.get_dMicroFlowdc( ) ) ); - BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMicroGradientFlowdCohesion, *R.get_dMicroGradientFlowdc( ) ) ); + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMicroGradientFlowdCohesion, *R.get_dMicroGradientFlowdc( ) ) ); - BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMacroFlowdDrivingStress, *R.get_dMacroFlowdDrivingStress( ) ) ); + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMacroFlowdDrivingStress, *R.get_dMacroFlowdDrivingStress( ) ) ); - BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMicroFlowdDrivingStress, *R.get_dMicroFlowdDrivingStress( ) ) ); + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMicroFlowdDrivingStress, *R.get_dMicroFlowdDrivingStress( ) ) ); - BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMicroGradientFlowdDrivingStress, *R.get_dMicroGradientFlowdDrivingStress( ) ) ); + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMicroGradientFlowdDrivingStress, *R.get_dMicroGradientFlowdDrivingStress( ) ) ); BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousdMacroFlowdCohesion, *R.get_previousdMacroFlowdc( ) ) ); @@ -2596,4 +2600,27 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives ){ BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousdMicroGradientFlowdDrivingStress, *R.get_previousdMicroGradientFlowdDrivingStress( ) ) ); + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMacroFlowdCohesion, *RJ.get_dMacroFlowdc( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMicroFlowdCohesion, *RJ.get_dMicroFlowdc( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMicroGradientFlowdCohesion, *RJ.get_dMicroGradientFlowdc( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMacroFlowdDrivingStress, *RJ.get_dMacroFlowdDrivingStress( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMicroFlowdDrivingStress, *RJ.get_dMicroFlowdDrivingStress( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( dMicroGradientFlowdDrivingStress, *RJ.get_dMicroGradientFlowdDrivingStress( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousdMacroFlowdCohesion, *RJ.get_previousdMacroFlowdc( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousdMicroFlowdCohesion, *RJ.get_previousdMicroFlowdc( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousdMicroGradientFlowdCohesion, *RJ.get_previousdMicroGradientFlowdc( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousdMacroFlowdDrivingStress, *RJ.get_previousdMacroFlowdDrivingStress( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousdMicroFlowdDrivingStress, *RJ.get_previousdMicroFlowdDrivingStress( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousdMicroGradientFlowdDrivingStress, *RJ.get_previousdMicroGradientFlowdDrivingStress( ) ) ); } From 6f1e820d4929230cd9bf30e29de4f043895860cc Mon Sep 17 00:00:00 2001 From: Nathan Miller Date: Mon, 8 Jan 2024 11:32:05 -0700 Subject: [PATCH 11/15] WIP: Added extraction of Jacobians of the driving stresses --- ...draMicromorphicDruckerPragerPlasticity.cpp | 54 +++++++++++++++++++ ...draMicromorphicDruckerPragerPlasticity.cpp | 6 +++ 2 files changed, 60 insertions(+) diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index 73aa8cbc..4e4ee79d 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -1836,6 +1836,24 @@ namespace tardigradeHydra{ const floatVector *microGradientDrivingStress; + const floatMatrix *dMacroDrivingStressdStress; + + const floatMatrix *dMacroDrivingStressdF; + + const floatMatrix *dMacroDrivingStressdFn; + + const floatMatrix *dMicroDrivingStressdStress; + + const floatMatrix *dMicroDrivingStressdF; + + const floatMatrix *dMicroDrivingStressdFn; + + const floatMatrix *dMicroGradientDrivingStressdStress; + + const floatMatrix *dMicroGradientDrivingStressdF; + + const floatMatrix *dMicroGradientDrivingStressdFn; + floatVector precedingDeformationGradient; const floatVector *macroFlowParameters = get_macroFlowParameters( ); @@ -1854,6 +1872,24 @@ namespace tardigradeHydra{ microGradientCohesion = get_previousMicroGradientCohesion( ); + dMacroDrivingStressdStress = get_previousdMacroDrivingStressdMacroStress( ); + + dMicroDrivingStressdStress = get_previousdSymmetricMicroDrivingStressdMicroStress( ); + + dMicroGradientDrivingStressdStress = get_previousdHigherOrderDrivingStressdHigherOrderStress( ); + + dMacroDrivingStressdF = get_previousdMacroDrivingStressdF( ); + + dMicroDrivingStressdF = get_previousdSymmetricMicroDrivingStressdF( ); + + dMicroGradientDrivingStressdF = get_previousdHigherOrderDrivingStressdF( ); + + dMacroDrivingStressdFn = get_previousdMacroDrivingStressdFn( ); + + dMicroDrivingStressdFn = get_previousdSymmetricMicroDrivingStressdFn( ); + + dMicroGradientDrivingStressdFn = get_previousdHigherOrderDrivingStressdFn( ); + macroDrivingStress = get_previousMacroDrivingStress( ); microDrivingStress = get_previousSymmetricMicroDrivingStress( ); @@ -1871,6 +1907,24 @@ namespace tardigradeHydra{ microGradientCohesion = get_microGradientCohesion( ); + dMacroDrivingStressdStress = get_dMacroDrivingStressdMacroStress( ); + + dMicroDrivingStressdStress = get_dSymmetricMicroDrivingStressdMicroStress( ); + + dMicroGradientDrivingStressdStress = get_dHigherOrderDrivingStressdHigherOrderStress( ); + + dMacroDrivingStressdF = get_dMacroDrivingStressdF( ); + + dMicroDrivingStressdF = get_dSymmetricMicroDrivingStressdF( ); + + dMicroGradientDrivingStressdF = get_dHigherOrderDrivingStressdF( ); + + dMacroDrivingStressdFn = get_dMacroDrivingStressdFn( ); + + dMicroDrivingStressdFn = get_dSymmetricMicroDrivingStressdFn( ); + + dMicroGradientDrivingStressdFn = get_dHigherOrderDrivingStressdFn( ); + macroDrivingStress = get_macroDrivingStress( ); microDrivingStress = get_symmetricMicroDrivingStress( ); diff --git a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index d6bfa18e..68b27810 100644 --- a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -2466,6 +2466,12 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives ){ } + virtual void setDrivingStressesJacobians( const bool isPrevious ) override{ + + setDrivingStresses( isPrevious ); + + } + }; class hydraBaseMicromorphicMock : public tardigradeHydra::hydraBaseMicromorphic{ From 651d14d50148c7ff06c7bc54123e1d7a8df3ab26 Mon Sep 17 00:00:00 2001 From: Nathan Miller Date: Mon, 8 Jan 2024 16:20:31 -0700 Subject: [PATCH 12/15] WIP: End of day check in --- ...draMicromorphicDruckerPragerPlasticity.cpp | 198 ++++++++- ...draMicromorphicDruckerPragerPlasticity.cpp | 383 ++++++++++++++++++ 2 files changed, 566 insertions(+), 15 deletions(-) diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index 4e4ee79d..799a5944 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -1824,6 +1824,8 @@ namespace tardigradeHydra{ * \param isPrevious: A flag for whether to set the current (false) or previous (true) derivatives */ + std::cout << "entering setFlowPotentialGradientsJacobians\n"; + const floatType *macroCohesion; const floatType *microCohesion; @@ -1856,6 +1858,12 @@ namespace tardigradeHydra{ floatVector precedingDeformationGradient; + floatMatrix dPrecedingFdSubFs; + + const floatMatrix *dF1dF; + + const floatMatrix *dF1dFn; + const floatVector *macroFlowParameters = get_macroFlowParameters( ); const floatVector *microFlowParameters = get_microFlowParameters( ); @@ -1866,11 +1874,17 @@ namespace tardigradeHydra{ precedingDeformationGradient = hydra->getPreviousPrecedingConfiguration( *getPlasticConfigurationIndex( ) ); - macroCohesion = get_previousMacroCohesion( ); + dPrecedingFdSubFs = hydra->getPreviousPrecedingConfigurationJacobian( *getPlasticConfigurationIndex( ) ); - microCohesion = get_previousMicroCohesion( ); + dF1dF = hydra->get_previousdF1dF( ); - microGradientCohesion = get_previousMicroGradientCohesion( ); + dF1dFn = hydra->get_previousdF1dFn( ); + + macroCohesion = get_previousMacroCohesion( ); + + microCohesion = get_previousMicroCohesion( ); + + microGradientCohesion = get_previousMicroGradientCohesion( ); dMacroDrivingStressdStress = get_previousdMacroDrivingStressdMacroStress( ); @@ -1890,22 +1904,28 @@ namespace tardigradeHydra{ dMicroGradientDrivingStressdFn = get_previousdHigherOrderDrivingStressdFn( ); - macroDrivingStress = get_previousMacroDrivingStress( ); + macroDrivingStress = get_previousMacroDrivingStress( ); - microDrivingStress = get_previousSymmetricMicroDrivingStress( ); + microDrivingStress = get_previousSymmetricMicroDrivingStress( ); - microGradientDrivingStress = get_previousHigherOrderDrivingStress( ); + microGradientDrivingStress = get_previousHigherOrderDrivingStress( ); } else{ precedingDeformationGradient = hydra->getPrecedingConfiguration( *getPlasticConfigurationIndex( ) ); - macroCohesion = get_macroCohesion( ); + dPrecedingFdSubFs = hydra->getPrecedingConfigurationJacobian( *getPlasticConfigurationIndex( ) ); - microCohesion = get_microCohesion( ); + dF1dF = hydra->get_dF1dF( ); - microGradientCohesion = get_microGradientCohesion( ); + dF1dFn = hydra->get_dF1dFn( ); + + macroCohesion = get_macroCohesion( ); + + microCohesion = get_microCohesion( ); + + microGradientCohesion = get_microGradientCohesion( ); dMacroDrivingStressdStress = get_dMacroDrivingStressdMacroStress( ); @@ -1925,11 +1945,43 @@ namespace tardigradeHydra{ dMicroGradientDrivingStressdFn = get_dHigherOrderDrivingStressdFn( ); - macroDrivingStress = get_macroDrivingStress( ); + macroDrivingStress = get_macroDrivingStress( ); - microDrivingStress = get_symmetricMicroDrivingStress( ); + microDrivingStress = get_symmetricMicroDrivingStress( ); - microGradientDrivingStress = get_higherOrderDrivingStress( ); + microGradientDrivingStress = get_higherOrderDrivingStress( ); + + } + + // Construct the derivatives of the preceding F + + floatMatrix dPrecedingFdF( precedingDeformationGradient.size( ), floatVector( hydra->getDeformationGradient( )->size( ), 0 ) ); + + floatMatrix dPrecedingFdFn( precedingDeformationGradient.size( ), floatVector( ( ( *hydra->getNumConfigurations( ) ) - 1 ) * hydra->getDeformationGradient( )->size( ), 0 ) ); + + for ( unsigned int i = 0; i < hydra->getDeformationGradient( )->size( ); i++ ){ + + for ( unsigned int j = 0; j < hydra->getDeformationGradient( )->size( ); j++ ){ + + for ( unsigned int k = 0; k < hydra->getDeformationGradient( )->size( ); k++ ){ + + dPrecedingFdF[ i ][ j ] += dPrecedingFdSubFs[ i ][ k ] * ( *dF1dF )[ k ][ j ]; + + } + + } + + for ( unsigned int j = 0; j < ( ( *hydra->getNumConfigurations( ) ) - 1 ) * hydra->getDeformationGradient( )->size( ); j++ ){ + + dPrecedingFdFn[ i ][ j ] = dPrecedingFdSubFs[ i ][ hydra->getDeformationGradient( )->size( ) + j ]; + + for ( unsigned int k = 0; k < hydra->getDeformationGradient( )->size( ); k++ ){ + + dPrecedingFdFn[ i ][ j ] += dPrecedingFdSubFs[ i ][ k ] * ( *dF1dFn )[ k ][ j ]; + + } + + } } @@ -1944,17 +1996,107 @@ namespace tardigradeHydra{ floatMatrix dMicroGradientFlowdDrivingStress, dMicroGradientFlowdCohesion, dMicroGradientFlowdPrecedingF; + floatMatrix d2MacroFlowdDrivingStress2, d2MacroFlowdDrivingStressdPrecedingF, + d2MicroFlowdDrivingStress2, d2MicroFlowdDrivingStressdPrecedingF, + d2MicroGradientFlowdDrivingStress2, d2MicroGradientFlowdDrivingStressdPrecedingF; + TARDIGRADE_ERROR_TOOLS_CATCH( computeSecondOrderDruckerPragerYieldEquation( *macroDrivingStress, *macroCohesion, precedingDeformationGradient, ( *macroFlowParameters )[ 0 ], ( *macroFlowParameters )[ 1 ], - tempYield, dMacroFlowdDrivingStress, dMacroFlowdCohesion, dMacroFlowdPrecedingF ) ); + tempYield, dMacroFlowdDrivingStress, dMacroFlowdCohesion, dMacroFlowdPrecedingF, + d2MacroFlowdDrivingStress2, d2MacroFlowdDrivingStressdPrecedingF ) ); TARDIGRADE_ERROR_TOOLS_CATCH( computeSecondOrderDruckerPragerYieldEquation( *microDrivingStress, *microCohesion, precedingDeformationGradient, ( *microFlowParameters )[ 0 ], ( *microFlowParameters )[ 1 ], - tempYield, dMicroFlowdDrivingStress, dMicroFlowdCohesion, dMicroFlowdPrecedingF ) ); + tempYield, dMicroFlowdDrivingStress, dMicroFlowdCohesion, dMicroFlowdPrecedingF, + d2MicroFlowdDrivingStress2, d2MicroFlowdDrivingStressdPrecedingF ) ); TARDIGRADE_ERROR_TOOLS_CATCH( computeHigherOrderDruckerPragerYieldEquation( *microGradientDrivingStress, *microGradientCohesion, precedingDeformationGradient, ( *microGradientFlowParameters )[ 0 ], ( *microGradientFlowParameters )[ 1 ], - tempVectorYield, dMicroGradientFlowdDrivingStress, dMicroGradientFlowdCohesion, dMicroGradientFlowdPrecedingF ) ); + tempVectorYield, dMicroGradientFlowdDrivingStress, dMicroGradientFlowdCohesion, dMicroGradientFlowdPrecedingF, + d2MicroGradientFlowdDrivingStress2, d2MicroGradientFlowdDrivingStressdPrecedingF ) ); + + floatMatrix d2MacroFlowdDrivingStressdMacroStress( dMacroFlowdDrivingStress.size( ), floatVector( macroDrivingStress->size( ) * macroDrivingStress->size( ), 0 ) ); + + floatMatrix d2MicroFlowdDrivingStressdMicroStress( dMicroFlowdDrivingStress.size( ), floatVector( microDrivingStress->size( ) * microDrivingStress->size( ), 0 ) ); + + floatMatrix d2MicroGradientFlowdDrivingStressdMicroGradientStress( dMicroGradientFlowdDrivingStress.size( ), floatVector( microGradientDrivingStress->size( ) * microGradientDrivingStress->size( ), 0 ) ); + + floatMatrix d2MacroFlowdDrivingStressdFn( dMacroFlowdDrivingStress.size( ), floatVector( macroDrivingStress->size( ) * ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ), 0 ) ); + + floatMatrix d2MicroFlowdDrivingStressdFn( dMicroFlowdDrivingStress.size( ), floatVector( microDrivingStress->size( ) * ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ), 0 ) ); + + floatMatrix d2MicroGradientFlowdDrivingStressdFn( dMicroGradientFlowdDrivingStress.size( ), floatVector( microGradientDrivingStress->size( ) * ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ), 0 ) ); + + for ( unsigned int I = 0; I < d2MicroFlowdDrivingStress2.size( ); I++ ){ + + for ( unsigned int J = 0; J < macroDrivingStress->size( ); J++ ){ + + for ( unsigned int K = 0; K < macroDrivingStress->size( ); K++ ){ + + for ( unsigned int L = 0; L < macroDrivingStress->size( ); L++ ){ + + d2MacroFlowdDrivingStressdMacroStress[ I ][ macroDrivingStress->size( ) * J + K ] + += d2MacroFlowdDrivingStress2[ I ][ macroDrivingStress->size( ) * J + L ] * ( *dMacroDrivingStressdStress )[ L ][ K ]; + + d2MicroFlowdDrivingStressdMicroStress[ I ][ microDrivingStress->size( ) * J + K ] + += d2MicroFlowdDrivingStress2[ I ][ microDrivingStress->size( ) * J + L ] * ( *dMicroDrivingStressdStress )[ L ][ K ]; + + } + + for ( unsigned int L = 0; L < ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ); L++ ){ + + d2MacroFlowdDrivingStressdFn[ I ][ ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ) * J + L ] + += d2MacroFlowdDrivingStressdPrecedingF[ I ][ macroDrivingStress->size( ) * J + K ] * dPrecedingFdFn[ K ][ L ] + + d2MacroFlowdDrivingStress2[ I ][ macroDrivingStress->size( ) * J + K ] * ( *dMacroDrivingStressdFn ) [ K ][ L ]; + + d2MicroFlowdDrivingStressdFn[ I ][ ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ) * J + L ] + += d2MicroFlowdDrivingStressdPrecedingF[ I ][ precedingDeformationGradient.size( ) * J + K ] * dPrecedingFdFn[ K ][ L ] + + d2MicroFlowdDrivingStress2[ I ][ microDrivingStress->size( ) * J + K ] * ( *dMicroDrivingStressdFn ) [ K ][ L ]; + + } + + } + + } + + } + + for ( unsigned int I = 0; I < d2MicroGradientFlowdDrivingStress2.size( ); I++ ){ + + for ( unsigned int J = 0; J < microGradientDrivingStress->size( ); J++ ){ + + for ( unsigned int K = 0; K < microGradientDrivingStress->size( ); K++ ){ + + for ( unsigned int L = 0; L < microGradientDrivingStress->size( ); L++ ){ + + d2MicroGradientFlowdDrivingStressdMicroGradientStress[ I ][ microGradientDrivingStress->size( ) * J + K ] + += d2MicroGradientFlowdDrivingStress2[ I ][ microGradientDrivingStress->size( ) * J + L ] * ( *dMicroGradientDrivingStressdStress )[ L ][ K ]; + + } + + for ( unsigned int L = 0; L < ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ); L++ ){ + + d2MicroGradientFlowdDrivingStressdFn[ I ][ ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ) * J + L ] + += d2MicroGradientFlowdDrivingStress2[ I ][ microGradientDrivingStress->size( ) * J + K ] * ( *dMicroGradientDrivingStressdFn ) [ K ][ L ]; + + } + + } + + for ( unsigned int K = 0; K < precedingDeformationGradient.size( ); K++ ){ + + for ( unsigned int L = 0; L < ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ); L++ ){ + + d2MicroGradientFlowdDrivingStressdFn[ I ][ ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ) * J + L ] + += d2MicroGradientFlowdDrivingStressdPrecedingF[ I ][ precedingDeformationGradient.size( ) * J + K ] * dPrecedingFdFn[ K ][ L ]; + + } + + } + + } + + } if ( isPrevious ){ @@ -1970,6 +2112,18 @@ namespace tardigradeHydra{ set_previousdMicroGradientFlowdDrivingStress( dMicroGradientFlowdDrivingStress ); + set_previousd2MacroFlowdDrivingStressdStress( d2MacroFlowdDrivingStressdMacroStress ); + + set_previousd2MicroFlowdDrivingStressdStress( d2MicroFlowdDrivingStressdMicroStress ); + + set_previousd2MicroGradientFlowdDrivingStressdStress( d2MicroGradientFlowdDrivingStressdMicroGradientStress ); + + set_previousd2MacroFlowdDrivingStressdFn( d2MacroFlowdDrivingStressdFn ); + + set_previousd2MicroFlowdDrivingStressdFn( d2MicroFlowdDrivingStressdFn ); + + set_previousd2MicroGradientFlowdDrivingStressdFn( d2MicroGradientFlowdDrivingStressdFn ); + } else{ @@ -1985,8 +2139,22 @@ namespace tardigradeHydra{ set_dMicroGradientFlowdDrivingStress( dMicroGradientFlowdDrivingStress ); + set_d2MacroFlowdDrivingStressdStress( d2MacroFlowdDrivingStressdMacroStress ); + + set_d2MicroFlowdDrivingStressdStress( d2MicroFlowdDrivingStressdMicroStress ); + + set_d2MicroGradientFlowdDrivingStressdStress( d2MicroGradientFlowdDrivingStressdMicroGradientStress ); + + set_d2MacroFlowdDrivingStressdFn( d2MacroFlowdDrivingStressdFn ); + + set_d2MicroFlowdDrivingStressdFn( d2MicroFlowdDrivingStressdFn ); + + set_d2MicroGradientFlowdDrivingStressdFn( d2MicroGradientFlowdDrivingStressdFn ); + } + std::cout << "exiting setFlowPotentialGradientsJacobians\n"; + } void residual::setMacroCohesion( ){ diff --git a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index 68b27810..ae171192 100644 --- a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -2396,6 +2396,24 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives ){ 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45 }; + floatMatrix dMacroDrivingStressdMacroStress = initialize( 9, 9 ); + + floatMatrix dMacroDrivingStressdF = initialize( 9, 9 ); + + floatMatrix dMacroDrivingStressdFn = initialize( 9, 9 ); + + floatMatrix dMicroDrivingStressdMicroStress = initialize( 9, 9 ); + + floatMatrix dMicroDrivingStressdF = initialize( 9, 9 ); + + floatMatrix dMicroDrivingStressdFn = initialize( 9, 9 ); + + floatMatrix dMicroGradientDrivingStressdHigherOrderStress = initialize( 27, 27 ); + + floatMatrix dMicroGradientDrivingStressdF = initialize( 27, 9 ); + + floatMatrix dMicroGradientDrivingStressdFn = initialize( 27, 9 ); + floatType previousMacroCohesion = -1.23; floatType previousMicroCohesion = -2.34; @@ -2420,6 +2438,14 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives ){ 2, 0.32982199, 0.60161431, 2, 0.58881096, 0.11473813 }; + floatMatrix initialize( unsigned int nrows, unsigned int ncols ){ + + floatMatrix value( nrows, floatVector( ncols, 0 ) ); + + return value; + + } + virtual void setCohesions( const bool isPrevious ) override{ if ( isPrevious ){ @@ -2470,6 +2496,24 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives ){ setDrivingStresses( isPrevious ); + set_dMacroDrivingStressdMacroStress( dMacroDrivingStressdMacroStress ); + + set_dMacroDrivingStressdF( dMacroDrivingStressdF ); + + set_dMacroDrivingStressdFn( dMacroDrivingStressdFn ); + + set_dSymmetricMicroDrivingStressdMicroStress( dMicroDrivingStressdMicroStress ); + + set_dSymmetricMicroDrivingStressdF( dMicroDrivingStressdF ); + + set_dSymmetricMicroDrivingStressdFn( dMicroDrivingStressdFn ); + + set_dHigherOrderDrivingStressdHigherOrderStress( dMicroGradientDrivingStressdHigherOrderStress ); + + set_dHigherOrderDrivingStressdF( dMicroGradientDrivingStressdF ); + + set_dHigherOrderDrivingStressdFn( dMicroGradientDrivingStressdFn ); + } }; @@ -2630,3 +2674,342 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives ){ BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousdMicroGradientFlowdDrivingStress, *RJ.get_previousdMicroGradientFlowdDrivingStress( ) ) ); } + +BOOST_AUTO_TEST_CASE( test_setFlowDerivatives2 ){ + /*! + * Test setting the cohesion values + */ + + // Form a hydra class and residual class + floatType time = 1.23; + + floatType deltaTime = 2.34; + + floatType temperature = 3.45; + + floatType previousTemperature = 4.56; + + variableVector deformationGradient = { 1.00157757, 0.00159138, 0.00672005, + 0.01747159, 1.01122277, 0.00555118, + 0.01112217, -0.00885205, 0.99308943 }; + + variableVector microDeformation = { 0.99460588, -0.0078411 , 0.01145249, + -0.00307139, 0.97798389, -0.00509779, + 0.01189977, -0.01587541, 0.98377259 }; + + variableVector gradientMicroDeformation = { -1.35868385e-02, -1.03142977e-02, 6.54880619e-03, -2.03947530e-02, + -3.31494137e-03, -3.45686183e-03, -3.15745117e-04, -3.70848549e-03, + -9.38693885e-03, -3.68243465e-03, 1.96694582e-02, 2.22080009e-02, + 9.18337942e-05, 6.19764759e-03, -1.92190802e-02, -9.13572591e-03, + -4.25868940e-03, 1.83154579e-02, -1.24772317e-02, -8.48286787e-04, + 2.42779893e-02, 9.74255963e-04, 5.64472629e-03, -1.89685667e-02, + 1.63170400e-02, 5.15300642e-03, 2.25340032e-03 }; + + floatVector previousDeformationGradient = { 9.94656270e-01, 4.82152400e-02, 3.31984800e-02, 2.81918700e-02, + 1.02086536e+00, -1.77592100e-02, -2.24798000e-03, -1.28410000e-04, + 9.77165250e-01 }; + + floatVector previousMicroDeformation = { 0.96917405, -0.01777599, 0.00870406, -0.02163002, 0.9998683 , + -0.01669352, 0.03355217, 0.04427456, 1.01778466 }; + + floatVector previousGradientMicroDeformation = { 0.05043761, 0.02160516, -0.0565408 , 0.01218304, -0.05851034, + 0.00485749, -0.00962607, -0.03455912, 0.04490067, 0.01552915, + -0.02878364, 0.00595866, 0.04750406, -0.02377005, -0.05041534, + -0.02922214, 0.06280788, 0.02850865, -0.00226005, 0.0146049 , + 0.01560184, 0.03224767, 0.05822091, -0.05294424, -0.03518206, + 0.01831308, 0.03774438 }; + + floatVector previousStateVariables = {-0.02495446, -0.00169657, 0.04855598, 0.00194851, 0.01128945, + -0.03793713, 0.03263408, 0.01030601, 0.0045068 , -0.01572362, + -0.01958792, -0.00829778, 0.01813008, 0.03754568, 0.00104223, + 0.01693138, 0.00859366, 0.01249035, 0.01746891, 0.03423424, + -0.0416805 , 0.02636828, -0.02563336, -0.0305777 , 0.0072457 , + -0.04042875, 0.03853268, 0.0127249 , 0.02234164, -0.04838708, + 0.00944319, 0.00567852, -0.03410404, -0.03469295, 0.01955295, + -0.01812336, 0.01919703, 0.00543832, -0.01110494, 0.04251325, + 0.034167 , -0.01426024, -0.04564085, -0.01952319, -0.01018143, + 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0 }; + + floatVector parameters = { 2, 0.53895133, 0.37172145, + 2, 0.37773052, 0.92739145, + 2, 0.53186824, 0.75454313, + 2, 0.95338442, 0.74042148, + 2, 0.38093104, 0.49241325, + 2, 0.82121039, 0.90566759, + 2, 0.01166325, 0.05331896, + 2, 0.32982199, 0.60161431, + 2, 0.58881096, 0.11473813 }; + + floatVector unknownVector = { 4.31202550e-01, -1.78960429e-01, -6.17986089e-01, 9.34988614e-01, + 3.01500733e-01, 7.30919703e-01, -9.49515284e-01, -4.66188370e-01, + 4.14220065e-03, -8.65102730e-01, 9.86066522e-01, -5.27075208e-01, + -2.51415635e-01, -5.71976170e-01, -7.89108268e-01, -5.35040429e-01, + -3.98779729e-01, 2.68884536e-01, -4.37530437e-01, -2.75446478e-01, + -9.88114313e-01, -2.68561748e-01, 6.77719634e-02, -6.75968326e-01, + 1.94866217e-01, -4.13695063e-01, 2.64100990e-01, -9.47606789e-01, + 7.75186921e-01, -9.67762739e-01, -7.46083938e-01, 5.54324923e-01, + -9.08209536e-01, 4.21997387e-01, 9.42092281e-01, 7.43365866e-01, + 4.20323303e-01, 9.17019486e-01, -1.40373324e-01, 7.45757829e-01, + -2.88084664e-01, 8.59527306e-01, -7.02444688e-01, 8.80058030e-01, + 6.65432395e-01, 1.01730274e+00, -1.88038495e-02, 4.82434492e-03, + -2.41803760e-02, 1.01105922e+00, -2.46131243e-02, -2.07588861e-02, + -1.37250795e-02, 1.01875623e+00, 9.93178816e-01, 1.99799676e-03, + 3.40516069e-03, -1.37268320e-02, 1.00360734e+00, 8.04758975e-03, + -1.00877303e-02, -4.06865705e-03, 9.97654446e-01, 2.16175331e-02, + 4.37468737e-03, 2.24126186e-02, 2.80173769e-03, 2.80710425e-05, + -2.48233895e-02, -9.55547806e-04, 2.13727499e-02, -1.50817155e-02, + -2.23954433e-02, -4.66105533e-03, -6.38017597e-03, 1.78576529e-02, + -2.36694442e-02, 2.10074615e-02, 9.04514995e-03, 2.02112997e-02, + 5.37645354e-03, 1.55976656e-02, -8.22280632e-03, -7.52168860e-03, + -5.50628848e-03, 1.27398541e-02, -6.53544128e-03, -1.28890097e-02, + 2.18834178e-02, 2.04005542e-02, + 0.01, 0.02, 0.03, 0.04, 0.05, + 0.06, 0.07, 0.08, 0.09, 0.10 }; + + unsigned int numConfigurations = 2; + + unsigned int numNonLinearSolveStateVariables = 10; + + std::vector< unsigned int > stateVariableIndices = { 0, 1, 2, 3, 4, 5, 6, 7, 8, 9 }; + + unsigned int dimension = 3; + + unsigned int configuration_unknown_count = 45; + + floatType tolr = 1e-2; + + floatType tola = 1e-3; + + unsigned int maxIterations = 24; + + unsigned int maxLSIterations = 45; + + floatType lsAlpha = 2.3; + + class stressMock : public tardigradeHydra::residualBaseMicromorphic{ + + public: + + using tardigradeHydra::residualBaseMicromorphic::residualBaseMicromorphic; + + }; + + class residualMock : public tardigradeHydra::micromorphicDruckerPragerPlasticity::residual{ + + public: + + using tardigradeHydra::micromorphicDruckerPragerPlasticity::residual::residual; + + floatMatrix initialize( unsigned int nrows, unsigned int ncols ){ + + floatMatrix value( nrows, floatVector( ncols, 0 ) ); + + return value; + + } + + }; + + class hydraBaseMicromorphicMock : public tardigradeHydra::hydraBaseMicromorphic{ + + public: + + using tardigradeHydra::hydraBaseMicromorphic::hydraBaseMicromorphic; + + stressMock elasticity; + + residualMock plasticity; + + std::vector< unsigned int > stateVariableIndices = { 0, 1, 2, 3, 4, 5, 6, 7, 8, 9 }; + + floatVector plasticParameters = { 2, 0.53895133, 0.37172145, + 2, 0.37773052, 0.92739145, + 2, 0.53186824, 0.75454313, + 2, 0.95338442, 0.74042148, + 2, 0.38093104, 0.49241325, + 2, 0.82121039, 0.90566759, + 2, 0.01166325, 0.05331896, + 2, 0.32982199, 0.60161431, + 2, 0.58881096, 0.11473813 }; + + private: + + using tardigradeHydra::hydraBaseMicromorphic::setResidualClasses; + + virtual void setResidualClasses( ) override{ + + std::vector< tardigradeHydra::residualBase* > residuals( 2 ); + + elasticity = stressMock( this, 45 ); + + plasticity = residualMock( this, 55, 1, stateVariableIndices, plasticParameters ); + + residuals[ 0 ] = &elasticity; + + residuals[ 1 ] = &plasticity; + + setResidualClasses( residuals ); + + } + + }; + + hydraBaseMicromorphicMock hydra( time, deltaTime, temperature, previousTemperature, deformationGradient, previousDeformationGradient, + microDeformation, previousMicroDeformation, gradientMicroDeformation, previousGradientMicroDeformation, + previousStateVariables, parameters, + numConfigurations, numNonLinearSolveStateVariables, + dimension, configuration_unknown_count, + tolr, tola, maxIterations, maxLSIterations, lsAlpha ); + + tardigradeHydra::unit_test::hydraBaseTester::updateUnknownVector( hydra, unknownVector ); + + residualMock R( &hydra, 55, 1, stateVariableIndices, parameters ); + + // Jacobian tests + + floatMatrix d2MacroFlowdDrivingStressdX( 9, floatVector( unknownVector.size( ), 0 ) ); + + floatMatrix d2MicroFlowdDrivingStressdX( 9, floatVector( unknownVector.size( ), 0 ) ); + + floatMatrix d2MicroGradientFlowdDrivingStressdX( 3 * 27, floatVector( unknownVector.size( ), 0 ) ); + + floatType eps = 1e-6; + + for ( unsigned int i = 0; i < unknownVector.size( ); i++ ){ + + floatVector delta( unknownVector.size( ), 0 ); + + delta[ i ] = eps * std::fabs( unknownVector[ i ] ) + eps; + + hydraBaseMicromorphicMock hydrap( time, deltaTime, temperature, previousTemperature, deformationGradient, previousDeformationGradient, + microDeformation, previousMicroDeformation, gradientMicroDeformation, previousGradientMicroDeformation, + previousStateVariables, parameters, + numConfigurations, numNonLinearSolveStateVariables, + dimension, configuration_unknown_count, + tolr, tola, maxIterations, maxLSIterations, lsAlpha ); + + hydraBaseMicromorphicMock hydram( time, deltaTime, temperature, previousTemperature, deformationGradient, previousDeformationGradient, + microDeformation, previousMicroDeformation, gradientMicroDeformation, previousGradientMicroDeformation, + previousStateVariables, parameters, + numConfigurations, numNonLinearSolveStateVariables, + dimension, configuration_unknown_count, + tolr, tola, maxIterations, maxLSIterations, lsAlpha ); + + tardigradeHydra::unit_test::hydraBaseTester::updateUnknownVector( hydrap, unknownVector + delta ); + + tardigradeHydra::unit_test::hydraBaseTester::updateUnknownVector( hydram, unknownVector - delta ); + + residualMock Rp( &hydrap, 55, 1, stateVariableIndices, parameters ); + + residualMock Rm( &hydram, 55, 1, stateVariableIndices, parameters ); + + floatVector vp = *Rp.get_dMacroFlowdDrivingStress( ); + + floatVector vm = *Rm.get_dMacroFlowdDrivingStress( ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + d2MacroFlowdDrivingStressdX[ j ][ i ] = ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ); + + } + + vp = *Rp.get_dMicroFlowdDrivingStress( ); + + vm = *Rm.get_dMicroFlowdDrivingStress( ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + d2MicroFlowdDrivingStressdX[ j ][ i ] = ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ); + + } + + vp = tardigradeVectorTools::appendVectors( *Rp.get_dMicroGradientFlowdDrivingStress( ) ); + + vm = tardigradeVectorTools::appendVectors( *Rm.get_dMicroGradientFlowdDrivingStress( ) ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + d2MicroGradientFlowdDrivingStressdX[ j ][ i ] = ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ); + + } + + } + + floatMatrix assembled_d2MacroFlowdDrivingStressdX( 9, floatVector( unknownVector.size( ), 0 ) ); + + floatMatrix assembled_d2MicroFlowdDrivingStressdX( 9, floatVector( unknownVector.size( ), 0 ) ); + + floatMatrix assembled_d2MicroGradientFlowdDrivingStressdX( 3 * 27, floatVector( unknownVector.size( ), 0 ) ); + + for ( unsigned int i = 0; i < 9; i++ ){ + + for ( unsigned int j = 0; j < 9; j++ ){ + + assembled_d2MacroFlowdDrivingStressdX[ i ][ j ] = ( *R.get_d2MacroFlowdDrivingStressdStress( ) )[ i ][ j ]; + + } + + for ( unsigned int j = 0; j < 9; j++ ){ + + assembled_d2MacroFlowdDrivingStressdX[ i ][ j + configuration_unknown_count ] = ( *R.get_d2MacroFlowdDrivingStressdFn( ) )[ i ][ j ]; + + } + + } + + for ( unsigned int i = 0; i < 9; i++ ){ + + for ( unsigned int j = 0; j < 9; j++ ){ + + assembled_d2MicroFlowdDrivingStressdX[ i ][ j + 9 ] = ( *R.get_d2MicroFlowdDrivingStressdStress( ) )[ i ][ j ]; + + } + + for ( unsigned int j = 0; j < 9; j++ ){ + + assembled_d2MicroFlowdDrivingStressdX[ i ][ j + configuration_unknown_count ] = ( *R.get_d2MicroFlowdDrivingStressdFn( ) )[ i ][ j ]; + + } + + } + + std::cout << "*R.get:\n"; tardigradeVectorTools::print( *R.get_d2MicroGradientFlowdDrivingStressdStress( ) ); + std::cout << "out\n"; + std::cout << "*R.get2:\n"; tardigradeVectorTools::print( *R.get_d2MicroGradientFlowdDrivingStressdFn( ) ); + std::cout << "out\n"; + for ( unsigned int i = 0; i < 3; i++ ){ + + for ( unsigned int j = 0; j < 27; j++ ){ + + for ( unsigned int k = 0; k < 27; k++ ){ + + assembled_d2MicroGradientFlowdDrivingStressdX[ 27 * i + j ][ k + 18 ] = ( *R.get_d2MicroGradientFlowdDrivingStressdStress( ) )[ i ][ 27 * j + k ]; + + } + + for ( unsigned int k = 0; k < 9; k++ ){ + + assembled_d2MicroGradientFlowdDrivingStressdX[ 27 * i + j ][ k + configuration_unknown_count ] = ( *R.get_d2MicroGradientFlowdDrivingStressdFn( ) )[ i ][ 9 * j + k ]; + + } + + } + + } + + std::cout << "assembled_d2MacroFlowdDrivingStressdX:\n"; tardigradeVectorTools::print( assembled_d2MacroFlowdDrivingStressdX ); + std::cout << "d2MacroFlowdDrivingStressdX:\n"; tardigradeVectorTools::print( d2MacroFlowdDrivingStressdX ); + + std::cout << "assembled_d2MicroFlowdDrivingStressdX:\n"; tardigradeVectorTools::print( assembled_d2MicroFlowdDrivingStressdX ); + std::cout << "d2MicroFlowdDrivingStressdX:\n"; tardigradeVectorTools::print( d2MicroFlowdDrivingStressdX ); + + std::cout << "assembled_d2MicroGradientFlowdDrivingStressdX:\n"; tardigradeVectorTools::print( assembled_d2MicroGradientFlowdDrivingStressdX ); + std::cout << "d2MicroGradientFlowdDrivingStressdX:\n"; tardigradeVectorTools::print( d2MicroGradientFlowdDrivingStressdX ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( assembled_d2MacroFlowdDrivingStressdX, d2MacroFlowdDrivingStressdX ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( assembled_d2MicroFlowdDrivingStressdX, d2MicroFlowdDrivingStressdX ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( assembled_d2MicroGradientFlowdDrivingStressdX, d2MicroGradientFlowdDrivingStressdX ) ); + +} From ec99b12fd84c1f83173afea32fdc169b337cd751 Mon Sep 17 00:00:00 2001 From: Nathan Miller Date: Mon, 8 Jan 2024 21:17:25 -0700 Subject: [PATCH 13/15] FEAT: Added Jacobians of the flow potential derivative w.r.t. the unknown vector --- ...draMicromorphicDruckerPragerPlasticity.cpp | 61 +++++++++++++++++-- ...hydraMicromorphicDruckerPragerPlasticity.h | 16 +++++ ...draMicromorphicDruckerPragerPlasticity.cpp | 25 ++++---- 3 files changed, 84 insertions(+), 18 deletions(-) diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index 799a5944..2ee1ed19 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -1736,6 +1736,24 @@ namespace tardigradeHydra{ } + void residual::setd2MicroGradientFlowdDrivingStressdChi( ){ + /*! + * Set the Jacobian of the derivative of the micro gradient flow potential w.r.t. the micro-gradient driving stress w.r.t. the micro deformation + */ + + setFlowPotentialGradientsJacobians( false ); + + } + + void residual::setd2MicroGradientFlowdDrivingStressdChin( ){ + /*! + * Set the Jacobian of the derivative of the micro gradient flow potential w.r.t. the micro-gradient driving stress w.r.t. the sub micro deformation + */ + + setFlowPotentialGradientsJacobians( false ); + + } + void residual::setPreviousd2MacroFlowdDrivingStressdStress( ){ /*! * Set the Jacobian of the previous derivative of the macro flow potential w.r.t. the macro driving stress w.r.t. the macro stress @@ -1817,6 +1835,24 @@ namespace tardigradeHydra{ } + void residual::setPreviousd2MicroGradientFlowdDrivingStressdChi( ){ + /*! + * Set the Jacobian of the previous derivative of the micro gradient flow potential w.r.t. the micro-gradient driving stress w.r.t. the micro deformation + */ + + setFlowPotentialGradientsJacobians( true ); + + } + + void residual::setPreviousd2MicroGradientFlowdDrivingStressdChin( ){ + /*! + * Set the Jacobian of the previous derivative of the micro gradient flow potential w.r.t. the micro-gradient driving stress w.r.t. the sub micro deformation + */ + + setFlowPotentialGradientsJacobians( true ); + + } + void residual::setFlowPotentialGradientsJacobians( const bool isPrevious ){ /*! * Set the Jacobians of the flow potential gradients @@ -1824,8 +1860,6 @@ namespace tardigradeHydra{ * \param isPrevious: A flag for whether to set the current (false) or previous (true) derivatives */ - std::cout << "entering setFlowPotentialGradientsJacobians\n"; - const floatType *macroCohesion; const floatType *microCohesion; @@ -1856,6 +1890,10 @@ namespace tardigradeHydra{ const floatMatrix *dMicroGradientDrivingStressdFn; + const floatMatrix *dMicroGradientDrivingStressdChi; + + const floatMatrix *dMicroGradientDrivingStressdChin; + floatVector precedingDeformationGradient; floatMatrix dPrecedingFdSubFs; @@ -1898,12 +1936,16 @@ namespace tardigradeHydra{ dMicroGradientDrivingStressdF = get_previousdHigherOrderDrivingStressdF( ); + dMicroGradientDrivingStressdChi = get_previousdHigherOrderDrivingStressdChi( ); + dMacroDrivingStressdFn = get_previousdMacroDrivingStressdFn( ); dMicroDrivingStressdFn = get_previousdSymmetricMicroDrivingStressdFn( ); dMicroGradientDrivingStressdFn = get_previousdHigherOrderDrivingStressdFn( ); + dMicroGradientDrivingStressdChin = get_previousdHigherOrderDrivingStressdChin( ); + macroDrivingStress = get_previousMacroDrivingStress( ); microDrivingStress = get_previousSymmetricMicroDrivingStress( ); @@ -1939,12 +1981,16 @@ namespace tardigradeHydra{ dMicroGradientDrivingStressdF = get_dHigherOrderDrivingStressdF( ); + dMicroGradientDrivingStressdChi = get_dHigherOrderDrivingStressdChi( ); + dMacroDrivingStressdFn = get_dMacroDrivingStressdFn( ); dMicroDrivingStressdFn = get_dSymmetricMicroDrivingStressdFn( ); dMicroGradientDrivingStressdFn = get_dHigherOrderDrivingStressdFn( ); + dMicroGradientDrivingStressdChin = get_dHigherOrderDrivingStressdChin( ); + macroDrivingStress = get_macroDrivingStress( ); microDrivingStress = get_symmetricMicroDrivingStress( ); @@ -2027,6 +2073,8 @@ namespace tardigradeHydra{ floatMatrix d2MicroGradientFlowdDrivingStressdFn( dMicroGradientFlowdDrivingStress.size( ), floatVector( microGradientDrivingStress->size( ) * ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ), 0 ) ); + floatMatrix d2MicroGradientFlowdDrivingStressdChin( dMicroGradientFlowdDrivingStress.size( ), floatVector( microGradientDrivingStress->size( ) * ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ), 0 ) ); + for ( unsigned int I = 0; I < d2MicroFlowdDrivingStress2.size( ); I++ ){ for ( unsigned int J = 0; J < macroDrivingStress->size( ); J++ ){ @@ -2079,6 +2127,9 @@ namespace tardigradeHydra{ d2MicroGradientFlowdDrivingStressdFn[ I ][ ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ) * J + L ] += d2MicroGradientFlowdDrivingStress2[ I ][ microGradientDrivingStress->size( ) * J + K ] * ( *dMicroGradientDrivingStressdFn ) [ K ][ L ]; + d2MicroGradientFlowdDrivingStressdChin[ I ][ ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ) * J + L ] + += d2MicroGradientFlowdDrivingStress2[ I ][ microGradientDrivingStress->size( ) * J + K ] * ( *dMicroGradientDrivingStressdChin ) [ K ][ L ]; + } } @@ -2124,6 +2175,8 @@ namespace tardigradeHydra{ set_previousd2MicroGradientFlowdDrivingStressdFn( d2MicroGradientFlowdDrivingStressdFn ); + set_previousd2MicroGradientFlowdDrivingStressdChin( d2MicroGradientFlowdDrivingStressdChin ); + } else{ @@ -2151,9 +2204,9 @@ namespace tardigradeHydra{ set_d2MicroGradientFlowdDrivingStressdFn( d2MicroGradientFlowdDrivingStressdFn ); - } + set_d2MicroGradientFlowdDrivingStressdChin( d2MicroGradientFlowdDrivingStressdChin ); - std::cout << "exiting setFlowPotentialGradientsJacobians\n"; + } } diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h index cc35289f..8c00efcc 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.h @@ -318,6 +318,10 @@ namespace tardigradeHydra{ virtual void setd2MicroGradientFlowdDrivingStressdFn( ); + virtual void setd2MicroGradientFlowdDrivingStressdChi( ); + + virtual void setd2MicroGradientFlowdDrivingStressdChin( ); + virtual void setPreviousd2MacroFlowdDrivingStressdStress( ); virtual void setPreviousd2MacroFlowdDrivingStressdF( ); @@ -336,6 +340,10 @@ namespace tardigradeHydra{ virtual void setPreviousd2MicroGradientFlowdDrivingStressdFn( ); + virtual void setPreviousd2MicroGradientFlowdDrivingStressdChi( ); + + virtual void setPreviousd2MicroGradientFlowdDrivingStressdChin( ); + virtual void setFlowPotentialGradientsJacobians( const bool isPrevious ); private: @@ -486,6 +494,10 @@ namespace tardigradeHydra{ TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, d2MicroGradientFlowdDrivingStressdFn, floatMatrix, setd2MicroGradientFlowdDrivingStressdFn ) + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, d2MicroGradientFlowdDrivingStressdChi, floatMatrix, setd2MicroGradientFlowdDrivingStressdChi ) + + TARDIGRADE_HYDRA_DECLARE_ITERATION_STORAGE( private, d2MicroGradientFlowdDrivingStressdChin, floatMatrix, setd2MicroGradientFlowdDrivingStressdChin ) + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousdMacroFlowdc, floatType, setPreviousdMacroFlowdc ) TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousdMacroFlowdDrivingStress, floatVector, setPreviousdMacroFlowdDrivingStress ) @@ -516,6 +528,10 @@ namespace tardigradeHydra{ TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousd2MicroGradientFlowdDrivingStressdFn, floatMatrix, setPreviousd2MicroGradientFlowdDrivingStressdFn ) + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousd2MicroGradientFlowdDrivingStressdChi, floatMatrix, setPreviousd2MicroGradientFlowdDrivingStressdChi ) + + TARDIGRADE_HYDRA_DECLARE_PREVIOUS_STORAGE( private, previousd2MicroGradientFlowdDrivingStressdChin, floatMatrix, setPreviousd2MicroGradientFlowdDrivingStressdChin ) + }; } diff --git a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index ae171192..a89cb764 100644 --- a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -2414,6 +2414,10 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives ){ floatMatrix dMicroGradientDrivingStressdFn = initialize( 27, 9 ); + floatMatrix dMicroGradientDrivingStressdChi = initialize( 27, 9 ); + + floatMatrix dMicroGradientDrivingStressdChin = initialize( 27, 9 ); + floatType previousMacroCohesion = -1.23; floatType previousMicroCohesion = -2.34; @@ -2514,6 +2518,10 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives ){ set_dHigherOrderDrivingStressdFn( dMicroGradientDrivingStressdFn ); + set_dHigherOrderDrivingStressdChi( dMicroGradientDrivingStressdChi ); + + set_dHigherOrderDrivingStressdChin( dMicroGradientDrivingStressdChin ); + } }; @@ -2973,10 +2981,6 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives2 ){ } - std::cout << "*R.get:\n"; tardigradeVectorTools::print( *R.get_d2MicroGradientFlowdDrivingStressdStress( ) ); - std::cout << "out\n"; - std::cout << "*R.get2:\n"; tardigradeVectorTools::print( *R.get_d2MicroGradientFlowdDrivingStressdFn( ) ); - std::cout << "out\n"; for ( unsigned int i = 0; i < 3; i++ ){ for ( unsigned int j = 0; j < 27; j++ ){ @@ -2989,7 +2993,9 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives2 ){ for ( unsigned int k = 0; k < 9; k++ ){ - assembled_d2MicroGradientFlowdDrivingStressdX[ 27 * i + j ][ k + configuration_unknown_count ] = ( *R.get_d2MicroGradientFlowdDrivingStressdFn( ) )[ i ][ 9 * j + k ]; + assembled_d2MicroGradientFlowdDrivingStressdX[ 27 * i + j ][ k + configuration_unknown_count ] = ( *R.get_d2MicroGradientFlowdDrivingStressdFn( ) )[ i ][ 9 * j + k ]; + + assembled_d2MicroGradientFlowdDrivingStressdX[ 27 * i + j ][ k + configuration_unknown_count + 9 ] = ( *R.get_d2MicroGradientFlowdDrivingStressdChin( ) )[ i ][ 9 * j + k ]; } @@ -2997,15 +3003,6 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives2 ){ } - std::cout << "assembled_d2MacroFlowdDrivingStressdX:\n"; tardigradeVectorTools::print( assembled_d2MacroFlowdDrivingStressdX ); - std::cout << "d2MacroFlowdDrivingStressdX:\n"; tardigradeVectorTools::print( d2MacroFlowdDrivingStressdX ); - - std::cout << "assembled_d2MicroFlowdDrivingStressdX:\n"; tardigradeVectorTools::print( assembled_d2MicroFlowdDrivingStressdX ); - std::cout << "d2MicroFlowdDrivingStressdX:\n"; tardigradeVectorTools::print( d2MicroFlowdDrivingStressdX ); - - std::cout << "assembled_d2MicroGradientFlowdDrivingStressdX:\n"; tardigradeVectorTools::print( assembled_d2MicroGradientFlowdDrivingStressdX ); - std::cout << "d2MicroGradientFlowdDrivingStressdX:\n"; tardigradeVectorTools::print( d2MicroGradientFlowdDrivingStressdX ); - BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( assembled_d2MacroFlowdDrivingStressdX, d2MacroFlowdDrivingStressdX ) ); BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( assembled_d2MicroFlowdDrivingStressdX, d2MicroFlowdDrivingStressdX ) ); From d12129d0ecf2621a9675b2f1888444bf0e5944a2 Mon Sep 17 00:00:00 2001 From: Nathan Miller Date: Mon, 8 Jan 2024 23:53:25 -0700 Subject: [PATCH 14/15] FEAT: Completed the calculation of the gradients of the flow potential --- ...draMicromorphicDruckerPragerPlasticity.cpp | 77 +- ...draMicromorphicDruckerPragerPlasticity.cpp | 688 +++++++++++++++++- 2 files changed, 738 insertions(+), 27 deletions(-) diff --git a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index 2ee1ed19..2e3eff94 100644 --- a/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -2061,12 +2061,20 @@ namespace tardigradeHydra{ tempVectorYield, dMicroGradientFlowdDrivingStress, dMicroGradientFlowdCohesion, dMicroGradientFlowdPrecedingF, d2MicroGradientFlowdDrivingStress2, d2MicroGradientFlowdDrivingStressdPrecedingF ) ); - floatMatrix d2MacroFlowdDrivingStressdMacroStress( dMacroFlowdDrivingStress.size( ), floatVector( macroDrivingStress->size( ) * macroDrivingStress->size( ), 0 ) ); + floatMatrix d2MacroFlowdDrivingStressdMacroStress( dMacroFlowdDrivingStress.size( ), floatVector( macroDrivingStress->size( ), 0 ) ); - floatMatrix d2MicroFlowdDrivingStressdMicroStress( dMicroFlowdDrivingStress.size( ), floatVector( microDrivingStress->size( ) * microDrivingStress->size( ), 0 ) ); + floatMatrix d2MicroFlowdDrivingStressdMicroStress( dMicroFlowdDrivingStress.size( ), floatVector( microDrivingStress->size( ), 0 ) ); floatMatrix d2MicroGradientFlowdDrivingStressdMicroGradientStress( dMicroGradientFlowdDrivingStress.size( ), floatVector( microGradientDrivingStress->size( ) * microGradientDrivingStress->size( ), 0 ) ); + floatMatrix d2MacroFlowdDrivingStressdF( dMacroFlowdDrivingStress.size( ), floatVector( precedingDeformationGradient.size( ), 0 ) ); + + floatMatrix d2MicroFlowdDrivingStressdF( dMicroFlowdDrivingStress.size( ), floatVector( precedingDeformationGradient.size( ), 0 ) ); + + floatMatrix d2MicroGradientFlowdDrivingStressdF( dMicroGradientFlowdDrivingStress.size( ), floatVector( microGradientDrivingStress->size( ) * precedingDeformationGradient.size( ), 0 ) ); + + floatMatrix d2MicroGradientFlowdDrivingStressdChi( dMicroGradientFlowdDrivingStress.size( ), floatVector( microGradientDrivingStress->size( ) * precedingDeformationGradient.size( ), 0 ) ); + floatMatrix d2MacroFlowdDrivingStressdFn( dMacroFlowdDrivingStress.size( ), floatVector( macroDrivingStress->size( ) * ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ), 0 ) ); floatMatrix d2MicroFlowdDrivingStressdFn( dMicroFlowdDrivingStress.size( ), floatVector( microDrivingStress->size( ) * ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ), 0 ) ); @@ -2081,27 +2089,31 @@ namespace tardigradeHydra{ for ( unsigned int K = 0; K < macroDrivingStress->size( ); K++ ){ - for ( unsigned int L = 0; L < macroDrivingStress->size( ); L++ ){ + d2MacroFlowdDrivingStressdMacroStress[ I ][ J ] + += d2MacroFlowdDrivingStress2[ I ][ K ] * ( *dMacroDrivingStressdStress )[ K ][ J ]; - d2MacroFlowdDrivingStressdMacroStress[ I ][ macroDrivingStress->size( ) * J + K ] - += d2MacroFlowdDrivingStress2[ I ][ macroDrivingStress->size( ) * J + L ] * ( *dMacroDrivingStressdStress )[ L ][ K ]; + d2MicroFlowdDrivingStressdMicroStress[ I ][ J ] + += d2MicroFlowdDrivingStress2[ I ][ K ] * ( *dMicroDrivingStressdStress )[ K ][ J ]; - d2MicroFlowdDrivingStressdMicroStress[ I ][ microDrivingStress->size( ) * J + K ] - += d2MicroFlowdDrivingStress2[ I ][ microDrivingStress->size( ) * J + L ] * ( *dMicroDrivingStressdStress )[ L ][ K ]; + d2MacroFlowdDrivingStressdF[ I ][ J ] + += d2MacroFlowdDrivingStress2[ I ][ K ] * ( *dMacroDrivingStressdF )[ K ][ J ] + + d2MacroFlowdDrivingStressdPrecedingF[ I ][ K ] * dPrecedingFdF[ K ][ J ]; - } + d2MicroFlowdDrivingStressdF[ I ][ J ] + += d2MicroFlowdDrivingStress2[ I ][ K ] * ( *dMicroDrivingStressdF )[ K ][ J ] + + d2MicroFlowdDrivingStressdPrecedingF[ I ][ K ] * dPrecedingFdF[ K ][ J ]; - for ( unsigned int L = 0; L < ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ); L++ ){ + } - d2MacroFlowdDrivingStressdFn[ I ][ ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ) * J + L ] - += d2MacroFlowdDrivingStressdPrecedingF[ I ][ macroDrivingStress->size( ) * J + K ] * dPrecedingFdFn[ K ][ L ] - + d2MacroFlowdDrivingStress2[ I ][ macroDrivingStress->size( ) * J + K ] * ( *dMacroDrivingStressdFn ) [ K ][ L ]; + for ( unsigned int K = 0; K < ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ); K++ ){ - d2MicroFlowdDrivingStressdFn[ I ][ ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ) * J + L ] - += d2MicroFlowdDrivingStressdPrecedingF[ I ][ precedingDeformationGradient.size( ) * J + K ] * dPrecedingFdFn[ K ][ L ] - + d2MicroFlowdDrivingStress2[ I ][ microDrivingStress->size( ) * J + K ] * ( *dMicroDrivingStressdFn ) [ K ][ L ]; + d2MacroFlowdDrivingStressdFn[ I ][ K ] + += d2MacroFlowdDrivingStressdPrecedingF[ I ][ J ] * dPrecedingFdFn[ J ][ K ] + + d2MacroFlowdDrivingStress2[ I ][ J ] * ( *dMacroDrivingStressdFn )[ J ][ K ]; - } + d2MicroFlowdDrivingStressdFn[ I ][ K ] + += d2MicroFlowdDrivingStressdPrecedingF[ I ][ J ] * dPrecedingFdFn[ J ][ K ] + + d2MicroFlowdDrivingStress2[ I ][ J ] * ( *dMicroDrivingStressdFn ) [ J ][ K ]; } @@ -2122,6 +2134,16 @@ namespace tardigradeHydra{ } + for ( unsigned int L = 0; L < precedingDeformationGradient.size( ); L++ ){ + + d2MicroGradientFlowdDrivingStressdF[ I ][ precedingDeformationGradient.size( ) * J + L ] + += d2MicroGradientFlowdDrivingStress2[ I ][ microGradientDrivingStress->size( ) * J + K ] * ( *dMicroGradientDrivingStressdF )[ K ][ L ]; + + d2MicroGradientFlowdDrivingStressdChi[ I ][ precedingDeformationGradient.size( ) * J + L ] + += d2MicroGradientFlowdDrivingStress2[ I ][ microGradientDrivingStress->size( ) * J + K ] * ( *dMicroGradientDrivingStressdChi )[ K ][ L ]; + + } + for ( unsigned int L = 0; L < ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ); L++ ){ d2MicroGradientFlowdDrivingStressdFn[ I ][ ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ) * J + L ] @@ -2136,6 +2158,13 @@ namespace tardigradeHydra{ for ( unsigned int K = 0; K < precedingDeformationGradient.size( ); K++ ){ + for ( unsigned int L = 0; L < precedingDeformationGradient.size( ); L++ ){ + + d2MicroGradientFlowdDrivingStressdF[ I ][ precedingDeformationGradient.size( ) * J + L ] + += d2MicroGradientFlowdDrivingStressdPrecedingF[ I ][ precedingDeformationGradient.size( ) * J + K ] * dPrecedingFdF[ K ][ L ]; + + } + for ( unsigned int L = 0; L < ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ); L++ ){ d2MicroGradientFlowdDrivingStressdFn[ I ][ ( ( *hydra->getNumConfigurations( ) ) - 1 ) * precedingDeformationGradient.size( ) * J + L ] @@ -2175,6 +2204,14 @@ namespace tardigradeHydra{ set_previousd2MicroGradientFlowdDrivingStressdFn( d2MicroGradientFlowdDrivingStressdFn ); + set_previousd2MacroFlowdDrivingStressdF( d2MacroFlowdDrivingStressdF ); + + set_previousd2MicroFlowdDrivingStressdF( d2MicroFlowdDrivingStressdF ); + + set_previousd2MicroGradientFlowdDrivingStressdF( d2MicroGradientFlowdDrivingStressdF ); + + set_previousd2MicroGradientFlowdDrivingStressdChi( d2MicroGradientFlowdDrivingStressdChi ); + set_previousd2MicroGradientFlowdDrivingStressdChin( d2MicroGradientFlowdDrivingStressdChin ); } @@ -2204,6 +2241,14 @@ namespace tardigradeHydra{ set_d2MicroGradientFlowdDrivingStressdFn( d2MicroGradientFlowdDrivingStressdFn ); + set_d2MacroFlowdDrivingStressdF( d2MacroFlowdDrivingStressdF ); + + set_d2MicroFlowdDrivingStressdF( d2MicroFlowdDrivingStressdF ); + + set_d2MicroGradientFlowdDrivingStressdF( d2MicroGradientFlowdDrivingStressdF ); + + set_d2MicroGradientFlowdDrivingStressdChi( d2MicroGradientFlowdDrivingStressdChi ); + set_d2MicroGradientFlowdDrivingStressdChin( d2MicroGradientFlowdDrivingStressdChin ); } diff --git a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp index a89cb764..3e243742 100644 --- a/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp +++ b/src/cpp/tests/test_tardigrade_hydraMicromorphicDruckerPragerPlasticity.cpp @@ -2500,27 +2500,56 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives ){ setDrivingStresses( isPrevious ); - set_dMacroDrivingStressdMacroStress( dMacroDrivingStressdMacroStress ); + if ( isPrevious ){ + + set_previousdMacroDrivingStressdMacroStress( dMacroDrivingStressdMacroStress ); + + set_previousdMacroDrivingStressdF( dMacroDrivingStressdF ); + + set_previousdMacroDrivingStressdFn( dMacroDrivingStressdFn ); + + set_previousdSymmetricMicroDrivingStressdMicroStress( dMicroDrivingStressdMicroStress ); + + set_previousdSymmetricMicroDrivingStressdF( dMicroDrivingStressdF ); + + set_previousdSymmetricMicroDrivingStressdFn( dMicroDrivingStressdFn ); + + set_previousdHigherOrderDrivingStressdHigherOrderStress( dMicroGradientDrivingStressdHigherOrderStress ); + + set_previousdHigherOrderDrivingStressdF( dMicroGradientDrivingStressdF ); + + set_previousdHigherOrderDrivingStressdFn( dMicroGradientDrivingStressdFn ); + + set_previousdHigherOrderDrivingStressdChi( dMicroGradientDrivingStressdChi ); + + set_previousdHigherOrderDrivingStressdChin( dMicroGradientDrivingStressdChin ); + + } + else{ + + set_dMacroDrivingStressdMacroStress( dMacroDrivingStressdMacroStress ); - set_dMacroDrivingStressdF( dMacroDrivingStressdF ); + set_dMacroDrivingStressdF( dMacroDrivingStressdF ); - set_dMacroDrivingStressdFn( dMacroDrivingStressdFn ); + set_dMacroDrivingStressdFn( dMacroDrivingStressdFn ); - set_dSymmetricMicroDrivingStressdMicroStress( dMicroDrivingStressdMicroStress ); + set_dSymmetricMicroDrivingStressdMicroStress( dMicroDrivingStressdMicroStress ); - set_dSymmetricMicroDrivingStressdF( dMicroDrivingStressdF ); + set_dSymmetricMicroDrivingStressdF( dMicroDrivingStressdF ); - set_dSymmetricMicroDrivingStressdFn( dMicroDrivingStressdFn ); + set_dSymmetricMicroDrivingStressdFn( dMicroDrivingStressdFn ); - set_dHigherOrderDrivingStressdHigherOrderStress( dMicroGradientDrivingStressdHigherOrderStress ); + set_dHigherOrderDrivingStressdHigherOrderStress( dMicroGradientDrivingStressdHigherOrderStress ); - set_dHigherOrderDrivingStressdF( dMicroGradientDrivingStressdF ); + set_dHigherOrderDrivingStressdF( dMicroGradientDrivingStressdF ); - set_dHigherOrderDrivingStressdFn( dMicroGradientDrivingStressdFn ); + set_dHigherOrderDrivingStressdFn( dMicroGradientDrivingStressdFn ); - set_dHigherOrderDrivingStressdChi( dMicroGradientDrivingStressdChi ); + set_dHigherOrderDrivingStressdChi( dMicroGradientDrivingStressdChi ); - set_dHigherOrderDrivingStressdChin( dMicroGradientDrivingStressdChin ); + set_dHigherOrderDrivingStressdChin( dMicroGradientDrivingStressdChin ); + + } } @@ -2800,6 +2829,24 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives2 ){ using tardigradeHydra::residualBaseMicromorphic::residualBaseMicromorphic; + floatVector previousPK2Stress = { 1, 2, 3, 4, 5, 6, 7, 8, 9 }; + + floatVector previousSIGMA = { 10, 11, 12, 13, 14, 15, 16, 17, 18 }; + + floatVector previousM = { 19, 20, 21, 22, 23, 24, 25, 26, 27, + 28, 29, 30, 31, 32, 33, 34, 35, 36, + 37, 38, 39, 40, 41, 42, 43, 44, 45 }; + + protected: + + using tardigradeHydra::residualBase::setPreviousStress; + + virtual void setPreviousStress( ) override{ + + setPreviousStress( tardigradeVectorTools::appendVectors( { previousPK2Stress, previousSIGMA, previousM } ) ); + + } + }; class residualMock : public tardigradeHydra::micromorphicDruckerPragerPlasticity::residual{ @@ -2840,6 +2887,28 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives2 ){ 2, 0.32982199, 0.60161431, 2, 0.58881096, 0.11473813 }; + floatVector _local_deltaPK2 = initialize( 9 ); + + floatVector _local_deltaSIGMA = initialize( 9 ); + + floatVector _local_deltaM = initialize( 27 ); + + floatVector initialize( unsigned int nrows ){ + + floatVector value( nrows, 0 ); + + return value; + + } + + floatMatrix initialize( unsigned int nrows, unsigned int ncols ){ + + floatMatrix value( nrows, floatVector( ncols, 0 ) ); + + return value; + + } + private: using tardigradeHydra::hydraBaseMicromorphic::setResidualClasses; @@ -2850,6 +2919,12 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives2 ){ elasticity = stressMock( this, 45 ); + elasticity.previousPK2Stress += _local_deltaPK2; + + elasticity.previousSIGMA += _local_deltaSIGMA; + + elasticity.previousM += _local_deltaM; + plasticity = residualMock( this, 55, 1, stateVariableIndices, plasticParameters ); residuals[ 0 ] = &elasticity; @@ -2881,8 +2956,39 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives2 ){ floatMatrix d2MicroGradientFlowdDrivingStressdX( 3 * 27, floatVector( unknownVector.size( ), 0 ) ); + floatMatrix d2MacroFlowdDrivingStressdF( 9, floatVector( 9, 0 ) ); + + floatMatrix d2MicroFlowdDrivingStressdF( 9, floatVector( 9, 0 ) ); + + floatMatrix d2MicroGradientFlowdDrivingStressdF( 3 * 27, floatVector( 9, 0 ) ); + + floatMatrix d2MacroFlowdDrivingStressdChi( 9, floatVector( 9, 0 ) ); + + floatMatrix d2MicroFlowdDrivingStressdChi( 9, floatVector( 9, 0 ) ); + + floatMatrix d2MicroGradientFlowdDrivingStressdChi( 3 * 27, floatVector( 9, 0 ) ); + + floatMatrix previousd2MacroFlowdDrivingStressdMacroStress( 9, floatVector( 9, 0 ) ); + + floatMatrix previousd2MicroFlowdDrivingStressdMicroStress( 9, floatVector( 9, 0 ) ); + + floatMatrix previousd2MicroGradientFlowdDrivingStressdHigherOrderStress( 3 * 27, floatVector( 27, 0 ) ); + + floatMatrix previousd2MacroFlowdDrivingStressdX( 9, floatVector( previousStateVariables.size( ), 0 ) ); + + floatMatrix previousd2MicroFlowdDrivingStressdX( 9, floatVector( previousStateVariables.size( ), 0 ) ); + + floatMatrix previousd2MicroGradientFlowdDrivingStressdX( 3 * 27, floatVector( previousStateVariables.size( ), 0 ) ); + + floatMatrix previousd2MacroFlowdDrivingStressdF( 9, floatVector( 9, 0 ) ); + + floatMatrix previousd2MicroFlowdDrivingStressdF( 9, floatVector( 9, 0 ) ); + + floatMatrix previousd2MicroGradientFlowdDrivingStressdF( 3 * 27, floatVector( 9, 0 ) ); + floatType eps = 1e-6; + // derivatives w.r.t. current stresses, Fn, Chin, gradchin, and ISVs for ( unsigned int i = 0; i < unknownVector.size( ); i++ ){ floatVector delta( unknownVector.size( ), 0 ); @@ -3009,4 +3115,564 @@ BOOST_AUTO_TEST_CASE( test_setFlowDerivatives2 ){ BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( assembled_d2MicroGradientFlowdDrivingStressdX, d2MicroGradientFlowdDrivingStressdX ) ); + // derivatives w.r.t. current F + for ( unsigned int i = 0; i < deformationGradient.size( ); i++ ){ + + floatVector delta( deformationGradient.size( ), 0 ); + + delta[ i ] = eps * std::fabs( deformationGradient[ i ] ) + eps; + + hydraBaseMicromorphicMock hydrap( time, deltaTime, temperature, previousTemperature, deformationGradient + delta, previousDeformationGradient, + microDeformation, previousMicroDeformation, gradientMicroDeformation, previousGradientMicroDeformation, + previousStateVariables, parameters, + numConfigurations, numNonLinearSolveStateVariables, + dimension, configuration_unknown_count, + tolr, tola, maxIterations, maxLSIterations, lsAlpha ); + + hydraBaseMicromorphicMock hydram( time, deltaTime, temperature, previousTemperature, deformationGradient - delta, previousDeformationGradient, + microDeformation, previousMicroDeformation, gradientMicroDeformation, previousGradientMicroDeformation, + previousStateVariables, parameters, + numConfigurations, numNonLinearSolveStateVariables, + dimension, configuration_unknown_count, + tolr, tola, maxIterations, maxLSIterations, lsAlpha ); + + tardigradeHydra::unit_test::hydraBaseTester::updateUnknownVector( hydrap, unknownVector ); + + tardigradeHydra::unit_test::hydraBaseTester::updateUnknownVector( hydram, unknownVector ); + + residualMock Rp( &hydrap, 55, 1, stateVariableIndices, parameters ); + + residualMock Rm( &hydram, 55, 1, stateVariableIndices, parameters ); + + floatVector vp = *Rp.get_dMacroFlowdDrivingStress( ); + + floatVector vm = *Rm.get_dMacroFlowdDrivingStress( ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + d2MacroFlowdDrivingStressdF[ j ][ i ] = ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ); + + } + + vp = *Rp.get_dMicroFlowdDrivingStress( ); + + vm = *Rm.get_dMicroFlowdDrivingStress( ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + d2MicroFlowdDrivingStressdF[ j ][ i ] = ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ); + + } + + vp = tardigradeVectorTools::appendVectors( *Rp.get_dMicroGradientFlowdDrivingStress( ) ); + + vm = tardigradeVectorTools::appendVectors( *Rm.get_dMicroGradientFlowdDrivingStress( ) ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + d2MicroGradientFlowdDrivingStressdF[ j ][ i ] = ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ); + + } + + } + + floatMatrix assembled_d2MicroGradientFlowdDrivingStressdF( 3 * 27, floatVector( 9, 0 ) ); + + for ( unsigned int i = 0; i < 3; i++ ){ + + for ( unsigned int j = 0; j < 27; j++ ){ + + for ( unsigned int k = 0; k < 9; k++ ){ + + assembled_d2MicroGradientFlowdDrivingStressdF[ 27 * i + j ][ k ] = ( *R.get_d2MicroGradientFlowdDrivingStressdF( ) )[ i ][ 9 * j + k ]; + + } + + } + + } + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( d2MacroFlowdDrivingStressdF, *R.get_d2MacroFlowdDrivingStressdF( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( d2MicroFlowdDrivingStressdF, *R.get_d2MicroFlowdDrivingStressdF( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( d2MicroGradientFlowdDrivingStressdF, assembled_d2MicroGradientFlowdDrivingStressdF ) ); + + // derivatives w.r.t. current chi + for ( unsigned int i = 0; i < microDeformation.size( ); i++ ){ + + floatVector delta( microDeformation.size( ), 0 ); + + delta[ i ] = eps * std::fabs( microDeformation[ i ] ) + eps; + + hydraBaseMicromorphicMock hydrap( time, deltaTime, temperature, previousTemperature, deformationGradient, previousDeformationGradient, + microDeformation + delta, previousMicroDeformation, gradientMicroDeformation, previousGradientMicroDeformation, + previousStateVariables, parameters, + numConfigurations, numNonLinearSolveStateVariables, + dimension, configuration_unknown_count, + tolr, tola, maxIterations, maxLSIterations, lsAlpha ); + + hydraBaseMicromorphicMock hydram( time, deltaTime, temperature, previousTemperature, deformationGradient, previousDeformationGradient, + microDeformation - delta, previousMicroDeformation, gradientMicroDeformation, previousGradientMicroDeformation, + previousStateVariables, parameters, + numConfigurations, numNonLinearSolveStateVariables, + dimension, configuration_unknown_count, + tolr, tola, maxIterations, maxLSIterations, lsAlpha ); + + tardigradeHydra::unit_test::hydraBaseTester::updateUnknownVector( hydrap, unknownVector ); + + tardigradeHydra::unit_test::hydraBaseTester::updateUnknownVector( hydram, unknownVector ); + + residualMock Rp( &hydrap, 55, 1, stateVariableIndices, parameters ); + + residualMock Rm( &hydram, 55, 1, stateVariableIndices, parameters ); + + floatVector vp = *Rp.get_dMacroFlowdDrivingStress( ); + + floatVector vm = *Rm.get_dMacroFlowdDrivingStress( ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + d2MacroFlowdDrivingStressdChi[ j ][ i ] = ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ); + + } + + vp = *Rp.get_dMicroFlowdDrivingStress( ); + + vm = *Rm.get_dMicroFlowdDrivingStress( ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + d2MicroFlowdDrivingStressdChi[ j ][ i ] = ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ); + + } + + vp = tardigradeVectorTools::appendVectors( *Rp.get_dMicroGradientFlowdDrivingStress( ) ); + + vm = tardigradeVectorTools::appendVectors( *Rm.get_dMicroGradientFlowdDrivingStress( ) ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + d2MicroGradientFlowdDrivingStressdChi[ j ][ i ] = ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ); + + } + + } + + floatMatrix assembled_d2MicroGradientFlowdDrivingStressdChi( 3 * 27, floatVector( 9, 0 ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( d2MacroFlowdDrivingStressdChi, floatMatrix( 9, floatVector( 9, 0 ) ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( d2MicroFlowdDrivingStressdChi, floatMatrix( 9, floatVector( 9, 0 ) ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( d2MicroGradientFlowdDrivingStressdChi, assembled_d2MicroGradientFlowdDrivingStressdChi ) ); + + // derivatives w.r.t. previous macro stress + for ( unsigned int i = 0; i < hydra.elasticity.previousPK2Stress.size( ); i++ ){ + + floatVector delta( hydra.elasticity.previousPK2Stress.size( ), 0 ); + + delta[ i ] = eps * std::fabs( hydra.elasticity.previousPK2Stress[ i ] ) + eps; + + hydraBaseMicromorphicMock hydrap( time, deltaTime, temperature, previousTemperature, deformationGradient, previousDeformationGradient, + microDeformation, previousMicroDeformation, gradientMicroDeformation, previousGradientMicroDeformation, + previousStateVariables, parameters, + numConfigurations, numNonLinearSolveStateVariables, + dimension, configuration_unknown_count, + tolr, tola, maxIterations, maxLSIterations, lsAlpha ); + + hydraBaseMicromorphicMock hydram( time, deltaTime, temperature, previousTemperature, deformationGradient, previousDeformationGradient, + microDeformation, previousMicroDeformation, gradientMicroDeformation, previousGradientMicroDeformation, + previousStateVariables, parameters, + numConfigurations, numNonLinearSolveStateVariables, + dimension, configuration_unknown_count, + tolr, tola, maxIterations, maxLSIterations, lsAlpha ); + + hydrap._local_deltaPK2 = delta; + + hydram._local_deltaPK2 = -delta; + + tardigradeHydra::unit_test::hydraBaseTester::updateUnknownVector( hydrap, unknownVector ); + + tardigradeHydra::unit_test::hydraBaseTester::updateUnknownVector( hydram, unknownVector ); + + residualMock Rp( &hydrap, 55, 1, stateVariableIndices, parameters ); + + residualMock Rm( &hydram, 55, 1, stateVariableIndices, parameters ); + + floatVector vp = *Rp.get_previousdMacroFlowdDrivingStress( ); + + floatVector vm = *Rm.get_previousdMacroFlowdDrivingStress( ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + previousd2MacroFlowdDrivingStressdMacroStress[ j ][ i ] = ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ); + + } + + vp = *Rp.get_previousdMicroFlowdDrivingStress( ); + + vm = *Rm.get_previousdMicroFlowdDrivingStress( ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( 0., ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ) ) ); + + } + + vp = tardigradeVectorTools::appendVectors( *Rp.get_previousdMicroGradientFlowdDrivingStress( ) ); + + vm = tardigradeVectorTools::appendVectors( *Rm.get_previousdMicroGradientFlowdDrivingStress( ) ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( 0., ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ) ) ); + + } + + } + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousd2MacroFlowdDrivingStressdMacroStress, *R.get_previousd2MacroFlowdDrivingStressdStress( ) ) ); + + // derivatives w.r.t. previous micro stress + for ( unsigned int i = 0; i < hydra.elasticity.previousSIGMA.size( ); i++ ){ + + floatVector delta( hydra.elasticity.previousSIGMA.size( ), 0 ); + + delta[ i ] = eps * std::fabs( hydra.elasticity.previousSIGMA[ i ] ) + eps; + + hydraBaseMicromorphicMock hydrap( time, deltaTime, temperature, previousTemperature, deformationGradient, previousDeformationGradient, + microDeformation, previousMicroDeformation, gradientMicroDeformation, previousGradientMicroDeformation, + previousStateVariables, parameters, + numConfigurations, numNonLinearSolveStateVariables, + dimension, configuration_unknown_count, + tolr, tola, maxIterations, maxLSIterations, lsAlpha ); + + hydraBaseMicromorphicMock hydram( time, deltaTime, temperature, previousTemperature, deformationGradient, previousDeformationGradient, + microDeformation, previousMicroDeformation, gradientMicroDeformation, previousGradientMicroDeformation, + previousStateVariables, parameters, + numConfigurations, numNonLinearSolveStateVariables, + dimension, configuration_unknown_count, + tolr, tola, maxIterations, maxLSIterations, lsAlpha ); + + hydrap._local_deltaSIGMA = delta; + + hydram._local_deltaSIGMA = -delta; + + tardigradeHydra::unit_test::hydraBaseTester::updateUnknownVector( hydrap, unknownVector ); + + tardigradeHydra::unit_test::hydraBaseTester::updateUnknownVector( hydram, unknownVector ); + + residualMock Rp( &hydrap, 55, 1, stateVariableIndices, parameters ); + + residualMock Rm( &hydram, 55, 1, stateVariableIndices, parameters ); + + floatVector vp = *Rp.get_previousdMacroFlowdDrivingStress( ); + + floatVector vm = *Rm.get_previousdMacroFlowdDrivingStress( ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( 0., ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ) ) ); + + } + + vp = *Rp.get_previousdMicroFlowdDrivingStress( ); + + vm = *Rm.get_previousdMicroFlowdDrivingStress( ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + previousd2MicroFlowdDrivingStressdMicroStress[ j ][ i ] = ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ); + + } + + vp = tardigradeVectorTools::appendVectors( *Rp.get_previousdMicroGradientFlowdDrivingStress( ) ); + + vm = tardigradeVectorTools::appendVectors( *Rm.get_previousdMicroGradientFlowdDrivingStress( ) ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( 0., ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ) ) ); + + } + + } + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousd2MicroFlowdDrivingStressdMicroStress, *R.get_previousd2MicroFlowdDrivingStressdStress( ) ) ); + + // derivatives w.r.t. previous micro gradient stress + for ( unsigned int i = 0; i < hydra.elasticity.previousM.size( ); i++ ){ + + floatVector delta( hydra.elasticity.previousM.size( ), 0 ); + + delta[ i ] = eps * std::fabs( hydra.elasticity.previousM[ i ] ) + eps; + + hydraBaseMicromorphicMock hydrap( time, deltaTime, temperature, previousTemperature, deformationGradient, previousDeformationGradient, + microDeformation, previousMicroDeformation, gradientMicroDeformation, previousGradientMicroDeformation, + previousStateVariables, parameters, + numConfigurations, numNonLinearSolveStateVariables, + dimension, configuration_unknown_count, + tolr, tola, maxIterations, maxLSIterations, lsAlpha ); + + hydraBaseMicromorphicMock hydram( time, deltaTime, temperature, previousTemperature, deformationGradient, previousDeformationGradient, + microDeformation, previousMicroDeformation, gradientMicroDeformation, previousGradientMicroDeformation, + previousStateVariables, parameters, + numConfigurations, numNonLinearSolveStateVariables, + dimension, configuration_unknown_count, + tolr, tola, maxIterations, maxLSIterations, lsAlpha ); + + hydrap._local_deltaM = delta; + + hydram._local_deltaM = -delta; + + tardigradeHydra::unit_test::hydraBaseTester::updateUnknownVector( hydrap, unknownVector ); + + tardigradeHydra::unit_test::hydraBaseTester::updateUnknownVector( hydram, unknownVector ); + + residualMock Rp( &hydrap, 55, 1, stateVariableIndices, parameters ); + + residualMock Rm( &hydram, 55, 1, stateVariableIndices, parameters ); + + floatVector vp = *Rp.get_previousdMacroFlowdDrivingStress( ); + + floatVector vm = *Rm.get_previousdMacroFlowdDrivingStress( ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( 0., ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ) ) ); + + } + + vp = *Rp.get_previousdMicroFlowdDrivingStress( ); + + vm = *Rm.get_previousdMicroFlowdDrivingStress( ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( 0., ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ) ) ); + + } + + vp = tardigradeVectorTools::appendVectors( *Rp.get_previousdMicroGradientFlowdDrivingStress( ) ); + + vm = tardigradeVectorTools::appendVectors( *Rm.get_previousdMicroGradientFlowdDrivingStress( ) ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + previousd2MicroGradientFlowdDrivingStressdHigherOrderStress[ j ][ i ] = ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ); + + } + + } + + floatMatrix assembled_previousd2MicroGradientFlowdDrivingStressdStress( 3 * 27, floatVector( 27, 0 ) ); + + for ( unsigned int i = 0; i < 3; i++ ){ + + for ( unsigned int j = 0; j < 27; j++ ){ + + for ( unsigned int k = 0; k < 27; k++ ){ + + assembled_previousd2MicroGradientFlowdDrivingStressdStress[ 27 * i + j ][ k ] += ( *R.get_previousd2MicroGradientFlowdDrivingStressdStress( ) )[ i ][ 27 * j + k ]; + + } + + } + + } + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousd2MicroGradientFlowdDrivingStressdHigherOrderStress, assembled_previousd2MicroGradientFlowdDrivingStressdStress ) ); + + // derivatives w.r.t. previous, Fn, Chin, gradchin, and ISVs + for ( unsigned int i = 0; i < previousStateVariables.size( ); i++ ){ + + floatVector delta( previousStateVariables.size( ), 0 ); + + delta[ i ] = eps * std::fabs( previousStateVariables[ i ] ) + eps; + + hydraBaseMicromorphicMock hydrap( time, deltaTime, temperature, previousTemperature, deformationGradient, previousDeformationGradient, + microDeformation, previousMicroDeformation, gradientMicroDeformation, previousGradientMicroDeformation, + previousStateVariables + delta, parameters, + numConfigurations, numNonLinearSolveStateVariables, + dimension, configuration_unknown_count, + tolr, tola, maxIterations, maxLSIterations, lsAlpha ); + + hydraBaseMicromorphicMock hydram( time, deltaTime, temperature, previousTemperature, deformationGradient, previousDeformationGradient, + microDeformation, previousMicroDeformation, gradientMicroDeformation, previousGradientMicroDeformation, + previousStateVariables - delta, parameters, + numConfigurations, numNonLinearSolveStateVariables, + dimension, configuration_unknown_count, + tolr, tola, maxIterations, maxLSIterations, lsAlpha ); + + tardigradeHydra::unit_test::hydraBaseTester::updateUnknownVector( hydrap, unknownVector ); + + tardigradeHydra::unit_test::hydraBaseTester::updateUnknownVector( hydram, unknownVector ); + + residualMock Rp( &hydrap, 55, 1, stateVariableIndices, parameters ); + + residualMock Rm( &hydram, 55, 1, stateVariableIndices, parameters ); + + floatVector vp = *Rp.get_previousdMacroFlowdDrivingStress( ); + + floatVector vm = *Rm.get_previousdMacroFlowdDrivingStress( ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + previousd2MacroFlowdDrivingStressdX[ j ][ i ] = ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ); + + } + + vp = *Rp.get_previousdMicroFlowdDrivingStress( ); + + vm = *Rm.get_previousdMicroFlowdDrivingStress( ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + previousd2MicroFlowdDrivingStressdX[ j ][ i ] = ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ); + + } + + vp = tardigradeVectorTools::appendVectors( *Rp.get_previousdMicroGradientFlowdDrivingStress( ) ); + + vm = tardigradeVectorTools::appendVectors( *Rm.get_previousdMicroGradientFlowdDrivingStress( ) ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + previousd2MicroGradientFlowdDrivingStressdX[ j ][ i ] = ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ); + + } + + } + + floatMatrix assembled_previousd2MacroFlowdDrivingStressdX( 9, floatVector( previousStateVariables.size( ), 0 ) ); + + floatMatrix assembled_previousd2MicroFlowdDrivingStressdX( 9, floatVector( previousStateVariables.size( ), 0 ) ); + + floatMatrix assembled_previousd2MicroGradientFlowdDrivingStressdX( 3 * 27, floatVector( previousStateVariables.size( ), 0 ) ); + + for ( unsigned int i = 0; i < 9; i++ ){ + + for ( unsigned int j = 0; j < 9; j++ ){ + + assembled_previousd2MacroFlowdDrivingStressdX[ i ][ j ] = ( *R.get_previousd2MacroFlowdDrivingStressdFn( ) )[ i ][ j ]; + + } + + } + + for ( unsigned int i = 0; i < 9; i++ ){ + + for ( unsigned int j = 0; j < 9; j++ ){ + + assembled_previousd2MicroFlowdDrivingStressdX[ i ][ j ] = ( *R.get_previousd2MicroFlowdDrivingStressdFn( ) )[ i ][ j ]; + + } + + } + + for ( unsigned int i = 0; i < 3; i++ ){ + + for ( unsigned int j = 0; j < 27; j++ ){ + + for ( unsigned int k = 0; k < 9; k++ ){ + + assembled_previousd2MicroGradientFlowdDrivingStressdX[ 27 * i + j ][ k ] = ( *R.get_previousd2MicroGradientFlowdDrivingStressdFn( ) )[ i ][ 9 * j + k ]; + + assembled_previousd2MicroGradientFlowdDrivingStressdX[ 27 * i + j ][ k + 9 ] = ( *R.get_previousd2MicroGradientFlowdDrivingStressdChin( ) )[ i ][ 9 * j + k ]; + + } + + } + + } + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( assembled_previousd2MacroFlowdDrivingStressdX, previousd2MacroFlowdDrivingStressdX ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( assembled_previousd2MicroFlowdDrivingStressdX, previousd2MicroFlowdDrivingStressdX ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( assembled_previousd2MicroGradientFlowdDrivingStressdX, previousd2MicroGradientFlowdDrivingStressdX ) ); + + // derivatives w.r.t. previous F + for ( unsigned int i = 0; i < deformationGradient.size( ); i++ ){ + + floatVector delta( deformationGradient.size( ), 0 ); + + delta[ i ] = eps * std::fabs( deformationGradient[ i ] ) + eps; + + hydraBaseMicromorphicMock hydrap( time, deltaTime, temperature, previousTemperature, deformationGradient, previousDeformationGradient + delta, + microDeformation, previousMicroDeformation, gradientMicroDeformation, previousGradientMicroDeformation, + previousStateVariables, parameters, + numConfigurations, numNonLinearSolveStateVariables, + dimension, configuration_unknown_count, + tolr, tola, maxIterations, maxLSIterations, lsAlpha ); + + hydraBaseMicromorphicMock hydram( time, deltaTime, temperature, previousTemperature, deformationGradient, previousDeformationGradient - delta, + microDeformation, previousMicroDeformation, gradientMicroDeformation, previousGradientMicroDeformation, + previousStateVariables, parameters, + numConfigurations, numNonLinearSolveStateVariables, + dimension, configuration_unknown_count, + tolr, tola, maxIterations, maxLSIterations, lsAlpha ); + + tardigradeHydra::unit_test::hydraBaseTester::updateUnknownVector( hydrap, unknownVector ); + + tardigradeHydra::unit_test::hydraBaseTester::updateUnknownVector( hydram, unknownVector ); + + residualMock Rp( &hydrap, 55, 1, stateVariableIndices, parameters ); + + residualMock Rm( &hydram, 55, 1, stateVariableIndices, parameters ); + + floatVector vp = *Rp.get_previousdMacroFlowdDrivingStress( ); + + floatVector vm = *Rm.get_previousdMacroFlowdDrivingStress( ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + previousd2MacroFlowdDrivingStressdF[ j ][ i ] = ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ); + + } + + vp = *Rp.get_previousdMicroFlowdDrivingStress( ); + + vm = *Rm.get_previousdMicroFlowdDrivingStress( ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + previousd2MicroFlowdDrivingStressdF[ j ][ i ] = ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ); + + } + + vp = tardigradeVectorTools::appendVectors( *Rp.get_previousdMicroGradientFlowdDrivingStress( ) ); + + vm = tardigradeVectorTools::appendVectors( *Rm.get_previousdMicroGradientFlowdDrivingStress( ) ); + + for ( unsigned int j = 0; j < vp.size( ); j++ ){ + + previousd2MicroGradientFlowdDrivingStressdF[ j ][ i ] = ( vp[ j ] - vm[ j ] ) / ( 2 * delta[ i ] ); + + } + + } + + floatMatrix assembled_previousd2MicroGradientFlowdDrivingStressdF( 3 * 27, floatVector( 9, 0 ) ); + + for ( unsigned int i = 0; i < 3; i++ ){ + + for ( unsigned int j = 0; j < 27; j++ ){ + + for ( unsigned int k = 0; k < 9; k++ ){ + + assembled_previousd2MicroGradientFlowdDrivingStressdF[ 27 * i + j ][ k ] = ( *R.get_previousd2MicroGradientFlowdDrivingStressdF( ) )[ i ][ 9 * j + k ]; + + } + + } + + } + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousd2MacroFlowdDrivingStressdF, *R.get_previousd2MacroFlowdDrivingStressdF( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousd2MicroFlowdDrivingStressdF, *R.get_previousd2MicroFlowdDrivingStressdF( ) ) ); + + BOOST_CHECK( tardigradeVectorTools::fuzzyEquals( previousd2MicroGradientFlowdDrivingStressdF, assembled_previousd2MicroGradientFlowdDrivingStressdF ) ); } From cae19f1a7aa2952ab4ffbd6d16a2734b310505ca Mon Sep 17 00:00:00 2001 From: Nathan Miller Date: Mon, 8 Jan 2024 23:54:36 -0700 Subject: [PATCH 15/15] DOC: Updated the changelog --- docs/sphinx/changelog.rst | 1 + 1 file changed, 1 insertion(+) diff --git a/docs/sphinx/changelog.rst b/docs/sphinx/changelog.rst index 97fa8e2d..a4101190 100644 --- a/docs/sphinx/changelog.rst +++ b/docs/sphinx/changelog.rst @@ -22,6 +22,7 @@ Internal Changes - Added the decomposition of the parameter vector (:pull:`39`). By `Nathan Miller`_. - Added the extraction of the nonlinear state variables (:pull:`40`). By `Nathan Miller`_. - Added the calculation of the cohesion (:pull:`41`). By `Nathan Miller`_. +- Added the calculation of the required quantities from the flow potential (:pull:`42`). By `Nathan Miller`_. ****************** 0.3.0 (01-03-2023)