From 3f8678c7c4e25f60ef2a8c0afcfc1b5d02bb3069 Mon Sep 17 00:00:00 2001
From: Pratheek Shanthraj
Date: Tue, 9 Dec 2014 18:23:48 +0000
Subject: [PATCH] plastic dot states evaluated at effective stress not damaged
stress
---
code/constitutive.f90 | 73 +++++++++++++++++++---------------
code/damage_anisoBrittle.f90 | 72 ---------------------------------
code/damage_isoDuctile.f90 | 11 ++---
code/plastic_dislokmc.f90 | 17 ++++----
code/plastic_dislotwin.f90 | 16 ++++----
code/plastic_j2.f90 | 14 ++++---
code/plastic_nonlocal.f90 | 20 ++++++----
code/plastic_phenopowerlaw.f90 | 19 +++++----
code/plastic_titanmod.f90 | 8 ++--
9 files changed, 102 insertions(+), 148 deletions(-)
diff --git a/code/constitutive.f90 b/code/constitutive.f90
index e292b09c5..b876c4e6d 100644
--- a/code/constitutive.f90
+++ b/code/constitutive.f90
@@ -705,37 +705,37 @@ subroutine constitutive_LpAndItsTangent(Lp, dLp_dTstar, Tstar_v, ipc, ip, el)
case (PLASTICITY_J2_ID)
nSlip = 1_pInt
call plastic_j2_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v, &
- constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el), &
- ipc,ip,el)
+ nSlip,constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el), &
+ ipc,ip,el)
case (PLASTICITY_PHENOPOWERLAW_ID)
nSlip = plastic_phenopowerlaw_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))
call plastic_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v, &
- constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el), &
- ipc,ip,el)
+ nSlip,constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el), &
+ ipc,ip,el)
case (PLASTICITY_NONLOCAL_ID)
nSlip = totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))
call plastic_nonlocal_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v, &
- constitutive_getTemperature(ipc,ip,el), &
- constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el), &
- ipc,ip,el)
+ constitutive_getTemperature(ipc,ip,el), &
+ nSlip,constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el), &
+ ipc,ip,el)
case (PLASTICITY_DISLOTWIN_ID)
nSlip = plastic_dislotwin_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))
call plastic_dislotwin_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v, &
- constitutive_getTemperature(ipc,ip,el), &
- constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el), &
- ipc,ip,el)
+ constitutive_getTemperature(ipc,ip,el), &
+ nSlip,constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el), &
+ ipc,ip,el)
case (PLASTICITY_DISLOKMC_ID)
nSlip = plastic_dislokmc_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))
call plastic_dislokmc_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v, &
- constitutive_getTemperature(ipc,ip,el), &
- constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el), &
- ipc,ip,el)
+ constitutive_getTemperature(ipc,ip,el), &
+ nSlip,constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el), &
+ ipc,ip,el)
case (PLASTICITY_TITANMOD_ID)
nSlip = plastic_titanmod_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))
call plastic_titanmod_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v, &
- constitutive_getTemperature(ipc,ip,el), &
- constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el), &
- ipc,ip,el)
+ constitutive_getTemperature(ipc,ip,el), &
+ nSlip,constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el), &
+ ipc,ip,el)
end select
@@ -1072,6 +1072,7 @@ subroutine constitutive_collectDotState(Tstar_v, Lp, FeArray, FpArray, subdt, su
mesh_maxNips
use material, only: &
phase_plasticity, &
+ phase_plasticityInstance, &
phase_damage, &
phase_vacancy, &
material_phase, &
@@ -1083,23 +1084,25 @@ subroutine constitutive_collectDotState(Tstar_v, Lp, FeArray, FpArray, subdt, su
PLASTICITY_dislokmc_ID, &
PLASTICITY_titanmod_ID, &
PLASTICITY_nonlocal_ID, &
- LOCAL_DAMAGE_anisoBrittle_ID, &
LOCAL_DAMAGE_gurson_ID, &
LOCAL_VACANCY_generation_ID
use plastic_j2, only: &
plastic_j2_dotState
use plastic_phenopowerlaw, only: &
- plastic_phenopowerlaw_dotState
+ plastic_phenopowerlaw_dotState, &
+ plastic_phenopowerlaw_totalNslip
use plastic_dislotwin, only: &
- plastic_dislotwin_dotState
+ plastic_dislotwin_dotState, &
+ plastic_dislotwin_totalNslip
use plastic_dislokmc, only: &
- plastic_dislokmc_dotState
+ plastic_dislokmc_dotState, &
+ plastic_dislokmc_totalNslip
use plastic_titanmod, only: &
- plastic_titanmod_dotState
+ plastic_titanmod_dotState, &
+ plastic_titanmod_totalNslip
use plastic_nonlocal, only: &
- plastic_nonlocal_dotState
- use damage_anisoBrittle, only: &
- damage_anisoBrittle_dotState
+ plastic_nonlocal_dotState, &
+ totalNslip
use damage_gurson, only: &
damage_gurson_dotState
use vacancy_generation, only: &
@@ -1135,23 +1138,29 @@ subroutine constitutive_collectDotState(Tstar_v, Lp, FeArray, FpArray, subdt, su
select case (phase_plasticity(material_phase(ipc,ip,el)))
case (PLASTICITY_J2_ID)
- call plastic_j2_dotState (Tstar_v,ipc,ip,el)
+ nSlip = 1_pInt
+ call plastic_j2_dotState (Tstar_v,nSlip,constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el),ipc,ip,el)
case (PLASTICITY_PHENOPOWERLAW_ID)
- call plastic_phenopowerlaw_dotState(Tstar_v,ipc,ip,el)
+ nSlip = plastic_phenopowerlaw_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))
+ call plastic_phenopowerlaw_dotState(Tstar_v,nSlip,constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el),ipc,ip,el)
case (PLASTICITY_DISLOTWIN_ID)
- call plastic_dislotwin_dotState (Tstar_v,constitutive_getTemperature(ipc,ip,el),ipc,ip,el)
+ nSlip = plastic_dislotwin_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))
+ call plastic_dislotwin_dotState (Tstar_v,constitutive_getTemperature(ipc,ip,el), &
+ nSlip,constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el),ipc,ip,el)
case (PLASTICITY_DISLOKMC_ID)
- call plastic_dislokmc_dotState (Tstar_v,constitutive_getTemperature(ipc,ip,el),ipc,ip,el)
+ nSlip = plastic_dislokmc_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))
+ call plastic_dislokmc_dotState (Tstar_v,constitutive_getTemperature(ipc,ip,el), &
+ nSlip,constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el),ipc,ip,el)
case (PLASTICITY_TITANMOD_ID)
- call plastic_titanmod_dotState (Tstar_v,constitutive_getTemperature(ipc,ip,el),ipc,ip,el)
+ call plastic_titanmod_dotState (Tstar_v,constitutive_getTemperature(ipc,ip,el), &
+ ipc,ip,el)
case (PLASTICITY_NONLOCAL_ID)
+ nSlip = totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))
call plastic_nonlocal_dotState (Tstar_v,FeArray,FpArray,constitutive_getTemperature(ipc,ip,el), &
- subdt,subfracArray,ip,el)
+ nSlip,constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el),subdt,subfracArray,ip,el)
end select
select case (phase_damage(material_phase(ipc,ip,el)))
- case (LOCAL_DAMAGE_anisoBrittle_ID)
- call damage_anisoBrittle_dotState(Tstar_v, ipc, ip, el)
case (LOCAL_DAMAGE_gurson_ID)
call damage_gurson_dotState(Tstar_v, Lp, ipc, ip, el)
end select
diff --git a/code/damage_anisoBrittle.f90 b/code/damage_anisoBrittle.f90
index 5b2521c78..56bffbb52 100644
--- a/code/damage_anisoBrittle.f90
+++ b/code/damage_anisoBrittle.f90
@@ -54,7 +54,6 @@ module damage_anisoBrittle
damage_anisoBrittle_init, &
damage_anisoBrittle_stateInit, &
damage_anisoBrittle_aTolState, &
- damage_anisoBrittle_dotState, &
damage_anisoBrittle_microstructure, &
damage_anisoBrittle_LdAndItsTangent, &
damage_anisoBrittle_getFd, &
@@ -340,77 +339,6 @@ subroutine damage_anisoBrittle_aTolState(phase,instance)
damageState(phase)%aTolState = tempTol
end subroutine damage_anisoBrittle_aTolState
-!--------------------------------------------------------------------------------------------------
-!> @brief calculates derived quantities from state
-!--------------------------------------------------------------------------------------------------
-subroutine damage_anisoBrittle_dotState(Tstar_v,ipc, ip, el)
- use material, only: &
- mappingConstitutive, &
- phase_damageInstance, &
- damageState
- use lattice, only: &
- lattice_Scleavage_v, &
- lattice_maxNcleavageFamily, &
- lattice_NcleavageSystem, &
- lattice_DamageMobility
-
- implicit none
- integer(pInt), intent(in) :: &
- ipc, & !< component-ID of integration point
- ip, & !< integration point
- el !< element
- real(pReal), intent(in), dimension(6) :: &
- Tstar_v !< 2nd Piola Kirchhoff stress tensor (Mandel)
- integer(pInt) :: &
- phase, &
- constituent, &
- instance, &
- f, i, index_d, index_o, index_myFamily
- real(pReal) :: &
- traction_d, traction_t, traction_n, traction_crit, &
- udotd, udott, udotn, &
- localDamage
-
- phase = mappingConstitutive(2,ipc,ip,el)
- constituent = mappingConstitutive(1,ipc,ip,el)
- instance = phase_damageInstance(phase)
-
- localDamage = minval(damageState(phase)%state(2+ damage_anisoBrittle_totalNcleavage(instance): &
- 1+2*damage_anisoBrittle_totalNcleavage(instance),constituent))
- damageState(phase)%dotState(1,constituent) = &
- (localDamage - damageState(phase)%state(1,constituent))/lattice_DamageMobility(phase)
-
- index_o = 2_pInt
- index_d = 2_pInt + damage_anisoBrittle_totalNcleavage(instance)
- do f = 1_pInt,lattice_maxNcleavageFamily
- index_myFamily = sum(lattice_NcleavageSystem(1:f-1_pInt,phase)) ! at which index starts my family
- do i = 1_pInt,damage_anisoBrittle_Ncleavage(f,instance) ! process each (active) cleavage system in family
- traction_d = dot_product(Tstar_v,lattice_Scleavage_v(1:6,1,index_myFamily+i,phase))
- traction_t = dot_product(Tstar_v,lattice_Scleavage_v(1:6,2,index_myFamily+i,phase))
- traction_n = dot_product(Tstar_v,lattice_Scleavage_v(1:6,3,index_myFamily+i,phase))
- traction_crit = damage_anisoBrittle_critLoad(f,instance)* &
- damageState(phase)%state0(index_d,constituent)* &
- damageState(phase)%state0(index_d,constituent)
-
- udotd = &
- damage_anisoBrittle_sdot_0(instance)* &
- (abs(traction_d)/traction_crit)**damage_anisoBrittle_N(instance)
- udott = &
- damage_anisoBrittle_sdot_0(instance)* &
- (abs(traction_t)/traction_crit)**damage_anisoBrittle_N(instance)
- udotn = &
- damage_anisoBrittle_sdot_0(instance)* &
- (max(0.0_pReal,traction_n)/traction_crit)**damage_anisoBrittle_N(instance)
-
- damageState(phase)%dotState(index_o,constituent) = &
- (udotd + udott + udotn)/damage_anisoBrittle_critDisp(f,instance)
-
- index_d = index_d + 1_pInt; index_o = index_o + 1_pInt
- enddo
- enddo
-
-end subroutine damage_anisoBrittle_dotState
-
!--------------------------------------------------------------------------------------------------
!> @brief calculates derived quantities from state
!--------------------------------------------------------------------------------------------------
diff --git a/code/damage_isoDuctile.f90 b/code/damage_isoDuctile.f90
index 12ea9c7f5..43b82b780 100644
--- a/code/damage_isoDuctile.f90
+++ b/code/damage_isoDuctile.f90
@@ -265,12 +265,12 @@ end subroutine damage_isoDuctile_aTolState
!> @brief calculates derived quantities from state
!--------------------------------------------------------------------------------------------------
subroutine damage_isoDuctile_microstructure(nSlip,accumulatedSlip,subdt,ipc, ip, el)
+ use numerics, only: &
+ residualStiffness
use material, only: &
phase_damageInstance, &
mappingConstitutive, &
damageState
- use math, only: &
- math_norm33
use lattice, only: &
lattice_DamageMobility
@@ -292,13 +292,14 @@ subroutine damage_isoDuctile_microstructure(nSlip,accumulatedSlip,subdt,ipc, ip,
instance = phase_damageInstance(phase)
damageState(phase)%state(2,constituent) = &
- min(damageState(phase)%state0(2,constituent), &
- damage_isoDuctile_critPlasticStrain(instance)/sum(accumulatedSlip))
+ max(residualStiffness, &
+ min(damageState(phase)%state0(2,constituent), &
+ damage_isoDuctile_critPlasticStrain(instance)/sum(accumulatedSlip)))
damageState(phase)%state(1,constituent) = &
damageState(phase)%state(2,constituent) + &
(damageState(phase)%subState0(1,constituent) - damageState(phase)%state(2,constituent))* &
- exp(-subdt/lattice_DamageMobility(phase))
+ exp(-subdt/(damageState(phase)%state(2,constituent)*lattice_DamageMobility(phase)))
end subroutine damage_isoDuctile_microstructure
diff --git a/code/plastic_dislokmc.f90 b/code/plastic_dislokmc.f90
index 0f52a7eee..86a629229 100644
--- a/code/plastic_dislokmc.f90
+++ b/code/plastic_dislokmc.f90
@@ -1111,7 +1111,8 @@ end subroutine plastic_dislokmc_microstructure
!--------------------------------------------------------------------------------------------------
!> @brief calculates plastic velocity gradient and its tangent
!--------------------------------------------------------------------------------------------------
-subroutine plastic_dislokmc_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,Temperature,slipDamage,ipc,ip,el)
+subroutine plastic_dislokmc_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,Temperature,nSlipDamage,slipDamage, &
+ ipc,ip,el)
use prec, only: &
tol_math_check
use math, only: &
@@ -1143,13 +1144,10 @@ subroutine plastic_dislokmc_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,Temperature,
LATTICE_fcc_ID
implicit none
- integer(pInt), intent(in) :: ipc,ip,el
+ integer(pInt), intent(in) :: nSlipDamage,ipc,ip,el
real(pReal), intent(in) :: Temperature
real(pReal), dimension(6), intent(in) :: Tstar_v
- real(pReal), &
- dimension(plastic_dislokmc_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))), &
- intent(in) :: &
- slipDamage
+ real(pReal), dimension(nSlipDamage), intent(in) :: slipDamage
real(pReal), dimension(3,3), intent(out) :: Lp
real(pReal), dimension(9,9), intent(out) :: dLp_dTstar99
@@ -1357,7 +1355,7 @@ end subroutine plastic_dislokmc_LpAndItsTangent
!--------------------------------------------------------------------------------------------------
!> @brief calculates the rate of change of microstructure
!--------------------------------------------------------------------------------------------------
-subroutine plastic_dislokmc_dotState(Tstar_v,Temperature,ipc,ip,el)
+subroutine plastic_dislokmc_dotState(Tstar_v,Temperature,nSlipDamage,slipDamage,ipc,ip,el)
use prec, only: &
tol_math_check
use math, only: &
@@ -1388,9 +1386,12 @@ subroutine plastic_dislokmc_dotState(Tstar_v,Temperature,ipc,ip,el)
real(pReal), intent(in) :: &
temperature !< temperature at integration point
integer(pInt), intent(in) :: &
+ nSlipDamage, &
ipc, & !< component-ID of integration point
ip, & !< integration point
el !< element
+ real(pReal), dimension(nSlipDamage), intent(in) :: &
+ slipDamage
integer(pInt) :: instance,ns,nt,f,i,j,k,index_myFamily,s1,s2, &
ph, &
@@ -1456,6 +1457,8 @@ subroutine plastic_dislokmc_dotState(Tstar_v,Temperature,ipc,ip,el)
tau_slip_neg = tau_slip_neg + plastic_dislokmc_nonSchmidCoeff(k,instance)* &
dot_product(Tstar_v,lattice_Sslip_v(1:6,2*k+1,index_myFamily+i,ph))
enddo nonSchmidSystems
+ tau_slip_pos = tau_slip_pos/slipDamage(j)
+ tau_slip_neg = tau_slip_pos/slipDamage(j)
significantPositiveStress: if((abs(tau_slip_pos)-plasticState(ph)%state(6*ns+4*nt+j, of)) > tol_math_check) then
!* Stress ratios
diff --git a/code/plastic_dislotwin.f90 b/code/plastic_dislotwin.f90
index 385bb7993..fed7fda46 100644
--- a/code/plastic_dislotwin.f90
+++ b/code/plastic_dislotwin.f90
@@ -1353,7 +1353,7 @@ end subroutine plastic_dislotwin_microstructure
!--------------------------------------------------------------------------------------------------
!> @brief calculates plastic velocity gradient and its tangent
!--------------------------------------------------------------------------------------------------
-subroutine plastic_dislotwin_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,Temperature,slipDamage,ipc,ip,el)
+subroutine plastic_dislotwin_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,Temperature,nSlipDamage,slipDamage,ipc,ip,el)
use prec, only: &
tol_math_check
use math, only: &
@@ -1389,13 +1389,10 @@ subroutine plastic_dislotwin_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,Temperature
LATTICE_fcc_ID
implicit none
- integer(pInt), intent(in) :: ipc,ip,el
+ integer(pInt), intent(in) :: nSlipDamage,ipc,ip,el
real(pReal), intent(in) :: Temperature
real(pReal), dimension(6), intent(in) :: Tstar_v
- real(pReal), &
- dimension(plastic_dislotwin_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))), &
- intent(in) :: &
- slipDamage
+ real(pReal), dimension(nSlipDamage), intent(in) :: slipDamage
real(pReal), dimension(3,3), intent(out) :: Lp
real(pReal), dimension(9,9), intent(out) :: dLp_dTstar99
@@ -1661,7 +1658,7 @@ end subroutine plastic_dislotwin_LpAndItsTangent
!--------------------------------------------------------------------------------------------------
!> @brief calculates the rate of change of microstructure
!--------------------------------------------------------------------------------------------------
-subroutine plastic_dislotwin_dotState(Tstar_v,Temperature,ipc,ip,el)
+subroutine plastic_dislotwin_dotState(Tstar_v,Temperature,nSlipDamage,slipDamage,ipc,ip,el)
use prec, only: &
tol_math_check
use math, only: &
@@ -1699,9 +1696,12 @@ subroutine plastic_dislotwin_dotState(Tstar_v,Temperature,ipc,ip,el)
real(pReal), intent(in) :: &
temperature !< temperature at integration point
integer(pInt), intent(in) :: &
+ nSlipDamage, &
ipc, & !< component-ID of integration point
ip, & !< integration point
el !< element
+ real(pReal), dimension(nSlipDamage), intent(in) :: &
+ slipDamage
integer(pInt) :: instance,ns,nt,nr,f,i,j,index_myFamily,s1,s2,a,b,sa,sb,ssa,ssb, &
ph, &
@@ -1741,7 +1741,7 @@ subroutine plastic_dislotwin_dotState(Tstar_v,Temperature,ipc,ip,el)
j = j+1_pInt
!* Resolved shear stress on slip system
- tau_slip(j) = dot_product(Tstar_v,lattice_Sslip_v(:,1,index_myFamily+i,ph))
+ tau_slip(j) = dot_product(Tstar_v,lattice_Sslip_v(:,1,index_myFamily+i,ph))/slipDamage(j)
if((abs(tau_slip(j))-plasticState(ph)%state(7*ns+4*nt+2*nr+j, of)) > tol_math_check) then
!* Stress ratios
diff --git a/code/plastic_j2.f90 b/code/plastic_j2.f90
index bb938cff6..5d9244d71 100644
--- a/code/plastic_j2.f90
+++ b/code/plastic_j2.f90
@@ -359,7 +359,7 @@ end subroutine plastic_j2_init
!--------------------------------------------------------------------------------------------------
!> @brief calculates plastic velocity gradient and its tangent
!--------------------------------------------------------------------------------------------------
-subroutine plastic_j2_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,slipDamage,ipc,ip,el)
+subroutine plastic_j2_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,nSlipDamage,slipDamage,ipc,ip,el)
use math, only: &
math_mul6x6, &
math_Mandel6to33, &
@@ -384,12 +384,13 @@ subroutine plastic_j2_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,slipDamage,ipc,ip,
real(pReal), dimension(6), intent(in) :: &
Tstar_v !< 2nd Piola Kirchhoff stress tensor in Mandel notation
- real(pReal), dimension(1), intent(in) :: &
- slipDamage
integer(pInt), intent(in) :: &
+ nSlipDamage, &
ipc, & !< component-ID of integration point
ip, & !< integration point
el !< element
+ real(pReal), dimension(nSlipDamage), intent(in) :: &
+ slipDamage
real(pReal), dimension(3,3) :: &
Tstar_dev_33 !< deviatoric part of the 2nd Piola Kirchhoff stress tensor as 2nd order tensor
@@ -437,7 +438,7 @@ end subroutine plastic_j2_LpAndItsTangent
!--------------------------------------------------------------------------------------------------
!> @brief calculates the rate of change of microstructure
!--------------------------------------------------------------------------------------------------
-subroutine plastic_j2_dotState(Tstar_v,ipc,ip,el)
+subroutine plastic_j2_dotState(Tstar_v,nSlipDamage,slipDamage,ipc,ip,el)
use math, only: &
math_mul6x6
use mesh, only: &
@@ -454,9 +455,12 @@ subroutine plastic_j2_dotState(Tstar_v,ipc,ip,el)
real(pReal), dimension(6), intent(in):: &
Tstar_v !< 2nd Piola Kirchhoff stress tensor in Mandel notation
integer(pInt), intent(in) :: &
+ nSlipDamage, &
ipc, & !< component-ID of integration point
ip, & !< integration point
el !< element
+ real(pReal), dimension(nSlipDamage), intent(in) :: &
+ slipDamage
real(pReal), dimension(6) :: &
Tstar_dev_v !< deviatoric part of the 2nd Piola Kirchhoff stress tensor in Mandel notation
real(pReal) :: &
@@ -483,7 +487,7 @@ subroutine plastic_j2_dotState(Tstar_v,ipc,ip,el)
! strain rate
gamma_dot = plastic_j2_gdot0(instance) * ( sqrt(1.5_pReal) * norm_Tstar_dev &
/ &!-----------------------------------------------------------------------------------
- (plastic_j2_fTaylor(instance)*plasticState(ph)%state(1,of)) )**plastic_j2_n(instance)
+ (slipDamage(1)*plastic_j2_fTaylor(instance)*plasticState(ph)%state(1,of)) )**plastic_j2_n(instance)
!--------------------------------------------------------------------------------------------------
! hardening coefficient
diff --git a/code/plastic_nonlocal.f90 b/code/plastic_nonlocal.f90
index bad91a8b9..6a360dd00 100644
--- a/code/plastic_nonlocal.f90
+++ b/code/plastic_nonlocal.f90
@@ -2016,7 +2016,8 @@ end subroutine plastic_nonlocal_kinetics
!--------------------------------------------------------------------------------------------------
!> @brief calculates plastic velocity gradient and its tangent
!--------------------------------------------------------------------------------------------------
-subroutine plastic_nonlocal_LpAndItsTangent(Lp, dLp_dTstar99, Tstar_v, Temperature, slipDamage, ipc, ip, el)
+subroutine plastic_nonlocal_LpAndItsTangent(Lp, dLp_dTstar99, Tstar_v, Temperature, nSlipDamage, slipDamage, &
+ ipc, ip, el)
use math, only: math_Plain3333to99, &
math_mul6x6, &
@@ -2042,13 +2043,14 @@ use mesh, only: mesh_ipVolume
implicit none
!*** input variables
-integer(pInt), intent(in) :: ipc, &
+integer(pInt), intent(in) :: nSlipDamage, &
+ ipc, &
ip, & !< current integration point
el !< current element number
real(pReal), intent(in) :: Temperature !< temperature
real(pReal), dimension(6), intent(in) :: Tstar_v !< 2nd Piola-Kirchhoff stress in Mandel notation
-real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))), intent(in) :: &
- slipDamage
+real(pReal), dimension(nSlipDamage), intent(in) :: &
+ slipDamage
!*** output variables
@@ -2401,7 +2403,8 @@ end subroutine plastic_nonlocal_deltaState
!---------------------------------------------------------------------------------------------------
!> @brief calculates the rate of change of microstructure
!---------------------------------------------------------------------------------------------------
-subroutine plastic_nonlocal_dotState(Tstar_v, Fe, Fp, Temperature, timestep,subfrac, ip,el)
+subroutine plastic_nonlocal_dotState(Tstar_v, Fe, Fp, Temperature,nSlipDamage,slipDamage, &
+ timestep,subfrac, ip,el)
use prec, only: DAMASK_NaN
use numerics, only: numerics_integrationMode, &
@@ -2455,7 +2458,8 @@ implicit none
!*** input variables
integer(pInt), intent(in) :: ip, & !< current integration point
- el !< current element number
+ el, & !< current element number
+ nSlipDamage
real(pReal), intent(in) :: Temperature, & !< temperature
timestep !< substepped crystallite time increment
real(pReal), dimension(6), intent(in) :: Tstar_v !< current 2nd Piola-Kirchhoff stress in Mandel notation
@@ -2464,6 +2468,8 @@ real(pReal), dimension(homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems), in
real(pReal), dimension(3,3,homogenization_maxNgrains,mesh_maxNips,mesh_NcpElems), intent(in) :: &
Fe, & !< elastic deformation gradient
Fp !< plastic deformation gradient
+real(pReal), dimension(nSlipDamage), intent(in) :: &
+ slipDamage
!*** local variables
@@ -2628,7 +2634,7 @@ forall (t = 1_pInt:4_pInt) &
do s = 1_pInt,ns ! loop over slip systems
sLattice = slipSystemLattice(s,instance)
- tau(s) = math_mul6x6(Tstar_v, lattice_Sslip_v(1:6,1,sLattice,ph)) + tauBack(s)
+ tau(s) = math_mul6x6(Tstar_v, lattice_Sslip_v(1:6,1,sLattice,ph))/slipDamage(s) + tauBack(s)
if (abs(tau(s)) < 1.0e-15_pReal) tau(s) = 1.0e-15_pReal
enddo
diff --git a/code/plastic_phenopowerlaw.f90 b/code/plastic_phenopowerlaw.f90
index 6148db3a1..06b783141 100644
--- a/code/plastic_phenopowerlaw.f90
+++ b/code/plastic_phenopowerlaw.f90
@@ -552,7 +552,7 @@ subroutine plastic_phenopowerlaw_init(fileUnit)
! allocate state arrays
sizeState = plastic_phenopowerlaw_totalNslip(instance) & ! s_slip
+ plastic_phenopowerlaw_totalNtwin(instance) & ! s_twin
- + 2_pInt & ! sum(gamma) + sum(f)
+ + 2_pInt & ! sum(gamma) + sum(f)
+ plastic_phenopowerlaw_totalNslip(instance) & ! accshear_slip
+ plastic_phenopowerlaw_totalNtwin(instance) ! accshear_twin
@@ -710,7 +710,8 @@ end subroutine plastic_phenopowerlaw_aTolState
!--------------------------------------------------------------------------------------------------
!> @brief calculates plastic velocity gradient and its tangent
!--------------------------------------------------------------------------------------------------
-subroutine plastic_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,slipDamage,ipc,ip,el)
+subroutine plastic_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,nSlipDamage,slipDamage, &
+ ipc,ip,el)
use math, only: &
math_Plain3333to99, &
math_Mandel6to33
@@ -737,14 +738,13 @@ subroutine plastic_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,slipDam
dLp_dTstar99 !< derivative of Lp with respect to 2nd Piola Kirchhoff stress
integer(pInt), intent(in) :: &
+ nSlipDamage, &
ipc, & !< component-ID of integration point
ip, & !< integration point
el !< element
real(pReal), dimension(6), intent(in) :: &
Tstar_v !< 2nd Piola Kirchhoff stress tensor in Mandel notation
- real(pReal), &
- dimension(plastic_phenopowerlaw_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))), &
- intent(in) :: &
+ real(pReal), dimension(nSlipDamage), intent(in) :: &
slipDamage
integer(pInt) :: &
@@ -864,7 +864,7 @@ end subroutine plastic_phenopowerlaw_LpAndItsTangent
!--------------------------------------------------------------------------------------------------
!> @brief calculates the rate of change of microstructure
!--------------------------------------------------------------------------------------------------
-subroutine plastic_phenopowerlaw_dotState(Tstar_v,ipc,ip,el)
+subroutine plastic_phenopowerlaw_dotState(Tstar_v,nSlipDamage,slipDamage,ipc,ip,el)
use lattice, only: &
lattice_Sslip_v, &
lattice_Stwin_v, &
@@ -884,9 +884,12 @@ subroutine plastic_phenopowerlaw_dotState(Tstar_v,ipc,ip,el)
real(pReal), dimension(6), intent(in) :: &
Tstar_v !< 2nd Piola Kirchhoff stress tensor in Mandel notation
integer(pInt), intent(in) :: &
+ nSlipDamage, &
ipc, & !< component-ID of integration point
ip, & !< integration point
el !< element !< microstructure state
+ real(pReal), dimension(nSlipDamage), intent(in) :: &
+ slipDamage
integer(pInt) :: &
instance,ph, &
@@ -958,8 +961,8 @@ subroutine plastic_phenopowerlaw_dotState(Tstar_v,ipc,ip,el)
dot_product(Tstar_v,lattice_Sslip_v(1:6,2*k+1,index_myFamily+i,ph))
enddo nonSchmidSystems
gdot_slip(j) = plastic_phenopowerlaw_gdot0_slip(instance)*0.5_pReal* &
- ((abs(tau_slip_pos)/plasticState(ph)%state(j,of))**plastic_phenopowerlaw_n_slip(instance) &
- +(abs(tau_slip_neg)/plasticState(ph)%state(j,of))**plastic_phenopowerlaw_n_slip(instance))&
+ ((abs(tau_slip_pos)/(slipDamage(j)*plasticState(ph)%state(j,of)))**plastic_phenopowerlaw_n_slip(instance) &
+ +(abs(tau_slip_neg)/(slipDamage(j)*plasticState(ph)%state(j,of)))**plastic_phenopowerlaw_n_slip(instance))&
*sign(1.0_pReal,tau_slip_pos)
enddo slipSystems1
enddo slipFamilies1
diff --git a/code/plastic_titanmod.f90 b/code/plastic_titanmod.f90
index 59b95b771..13e685428 100644
--- a/code/plastic_titanmod.f90
+++ b/code/plastic_titanmod.f90
@@ -1326,7 +1326,8 @@ end subroutine plastic_titanmod_microstructure
!--------------------------------------------------------------------------------------------------
!> @brief calculates plastic velocity gradient and its tangent
!--------------------------------------------------------------------------------------------------
-subroutine plastic_titanmod_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,temperature,slipDamage,ipc,ip,el)
+subroutine plastic_titanmod_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,temperature,nSlipDamage,slipDamage, &
+ ipc,ip,el)
use math, only: &
math_Plain3333to99, &
math_Mandel6to33
@@ -1358,6 +1359,7 @@ subroutine plastic_titanmod_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,temperature,
dLp_dTstar99 !< derivative of Lp with respect to 2nd Piola Kirchhoff stress
integer(pInt), intent(in) :: &
+ nSlipDamage, &
ipc, & !< component-ID of integration point
ip, & !< integration point
el !< element
@@ -1365,9 +1367,7 @@ subroutine plastic_titanmod_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,temperature,
Tstar_v !< 2nd Piola Kirchhoff stress tensor in Mandel notation
real(pReal), intent(in) :: &
temperature !< temperature at IP
- real(pReal), &
- dimension(plastic_titanmod_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))), &
- intent(in) :: &
+ real(pReal), dimension(nSlipDamage), intent(in) :: &
slipDamage
integer(pInt) :: &
index_myFamily, instance, &