From 14d71eb35b47ffd18ceb16f4529a4ee9e1abd940 Mon Sep 17 00:00:00 2001 From: Pratheek Shanthraj Date: Tue, 28 Oct 2014 02:42:25 +0000 Subject: [PATCH] cleaner interface for ductile damage models --- code/constitutive.f90 | 114 ++++++++++++++++++++++------ code/constitutive_dislokmc.f90 | 10 ++- code/constitutive_dislotwin.f90 | 10 ++- code/constitutive_j2.f90 | 6 +- code/constitutive_nonlocal.f90 | 11 ++- code/constitutive_phenopowerlaw.f90 | 20 +++-- code/constitutive_titanmod.f90 | 18 +++-- code/damage_ductile.f90 | 18 +++++ code/damage_gurson.f90 | 20 +++++ 9 files changed, 175 insertions(+), 52 deletions(-) diff --git a/code/constitutive.f90 b/code/constitutive.f90 index e74a27f77..26cd5cffc 100644 --- a/code/constitutive.f90 +++ b/code/constitutive.f90 @@ -33,6 +33,7 @@ module constitutive constitutive_getLocalDamage, & constitutive_putLocalDamage, & constitutive_getDamage, & + constitutive_getSlipDamage, & constitutive_getDamageStrain, & constitutive_getDamageDiffusion33, & constitutive_getAdiabaticTemperature, & @@ -612,6 +613,7 @@ subroutine constitutive_LpAndItsTangent(Lp, dLp_dTstar, Tstar_v, ipc, ip, el) math_identity2nd use material, only: & phase_plasticity, & + phase_plasticityInstance, & material_phase, & plasticState,& mappingConstitutive, & @@ -625,15 +627,20 @@ subroutine constitutive_LpAndItsTangent(Lp, dLp_dTstar, Tstar_v, ipc, ip, el) use constitutive_j2, only: & constitutive_j2_LpAndItsTangent use constitutive_phenopowerlaw, only: & - constitutive_phenopowerlaw_LpAndItsTangent + constitutive_phenopowerlaw_LpAndItsTangent, & + constitutive_phenopowerlaw_totalNslip use constitutive_dislotwin, only: & - constitutive_dislotwin_LpAndItsTangent + constitutive_dislotwin_LpAndItsTangent, & + constitutive_dislotwin_totalNslip use constitutive_dislokmc, only: & - constitutive_dislokmc_LpAndItsTangent + constitutive_dislokmc_LpAndItsTangent, & + constitutive_dislokmc_totalNslip use constitutive_titanmod, only: & - constitutive_titanmod_LpAndItsTangent + constitutive_titanmod_LpAndItsTangent, & + constitutive_titanmod_totalNslip use constitutive_nonlocal, only: & - constitutive_nonlocal_LpAndItsTangent + constitutive_nonlocal_LpAndItsTangent, & + totalNslip implicit none integer(pInt), intent(in) :: & @@ -646,27 +653,48 @@ subroutine constitutive_LpAndItsTangent(Lp, dLp_dTstar, Tstar_v, ipc, ip, el) Lp !< plastic velocity gradient real(pReal), intent(out), dimension(9,9) :: & dLp_dTstar !< derivative of Lp with respect to Tstar (4th-order tensor) - real(pReal) :: damage, Tstar_v_effective(6) - - damage = constitutive_getDamage(ipc,ip,el) - Tstar_v_effective = Tstar_v/(damage*damage) + integer(pInt) :: & + nSlip + select case (phase_plasticity(material_phase(ipc,ip,el))) case (PLASTICITY_NONE_ID) Lp = 0.0_pReal dLp_dTstar = 0.0_pReal case (PLASTICITY_J2_ID) - call constitutive_j2_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v_effective,ipc,ip,el) + nSlip = 1_pInt + call constitutive_j2_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v, & + constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el), & + ipc,ip,el) case (PLASTICITY_PHENOPOWERLAW_ID) - call constitutive_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v_effective,ipc,ip,el) + nSlip = constitutive_phenopowerlaw_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el))) + call constitutive_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v, & + constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el), & + ipc,ip,el) case (PLASTICITY_NONLOCAL_ID) - call constitutive_nonlocal_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v_effective,constitutive_getTemperature(ipc,ip,el),ip,el) + nSlip = totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el))) + call constitutive_nonlocal_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v, & + constitutive_getTemperature(ipc,ip,el), & + constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el), & + ipc,ip,el) case (PLASTICITY_DISLOTWIN_ID) - call constitutive_dislotwin_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v_effective,constitutive_getTemperature(ipc,ip,el),ipc,ip,el) + nSlip = constitutive_dislotwin_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el))) + call constitutive_dislotwin_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v, & + constitutive_getTemperature(ipc,ip,el), & + constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el), & + ipc,ip,el) case (PLASTICITY_DISLOKMC_ID) - call constitutive_dislokmc_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v_effective,constitutive_getTemperature(ipc,ip,el),ipc,ip,el) + nSlip = constitutive_dislokmc_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el))) + call constitutive_dislokmc_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v, & + constitutive_getTemperature(ipc,ip,el), & + constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el), & + ipc,ip,el) case (PLASTICITY_TITANMOD_ID) - call constitutive_titanmod_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v_effective,constitutive_getTemperature(ipc,ip,el),ipc,ip,el) + nSlip = constitutive_titanmod_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el))) + call constitutive_titanmod_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v, & + constitutive_getTemperature(ipc,ip,el), & + constitutive_getSlipDamage(nSlip,Tstar_v,ipc,ip,el), & + ipc,ip,el) end select @@ -1093,6 +1121,45 @@ function constitutive_getDamage(ipc, ip, el) end function constitutive_getDamage +!-------------------------------------------------------------------------------------------------- +!> @brief Returns the damage on each slip system +!-------------------------------------------------------------------------------------------------- +function constitutive_getSlipDamage(nSlip, Tstar_v, ipc, ip, el) + use prec, only: & + pReal + use material, only: & + material_phase, & + LOCAL_DAMAGE_ductile_ID, & + LOCAL_DAMAGE_gurson_ID, & + phase_damage + use damage_ductile, only: & + damage_ductile_getSlipDamage + use damage_gurson, only: & + damage_gurson_getSlipDamage + + implicit none + integer(pInt), intent(in) :: & + nSlip, & + ipc, & !< grain number + ip, & !< integration point number + el !< element number + real(pReal), intent(in), dimension(6) :: & + Tstar_v !< 2nd Piola-Kirchhoff stress + real(pReal) :: constitutive_getSlipDamage(nSlip) + + select case (phase_damage(material_phase(ipc,ip,el))) + case (LOCAL_DAMAGE_ductile_ID) + constitutive_getSlipDamage = damage_ductile_getSlipDamage(ipc, ip, el) + + case (LOCAL_DAMAGE_gurson_ID) + constitutive_getSlipDamage = damage_gurson_getSlipDamage(Tstar_v, ipc, ip, el) + + case default + constitutive_getSlipDamage = 1.0_pReal + end select + +end function constitutive_getSlipDamage + !-------------------------------------------------------------------------------------------------- !> @brief returns damage deformation gradient !-------------------------------------------------------------------------------------------------- @@ -1644,12 +1711,9 @@ function constitutive_postResults(Tstar_v, FeArray, ipc, ip, el) FeArray !< elastic deformation gradient real(pReal), intent(in), dimension(6) :: & Tstar_v !< 2nd Piola Kirchhoff stress tensor (Mandel) - real(pReal) :: damage, Tstar_v_effective(6) - integer(pInt) :: startPos, endPos + integer(pInt) :: & + startPos, endPos - damage = constitutive_getDamage(ipc,ip,el) - Tstar_v_effective = damage*damage*Tstar_v - constitutive_postResults = 0.0_pReal startPos = 1_pInt @@ -1658,19 +1722,19 @@ function constitutive_postResults(Tstar_v, FeArray, ipc, ip, el) case (PLASTICITY_TITANMOD_ID) constitutive_postResults(startPos:endPos) = constitutive_titanmod_postResults(ipc,ip,el) case (PLASTICITY_J2_ID) - constitutive_postResults(startPos:endPos) = constitutive_j2_postResults(Tstar_v_effective,ipc,ip,el) + constitutive_postResults(startPos:endPos) = constitutive_j2_postResults(Tstar_v,ipc,ip,el) case (PLASTICITY_PHENOPOWERLAW_ID) constitutive_postResults(startPos:endPos) = & - constitutive_phenopowerlaw_postResults(Tstar_v_effective,ipc,ip,el) + constitutive_phenopowerlaw_postResults(Tstar_v,ipc,ip,el) case (PLASTICITY_DISLOTWIN_ID) constitutive_postResults(startPos:endPos) = & - constitutive_dislotwin_postResults(Tstar_v_effective,constitutive_getTemperature(ipc,ip,el),ipc,ip,el) + constitutive_dislotwin_postResults(Tstar_v,constitutive_getTemperature(ipc,ip,el),ipc,ip,el) case (PLASTICITY_DISLOKMC_ID) constitutive_postResults(startPos:endPos) = & - constitutive_dislokmc_postResults(Tstar_v_effective,constitutive_getTemperature(ipc,ip,el),ipc,ip,el) + constitutive_dislokmc_postResults(Tstar_v,constitutive_getTemperature(ipc,ip,el),ipc,ip,el) case (PLASTICITY_NONLOCAL_ID) constitutive_postResults(startPos:endPos) = & - constitutive_nonlocal_postResults (Tstar_v_effective,FeArray,ip,el) + constitutive_nonlocal_postResults (Tstar_v,FeArray,ip,el) end select #ifdef multiphysicsOut diff --git a/code/constitutive_dislokmc.f90 b/code/constitutive_dislokmc.f90 index e86d97f51..221c1a108 100644 --- a/code/constitutive_dislokmc.f90 +++ b/code/constitutive_dislokmc.f90 @@ -46,7 +46,7 @@ module constitutive_dislokmc integer(pInt), dimension(:), allocatable, target, public :: & constitutive_dislokmc_Noutput !< number of outputs per instance of this plasticity - integer(pInt), dimension(:), allocatable, private :: & + integer(pInt), dimension(:), allocatable, public, protected :: & constitutive_dislokmc_totalNslip, & !< total number of active slip systems for each instance constitutive_dislokmc_totalNtwin !< total number of active twin systems for each instance @@ -1193,7 +1193,7 @@ end subroutine constitutive_dislokmc_microstructure !-------------------------------------------------------------------------------------------------- !> @brief calculates plastic velocity gradient and its tangent !-------------------------------------------------------------------------------------------------- -subroutine constitutive_dislokmc_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v,Temperature,ipc,ip,el) +subroutine constitutive_dislokmc_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v,Temperature,slipDamage,ipc,ip,el) use prec, only: & tol_math_check use math, only: & @@ -1232,6 +1232,10 @@ subroutine constitutive_dislokmc_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v,Temperatu integer(pInt), intent(in) :: ipc,ip,el real(pReal), intent(in) :: Temperature real(pReal), dimension(6), intent(in) :: Tstar_v + real(pReal), & + dimension(constitutive_dislokmc_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))), & + intent(in) :: & + slipDamage real(pReal), dimension(3,3), intent(out) :: Lp real(pReal), dimension(9,9), intent(out) :: dLp_dTstar @@ -1273,7 +1277,7 @@ subroutine constitutive_dislokmc_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v,Temperatu !-------------------------------------------------------------------------------------------------- ! Calculation of Lp - tau_slip_pos(j) = dot_product(Tstar_v,lattice_Sslip_v(1:6,1,index_myFamily+i,ph)) + tau_slip_pos(j) = dot_product(Tstar_v,lattice_Sslip_v(1:6,1,index_myFamily+i,ph))/slipDamage(j) tau_slip_neg(j) = tau_slip_pos(j) nonSchmid_tensor(1:3,1:3,1) = lattice_Sslip(1:3,1:3,1,index_myFamily+i,ph) nonSchmid_tensor(1:3,1:3,2) = nonSchmid_tensor(1:3,1:3,1) diff --git a/code/constitutive_dislotwin.f90 b/code/constitutive_dislotwin.f90 index b0c60120a..1056e66ac 100644 --- a/code/constitutive_dislotwin.f90 +++ b/code/constitutive_dislotwin.f90 @@ -48,7 +48,7 @@ module constitutive_dislotwin integer(pInt), dimension(:), allocatable, target, public :: & constitutive_dislotwin_Noutput !< number of outputs per instance of this plasticity - integer(pInt), dimension(:), allocatable, private :: & + integer(pInt), dimension(:), allocatable, public, protected :: & constitutive_dislotwin_totalNslip, & !< total number of active slip systems for each instance constitutive_dislotwin_totalNtwin, & !< total number of active twin systems for each instance constitutive_dislotwin_totalNtrans !< number of active transformation systems @@ -1355,7 +1355,7 @@ end subroutine constitutive_dislotwin_microstructure !-------------------------------------------------------------------------------------------------- !> @brief calculates plastic velocity gradient and its tangent !-------------------------------------------------------------------------------------------------- -subroutine constitutive_dislotwin_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v,Temperature,ipc,ip,el) +subroutine constitutive_dislotwin_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v,Temperature,slipDamage,ipc,ip,el) use prec, only: & tol_math_check use math, only: & @@ -1398,6 +1398,10 @@ subroutine constitutive_dislotwin_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v,Temperat integer(pInt), intent(in) :: ipc,ip,el real(pReal), intent(in) :: Temperature real(pReal), dimension(6), intent(in) :: Tstar_v + real(pReal), & + dimension(constitutive_dislotwin_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))), & + intent(in) :: & + slipDamage real(pReal), dimension(3,3), intent(out) :: Lp real(pReal), dimension(9,9), intent(out) :: dLp_dTstar @@ -1463,7 +1467,7 @@ subroutine constitutive_dislotwin_LpAndItsTangent(Lp,dLp_dTstar,Tstar_v,Temperat !* Calculation of Lp !* 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/constitutive_j2.f90 b/code/constitutive_j2.f90 index 6aef32521..80acbf737 100644 --- a/code/constitutive_j2.f90 +++ b/code/constitutive_j2.f90 @@ -346,7 +346,7 @@ end subroutine constitutive_j2_init !-------------------------------------------------------------------------------------------------- !> @brief calculates plastic velocity gradient and its tangent !-------------------------------------------------------------------------------------------------- -subroutine constitutive_j2_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,ipc,ip,el) +subroutine constitutive_j2_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,slipDamage,ipc,ip,el) use math, only: & math_mul6x6, & math_Mandel6to33, & @@ -371,6 +371,8 @@ subroutine constitutive_j2_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,ipc,ip,el) 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) :: & ipc, & !< component-ID of integration point ip, & !< integration point @@ -398,7 +400,7 @@ subroutine constitutive_j2_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,ipc,ip,el) dLp_dTstar99 = 0.0_pReal else gamma_dot = constitutive_j2_gdot0(instance) & - * (sqrt(1.5_pReal) * norm_Tstar_dev / (constitutive_j2_fTaylor(instance) * & + * (sqrt(1.5_pReal) * norm_Tstar_dev / (slipDamage(1)*constitutive_j2_fTaylor(instance) * & plasticState(mappingConstitutive(2,ipc,ip,el))%state(1,mappingConstitutive(1,ipc,ip,el)))) & **constitutive_j2_n(instance) diff --git a/code/constitutive_nonlocal.f90 b/code/constitutive_nonlocal.f90 index ca0cf1cca..b48b48f75 100644 --- a/code/constitutive_nonlocal.f90 +++ b/code/constitutive_nonlocal.f90 @@ -69,7 +69,7 @@ module constitutive_nonlocal iV, & !< state indices for dislcation velocities iD !< state indices for stable dipole height - integer(pInt), dimension(:), allocatable, private :: & + integer(pInt), dimension(:), allocatable, public, protected :: & totalNslip !< total number of active slip systems for each instance integer(pInt), dimension(:,:), allocatable, private :: & @@ -2016,7 +2016,7 @@ end subroutine constitutive_nonlocal_kinetics !-------------------------------------------------------------------------------------------------- !> @brief calculates plastic velocity gradient and its tangent !-------------------------------------------------------------------------------------------------- -subroutine constitutive_nonlocal_LpAndItsTangent(Lp, dLp_dTstar99, Tstar_v, Temperature, ip, el) +subroutine constitutive_nonlocal_LpAndItsTangent(Lp, dLp_dTstar99, Tstar_v, Temperature, slipDamage, ipc, ip, el) use math, only: math_Plain3333to99, & math_mul6x6, & @@ -2042,10 +2042,13 @@ use mesh, only: mesh_ipVolume implicit none !*** input variables -integer(pInt), intent(in) :: ip, & !< current integration point +integer(pInt), intent(in) :: 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 !*** output variables @@ -2111,7 +2114,7 @@ tauThreshold = plasticState(ph)%state(iTauF(1:ns,instance),of) do s = 1_pInt,ns sLattice = slipSystemLattice(s,instance) - tau(s) = math_mul6x6(Tstar_v, lattice_Sslip_v(1:6,1,sLattice,ph)) + tau(s) = math_mul6x6(Tstar_v, lattice_Sslip_v(1:6,1,sLattice,ph))/slipDamage(s) tauNS(s,1) = tau(s) tauNS(s,2) = tau(s) if (tau(s) > 0.0_pReal) then diff --git a/code/constitutive_phenopowerlaw.f90 b/code/constitutive_phenopowerlaw.f90 index 27e9feb49..e2994cdc3 100644 --- a/code/constitutive_phenopowerlaw.f90 +++ b/code/constitutive_phenopowerlaw.f90 @@ -25,7 +25,7 @@ module constitutive_phenopowerlaw integer(pInt), dimension(:), allocatable, target, public :: & constitutive_phenopowerlaw_Noutput !< number of outputs per instance of this constitution - integer(pInt), dimension(:), allocatable, private :: & + integer(pInt), dimension(:), allocatable, public, protected :: & constitutive_phenopowerlaw_totalNslip, & !< no. of slip system used in simulation constitutive_phenopowerlaw_totalNtwin !< no. of twin system used in simulation @@ -669,7 +669,7 @@ end subroutine constitutive_phenopowerlaw_aTolState !-------------------------------------------------------------------------------------------------- !> @brief calculates plastic velocity gradient and its tangent !-------------------------------------------------------------------------------------------------- -subroutine constitutive_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,ipc,ip,el) +subroutine constitutive_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,slipDamage,ipc,ip,el) use math, only: & math_Plain3333to99, & math_Mandel6to33 @@ -699,12 +699,16 @@ subroutine constitutive_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,ip real(pReal), dimension(9,9), intent(out) :: & dLp_dTstar99 !< derivative of Lp with respect to 2nd Piola Kirchhoff stress - real(pReal), dimension(6), intent(in) :: & - Tstar_v !< 2nd Piola Kirchhoff stress tensor in Mandel notation integer(pInt), intent(in) :: & 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(constitutive_phenopowerlaw_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))), & + intent(in) :: & + slipDamage integer(pInt) :: & instance, & @@ -758,12 +762,12 @@ subroutine constitutive_phenopowerlaw_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,ip lattice_Sslip(1:3,1:3,2*k+1,index_myFamily+i,ph) enddo gdot_slip_pos(j) = 0.5_pReal*constitutive_phenopowerlaw_gdot0_slip(instance)* & - ((abs(tau_slip_pos(j))/plasticState(ph)%state(j,of))**constitutive_phenopowerlaw_n_slip(instance))*& - sign(1.0_pReal,tau_slip_pos(j)) + ((abs(tau_slip_pos(j))/(slipDamage(j)*plasticState(ph)%state(j,of)))** & + constitutive_phenopowerlaw_n_slip(instance))*sign(1.0_pReal,tau_slip_pos(j)) gdot_slip_neg(j) = 0.5_pReal*constitutive_phenopowerlaw_gdot0_slip(instance)* & - ((abs(tau_slip_neg(j))/plasticState(ph)%state(j,of))**constitutive_phenopowerlaw_n_slip(instance))*& - sign(1.0_pReal,tau_slip_neg(j)) + ((abs(tau_slip_neg(j))/(slipDamage(j)*plasticState(ph)%state(j,of)))**& + constitutive_phenopowerlaw_n_slip(instance))*sign(1.0_pReal,tau_slip_neg(j)) Lp = Lp + (1.0_pReal-plasticState(ph)%state(index_F,of))*& ! 1-F (gdot_slip_pos(j)+gdot_slip_neg(j))*lattice_Sslip(1:3,1:3,1,index_myFamily+i,ph) diff --git a/code/constitutive_titanmod.f90 b/code/constitutive_titanmod.f90 index d0dacccf6..c389008ff 100644 --- a/code/constitutive_titanmod.f90 +++ b/code/constitutive_titanmod.f90 @@ -45,7 +45,7 @@ module constitutive_titanmod integer(pInt), dimension(:), allocatable, target, public :: & constitutive_titanmod_Noutput !< number of outputs per instance of this plasticity !< ID of the lattice structure - integer(pInt), dimension(:), allocatable, private :: & + integer(pInt), dimension(:), allocatable, public, protected :: & constitutive_titanmod_totalNslip, & !< total number of active slip systems for each instance constitutive_titanmod_totalNtwin !< total number of active twin systems for each instance @@ -1326,7 +1326,7 @@ end subroutine constitutive_titanmod_microstructure !-------------------------------------------------------------------------------------------------- !> @brief calculates plastic velocity gradient and its tangent !-------------------------------------------------------------------------------------------------- -subroutine constitutive_titanmod_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,temperature,ipc,ip,el) +subroutine constitutive_titanmod_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,temperature,slipDamage,ipc,ip,el) use math, only: & math_Plain3333to99, & math_Mandel6to33 @@ -1357,14 +1357,18 @@ subroutine constitutive_titanmod_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,tempera real(pReal), dimension(9,9), intent(out) :: & dLp_dTstar99 !< derivative of Lp with respect to 2nd Piola Kirchhoff stress - real(pReal), dimension(6), intent(in) :: & - Tstar_v !< 2nd Piola Kirchhoff stress tensor in Mandel notation - real(pReal), intent(in) :: & - temperature !< temperature at IP integer(pInt), intent(in) :: & 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), intent(in) :: & + temperature !< temperature at IP + real(pReal), & + dimension(constitutive_titanmod_totalNslip(phase_plasticityInstance(material_phase(ipc,ip,el)))), & + intent(in) :: & + slipDamage integer(pInt) :: & index_myFamily, instance, & ns,nt, & @@ -1422,7 +1426,7 @@ subroutine constitutive_titanmod_LpAndItsTangent(Lp,dLp_dTstar99,Tstar_v,tempera !* Calculation of Lp !* 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(lattice_structure(ph)==LATTICE_hex_ID) then ! only for prismatic and pyr systems in hex screwvelocity_prefactor=constitutive_titanmod_debyefrequency(instance)* & plasticState(ph)%state(4_pInt*ns+nt+j, of)*(constitutive_titanmod_burgersPerSlipSys(j,instance)/ & diff --git a/code/damage_ductile.f90 b/code/damage_ductile.f90 index b057ce3ab..74c88aff6 100644 --- a/code/damage_ductile.f90 +++ b/code/damage_ductile.f90 @@ -45,6 +45,7 @@ module damage_ductile damage_ductile_dotState, & damage_ductile_microstructure, & damage_ductile_getDamage, & + damage_ductile_getSlipDamage, & damage_ductile_putLocalDamage, & damage_ductile_getLocalDamage, & damage_ductile_postResults @@ -352,6 +353,23 @@ function damage_ductile_getDamage(ipc, ip, el) end function damage_ductile_getDamage +!-------------------------------------------------------------------------------------------------- +!> @brief returns slip damage +!-------------------------------------------------------------------------------------------------- +function damage_ductile_getSlipDamage(ipc, ip, el) + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< grain number + ip, & !< integration point number + el !< element number + real(pReal) :: damage_ductile_getSlipDamage, damage + + damage = damage_ductile_getDamage(ipc, ip, el) + damage_ductile_getSlipDamage = damage*damage + +end function damage_ductile_getSlipDamage + !-------------------------------------------------------------------------------------------------- !> @brief puts local damage !-------------------------------------------------------------------------------------------------- diff --git a/code/damage_gurson.f90 b/code/damage_gurson.f90 index cc39e8347..0f849143f 100644 --- a/code/damage_gurson.f90 +++ b/code/damage_gurson.f90 @@ -51,6 +51,7 @@ module damage_gurson damage_gurson_dotState, & damage_gurson_microstructure, & damage_gurson_getDamage, & + damage_gurson_getSlipDamage, & damage_gurson_putLocalDamage, & damage_gurson_getLocalDamage, & damage_gurson_postResults @@ -410,6 +411,25 @@ function damage_gurson_getDamage(ipc, ip, el) end function damage_gurson_getDamage +!-------------------------------------------------------------------------------------------------- +!> @brief returns slip damage +!-------------------------------------------------------------------------------------------------- +function damage_gurson_getSlipDamage(Tstar_v, ipc, ip, el) + + implicit none + integer(pInt), intent(in) :: & + ipc, & !< grain number + ip, & !< integration point number + el !< element number + real(pReal), dimension(6), intent(in) :: & + Tstar_v !< 2nd Piola Kirchhoff stress tensor in Mandel notation + real(pReal) :: damage_gurson_getSlipDamage, porosity + + porosity = damage_gurson_getDamage(ipc, ip, el) + damage_gurson_getSlipDamage = porosity*porosity ! Gurson yield function should go here + +end function damage_gurson_getSlipDamage + !-------------------------------------------------------------------------------------------------- !> @brief puts local damage !--------------------------------------------------------------------------------------------------