feat: plastic strain output #3384
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
|
|
@@ -115,6 +115,11 @@ class ElasticIsotropicPressureDependentUpdates : public SolidBaseUpdates | |
| localIndex const q, | ||
| real64 ( &elasticStrain )[6] ) const override final; | ||
|
|
||
| GEOS_HOST_DEVICE | ||
| virtual void getElasticStrainInc( localIndex const k, | ||
| localIndex const q, | ||
| real64 ( &elasticStrainInc )[6] ) const override final; | ||
|
|
||
| GEOS_HOST_DEVICE | ||
| virtual void viscousStateUpdate( localIndex const k, | ||
| localIndex const q, | ||
|
|
@@ -223,6 +228,67 @@ void ElasticIsotropicPressureDependentUpdates::getElasticStrain( localIndex cons | |
| } | ||
|
|
||
|
|
||
| GEOS_HOST_DEVICE | ||
| inline | ||
| void ElasticIsotropicPressureDependentUpdates::getElasticStrainInc( localIndex const k, | ||
| localIndex const q, | ||
| real64 ( & elasticStrainInc)[6] ) const | ||
| { | ||
| real64 const mu = m_shearModulus[k]; | ||
| real64 const p0 = m_refPressure; | ||
| real64 const eps_v0 = m_refStrainVol; | ||
| real64 const Cr = m_recompressionIndex[k]; | ||
| real64 deviator[6]; | ||
| real64 stress[6]; | ||
ryar9534 marked this conversation as resolved.
Show resolved
Hide resolved
|
||
| real64 P; | ||
| real64 Q; | ||
| real64 elasticStrainVol; | ||
ryar9534 marked this conversation as resolved.
Show resolved
Hide resolved
|
||
| real64 elasticStrainDev; | ||
ryar9534 marked this conversation as resolved.
Show resolved
Hide resolved
|
||
|
|
||
| for( localIndex i=0; i<6; ++i ) | ||
| { | ||
| stress[i] = m_newStress[k][q][i]; | ||
| } | ||
|
|
||
| twoInvariant::stressDecomposition( stress, | ||
| P, | ||
| Q, | ||
| deviator ); | ||
|
|
||
| elasticStrainVol = std::log( P/p0 ) * Cr * (-1.0) + eps_v0; | ||
ryar9534 marked this conversation as resolved.
Show resolved
Hide resolved
|
||
| elasticStrainDev = Q/3./mu; | ||
ryar9534 marked this conversation as resolved.
Show resolved
Hide resolved
|
||
|
|
||
| twoInvariant::strainRecomposition( elasticStrainVol, | ||
| elasticStrainDev, | ||
| deviator, | ||
| elasticStrainInc ); | ||
|
|
||
| real64 oldStrain[6]; | ||
ryar9534 marked this conversation as resolved.
Show resolved
Hide resolved
|
||
| for( localIndex i=0; i<6; ++i ) | ||
| { | ||
| stress[i] = m_oldStress[k][q][i]; | ||
| } | ||
ryar9534 marked this conversation as resolved.
Show resolved
Hide resolved
|
||
|
|
||
| twoInvariant::stressDecomposition( stress, | ||
| P, | ||
| Q, | ||
| deviator ); | ||
|
|
||
| elasticStrainVol = std::log( P/p0 ) * Cr * (-1.0) + eps_v0; | ||
ryar9534 marked this conversation as resolved.
Show resolved
Hide resolved
|
||
| elasticStrainDev = Q/3./mu; | ||
ryar9534 marked this conversation as resolved.
Show resolved
Hide resolved
|
||
|
|
||
| twoInvariant::strainRecomposition( elasticStrainVol, | ||
| elasticStrainDev, | ||
| deviator, | ||
| oldStrain ); | ||
|
|
||
| for (localIndex i = 0; i<6; ++i) | ||
| { | ||
| elasticStrainInc[i] -= oldStrain[i]; | ||
| } | ||
| } | ||
|
|
||
|
|
||
| GEOS_HOST_DEVICE | ||
| inline | ||
| void ElasticIsotropicPressureDependentUpdates::smallStrainUpdate( localIndex const k, | ||
|
|
||
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
|
|
@@ -337,6 +337,24 @@ class SolidBaseUpdates | |
| GEOS_ERROR( "getElasticStrain() not implemented for this model" ); | ||
| } | ||
|
|
||
| /** | ||
| * @brief Return the current elastic strain increment at a given material point (small-strain interface) | ||
| * | ||
| * @param k the element inex | ||
| * @param q the quadrature index | ||
| * @param elasticStrainInc Current elastic strain increment | ||
| */ | ||
| GEOS_HOST_DEVICE | ||
| virtual void getElasticStrainInc( localIndex const k, | ||
| localIndex const q, | ||
| real64 ( & elasticStrainInc )[6] ) const | ||
| { | ||
| GEOS_UNUSED_VAR( k ); | ||
| GEOS_UNUSED_VAR( q ); | ||
| GEOS_UNUSED_VAR( elasticStrainInc ); | ||
| GEOS_ERROR( "getElasticStrainInc() not implemented for this model (called when computing plastic strains)" ); | ||
ryar9534 marked this conversation as resolved.
Show resolved
Hide resolved
|
||
| } | ||
|
|
||
| /** | ||
| * @brief Perform a viscous (rate-dependent) state update | ||
| * | ||
|
|
||
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
|
|
@@ -87,6 +87,14 @@ DECLARE_FIELD( strain, | |
| WRITE_AND_READ, | ||
| "Average strain in cell" ); | ||
|
|
||
| DECLARE_FIELD( plasticStrain, | ||
| "plasticStrain", | ||
| array2dLayoutStrain, | ||
| 0, | ||
| LEVEL_0, | ||
| WRITE_AND_READ, | ||
| "Average plastic strain in cell" ); | ||
|
|
||
| DECLARE_FIELD( incrementalBubbleDisplacement, | ||
| "incrementalBubbleDisplacement", | ||
| array2d< real64 >, | ||
|
|
||
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
|
|
@@ -161,6 +161,8 @@ void SolidMechanicsLagrangianFEM::registerDataOnMesh( Group & meshBodies ) | |
| setConstitutiveNamesCallSuper( subRegion ); | ||
|
|
||
| subRegion.registerField< solidMechanics::strain >( getName() ).setDimLabels( 1, voightLabels ).reference().resizeDimension< 1 >( 6 ); | ||
| subRegion.registerField< solidMechanics::plasticStrain >( getName() ).setDimLabels( 1, voightLabels ).reference().resizeDimension< 1 >( 6 ); | ||
|
|
||
| } ); | ||
|
|
||
| NodeManager & nodes = meshLevel.getNodeManager(); | ||
|
|
@@ -454,6 +456,7 @@ void SolidMechanicsLagrangianFEM::initializePostInitialConditionsPreSubGroups() | |
| m_targetNodes.insert( m_nonSendOrReceiveNodes.begin(), | ||
| m_nonSendOrReceiveNodes.end() ); | ||
|
|
||
|
|
||
| } ); | ||
| } ); | ||
|
|
||
|
|
@@ -944,24 +947,37 @@ void SolidMechanicsLagrangianFEM::implicitStepComplete( real64 const & GEOS_UNUS | |
| { | ||
| string const & solidMaterialName = subRegion.template getReference< string >( viewKeyStruct::solidMaterialNamesString() ); | ||
| SolidBase & constitutiveRelation = getConstitutiveModel< SolidBase >( subRegion, solidMaterialName ); | ||
|
|
||
|
|
||
| solidMechanics::arrayView2dLayoutStrain strain = subRegion.getField< solidMechanics::strain >(); | ||
| solidMechanics::arrayView2dLayoutStrain plasticStrain = subRegion.getField< solidMechanics::plasticStrain >(); | ||
|
|
||
| constitutive::ConstitutivePassThru< SolidBase >::execute(constitutiveRelation, [&] (auto & solidModel) | ||
| { | ||
|
|
||
| using SOLID_TYPE = TYPEOFREF( solidModel ); | ||
|
|
||
| finiteElement::FiniteElementBase & subRegionFE = subRegion.template getReference< finiteElement::FiniteElementBase >( this->getDiscretizationName()); | ||
| finiteElement::FiniteElementDispatchHandler< BASE_FE_TYPES >::dispatch3D( subRegionFE, [&] ( auto const finiteElement ) | ||
| { | ||
| using FE_TYPE = decltype( finiteElement ); | ||
| AverageStrainOverQuadraturePointsKernelFactory::createAndLaunch< FE_TYPE, SOLID_TYPE, parallelDevicePolicy<> >( nodeManager, | ||
| mesh.getEdgeManager(), | ||
| mesh.getFaceManager(), | ||
| subRegion, | ||
| finiteElement, | ||
| solidModel, | ||
| disp, | ||
| uhat, | ||
| strain, | ||
| plasticStrain ); | ||
| } ); | ||
|
|
||
|
|
||
| } ); | ||
|
|
||
| constitutiveRelation.saveConvergedState(); | ||
|
|
||
| } ); | ||
| } ); | ||
|
|
||
|
|
||
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
|
|
@@ -23,6 +23,7 @@ | |
| #include "common/DataTypes.hpp" | ||
| #include "common/GEOS_RAJA_Interface.hpp" | ||
| #include "finiteElement/FiniteElementDispatch.hpp" | ||
| #include "constitutive/ConstitutivePassThru.hpp" | ||
| #include "mesh/CellElementSubRegion.hpp" | ||
| #include "mesh/utilities/AverageOverQuadraturePointsKernel.hpp" | ||
| #include "physicsSolvers/solidMechanics/SolidMechanicsFields.hpp" | ||
|
|
@@ -33,17 +34,18 @@ namespace geos | |
| * @class AverageStrainOverQuadraturePoints | ||
| * @tparam SUBREGION_TYPE the subRegion type | ||
| * @tparam FE_TYPE the finite element type | ||
| * @tparam SOLID_TYPE the solid mechanics constitutuve type | ||
| */ | ||
| template< typename FE_TYPE, | ||
| typename SOLID_TYPE > | ||
| class AverageStrainOverQuadraturePoints : | ||
| public AverageOverQuadraturePointsBase< CellElementSubRegion, | ||
| FE_TYPE > | ||
| { | ||
| public: | ||
|
|
||
| /// Alias for the base class; | ||
| using Base = AverageOverQuadraturePointsBase< CellElementSubRegion, | ||
| FE_TYPE >; | ||
|
|
||
| using Base::m_elementVolume; | ||
|
|
@@ -63,24 +65,31 @@ class AverageStrainOverQuadraturePoints : | |
| AverageStrainOverQuadraturePoints( NodeManager & nodeManager, | ||
| EdgeManager const & edgeManager, | ||
| FaceManager const & faceManager, | ||
| CellElementSubRegion const & elementSubRegion, | ||
| FE_TYPE const & finiteElementSpace, | ||
| SOLID_TYPE const & solidModel, | ||
| fields::solidMechanics::arrayViewConst2dLayoutTotalDisplacement const displacement, | ||
| fields::solidMechanics::arrayViewConst2dLayoutIncrDisplacement const displacementInc, | ||
| fields::solidMechanics::arrayView2dLayoutStrain const avgStrain, | ||
| fields::solidMechanics::arrayView2dLayoutStrain const avgPlasticStrain): | ||
| Base( nodeManager, | ||
| edgeManager, | ||
| faceManager, | ||
| elementSubRegion, | ||
| finiteElementSpace ), | ||
| m_solidUpdate(solidModel.createKernelUpdates()), | ||
| m_displacement( displacement ), | ||
| m_displacementInc( displacementInc ), | ||
| m_avgStrain( avgStrain ), | ||
| m_avgPlasticStrain( avgPlasticStrain ) | ||
| {} | ||
|
|
||
| /** | ||
| * @copydoc finiteElement::KernelBase::StackVariables | ||
| */ | ||
| struct StackVariables : Base::StackVariables | ||
| {real64 uLocal[FE_TYPE::maxSupportPoints][3]; | ||
| real64 uHatLocal[FE_TYPE::maxSupportPoints][3]; }; | ||
|
|
||
| /** | ||
| * @brief Performs the setup phase for the kernel. | ||
|
|
@@ -99,12 +108,14 @@ class AverageStrainOverQuadraturePoints : | |
| for( int i = 0; i < 3; ++i ) | ||
| { | ||
| stack.uLocal[a][i] = m_displacement[localNodeIndex][i]; | ||
| stack.uHatLocal[a][i] = m_displacementInc[localNodeIndex][i]; | ||
| } | ||
| } | ||
|
|
||
| for( int icomp = 0; icomp < 6; ++icomp ) | ||
| { | ||
| m_avgStrain[k][icomp] = 0.0; | ||
| //m_avgPlasticStrain[k][icomp] = 0.0; | ||
| } | ||
| } | ||
|
|
||
|
|
@@ -124,11 +135,18 @@ class AverageStrainOverQuadraturePoints : | |
| real64 dNdX[ FE_TYPE::maxSupportPoints ][3]; | ||
| real64 const detJxW = m_finiteElementSpace.template getGradN< FE_TYPE >( k, q, stack.xLocal, stack.feStack, dNdX ); | ||
| real64 strain[6] = {0.0}; | ||
| real64 strainInc[6] = {0.0}; | ||
| FE_TYPE::symmetricGradient( dNdX, stack.uLocal, strain ); | ||
| FE_TYPE::symmetricGradient( dNdX, stack.uHatLocal, strainInc ); | ||
|
|
||
| real64 elasticStrainInc[6] = {0.0}; | ||
| m_solidUpdate.getElasticStrainInc(k, q, elasticStrainInc); | ||
|
|
||
| for( int icomp = 0; icomp < 6; ++icomp ) | ||
| { | ||
| m_avgStrain[k][icomp] += detJxW*strain[icomp]/m_elementVolume[k]; | ||
| m_avgPlasticStrain[k][icomp] += detJxW*(strainInc[icomp] - elasticStrainInc[icomp])/m_elementVolume[k]; | ||
| //m_avgPlasticStrain[k][icomp] += detJxW*(elasticStrainInc[icomp])/m_elementVolume[k]; | ||
| } | ||
| } | ||
|
|
||
|
|
@@ -160,12 +178,21 @@ class AverageStrainOverQuadraturePoints : | |
|
|
||
| protected: | ||
|
|
||
| /// The material | ||
| typename SOLID_TYPE::KernelWrapper const m_solidUpdate; | ||
|
|
||
| /// The displacement solution | ||
| fields::solidMechanics::arrayViewConst2dLayoutTotalDisplacement const m_displacement; | ||
|
|
||
| /// The displacement increment | ||
| fields::solidMechanics::arrayViewConst2dLayoutIncrDisplacement const m_displacementInc; | ||
|
|
||
| /// The average strain | ||
| fields::solidMechanics::arrayView2dLayoutStrain const m_avgStrain; | ||
|
|
||
| /// The average plastic strain | ||
| fields::solidMechanics::arrayView2dLayoutStrain const m_avgPlasticStrain; | ||
|
|
||
| }; | ||
|
|
||
|
|
||
|
|
@@ -182,6 +209,7 @@ class AverageStrainOverQuadraturePointsKernelFactory | |
| * @brief Create a new kernel and launch | ||
| * @tparam SUBREGION_TYPE the subRegion type | ||
| * @tparam FE_TYPE the finite element type | ||
| * @tparam SOLID_TYPE the constitutive type | ||
| * @tparam POLICY the kernel policy | ||
| * @param nodeManager the node manager | ||
| * @param edgeManager the edge manager | ||
|
|
@@ -191,23 +219,26 @@ class AverageStrainOverQuadraturePointsKernelFactory | |
| * @param property the property at quadrature points | ||
| * @param averageProperty the property averaged over quadrature points | ||
| */ | ||
| template< typename FE_TYPE, | ||
| typename SOLID_TYPE, | ||
| typename POLICY > | ||
| static void | ||
| createAndLaunch( NodeManager & nodeManager, | ||
| EdgeManager const & edgeManager, | ||
| FaceManager const & faceManager, | ||
| CellElementSubRegion const & elementSubRegion, | ||
| FE_TYPE const & finiteElementSpace, | ||
| SOLID_TYPE const & solidModel, | ||
| fields::solidMechanics::arrayViewConst2dLayoutTotalDisplacement const displacement, | ||
| fields::solidMechanics::arrayViewConst2dLayoutIncrDisplacement const displacementInc, | ||
| fields::solidMechanics::arrayView2dLayoutStrain const avgStrain, | ||
| fields::solidMechanics::arrayView2dLayoutStrain const avgPlasticStrain) | ||
| { | ||
| AverageStrainOverQuadraturePoints< FE_TYPE, SOLID_TYPE > | ||
| kernel( nodeManager, edgeManager, faceManager, elementSubRegion, finiteElementSpace, | ||
| solidModel, displacement, displacementInc, avgStrain, avgPlasticStrain ); | ||
|
|
||
| AverageStrainOverQuadraturePoints< FE_TYPE, SOLID_TYPE >::template | ||
| kernelLaunch< POLICY >( elementSubRegion.size(), kernel ); | ||
| } | ||
| }; | ||
|
|
||
Oops, something went wrong.
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Since this isn't really a "getter"