simplified access pattern

This commit is contained in:
Martin Diehl 2021-02-13 11:01:08 +01:00
parent 775a51faa1
commit d9699b0f2e
5 changed files with 71 additions and 118 deletions

View File

@ -43,82 +43,65 @@ submodule(phase) damagee
end function isoductile_init end function isoductile_init
module subroutine source_damage_isoBrittle_deltaState(C, Fe, ph, me) module subroutine isobrittle_deltaState(C, Fe, ph, me)
integer, intent(in) :: ph,me integer, intent(in) :: ph,me
real(pReal), intent(in), dimension(3,3) :: & real(pReal), intent(in), dimension(3,3) :: &
Fe Fe
real(pReal), intent(in), dimension(6,6) :: & real(pReal), intent(in), dimension(6,6) :: &
C C
end subroutine source_damage_isoBrittle_deltaState end subroutine isobrittle_deltaState
module subroutine anisobrittle_dotState(S, co, ip, el) module subroutine anisobrittle_dotState(S, ph, me)
integer, intent(in) :: & integer, intent(in) :: ph,me
co, & !< component-ID of integration point
ip, & !< integration point
el !< element
real(pReal), intent(in), dimension(3,3) :: & real(pReal), intent(in), dimension(3,3) :: &
S S
end subroutine anisobrittle_dotState end subroutine anisobrittle_dotState
module subroutine anisoductile_dotState(co, ip, el) module subroutine anisoductile_dotState(ph,me)
integer, intent(in) :: & integer, intent(in) :: ph,me
co, & !< component-ID of integration point
ip, & !< integration point
el !< element
end subroutine anisoductile_dotState end subroutine anisoductile_dotState
module subroutine isoductile_dotState(co, ip, el) module subroutine isoductile_dotState(ph,me)
integer, intent(in) :: & integer, intent(in) :: ph,me
co, & !< component-ID of integration point
ip, & !< integration point
el !< element
end subroutine isoductile_dotState end subroutine isoductile_dotState
module subroutine source_damage_anisobrittle_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, phase, constituent) module subroutine anisobrittle_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, ph, me)
integer, intent(in) :: & integer, intent(in) :: ph,me
phase, & !< phase ID of element
constituent !< position of element within its phase instance
real(pReal), intent(in) :: & real(pReal), intent(in) :: &
phi !< damage parameter phi !< damage parameter
real(pReal), intent(out) :: & real(pReal), intent(out) :: &
localphiDot, & localphiDot, &
dLocalphiDot_dPhi dLocalphiDot_dPhi
end subroutine source_damage_anisoBrittle_getRateAndItsTangent end subroutine anisobrittle_getRateAndItsTangent
module subroutine source_damage_anisoDuctile_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, phase, constituent) module subroutine anisoductile_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, ph,me)
integer, intent(in) :: & integer, intent(in) :: ph,me
phase, & !< phase ID of element
constituent !< position of element within its phase instance
real(pReal), intent(in) :: & real(pReal), intent(in) :: &
phi !< damage parameter phi !< damage parameter
real(pReal), intent(out) :: & real(pReal), intent(out) :: &
localphiDot, & localphiDot, &
dLocalphiDot_dPhi dLocalphiDot_dPhi
end subroutine source_damage_anisoDuctile_getRateAndItsTangent end subroutine anisoductile_getRateAndItsTangent
module subroutine source_damage_isoBrittle_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, phase, constituent) module subroutine isobrittle_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, ph,me)
integer, intent(in) :: & integer, intent(in) :: ph,me
phase, & !< phase ID of element
constituent !< position of element within its phase instance
real(pReal), intent(in) :: & real(pReal), intent(in) :: &
phi !< damage parameter phi !< damage parameter
real(pReal), intent(out) :: & real(pReal), intent(out) :: &
localphiDot, & localphiDot, &
dLocalphiDot_dPhi dLocalphiDot_dPhi
end subroutine source_damage_isoBrittle_getRateAndItsTangent end subroutine isobrittle_getRateAndItsTangent
module subroutine source_damage_isoDuctile_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, phase, constituent) module subroutine isoductile_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, ph,me)
integer, intent(in) :: & integer, intent(in) :: ph,me
phase, & !< phase ID of element
constituent !< position of element within its phase instance
real(pReal), intent(in) :: & real(pReal), intent(in) :: &
phi !< damage parameter phi !< damage parameter
real(pReal), intent(out) :: & real(pReal), intent(out) :: &
localphiDot, & localphiDot, &
dLocalphiDot_dPhi dLocalphiDot_dPhi
end subroutine source_damage_isoDuctile_getRateAndItsTangent end subroutine isoductile_getRateAndItsTangent
module subroutine anisobrittle_results(phase,group) module subroutine anisobrittle_results(phase,group)
integer, intent(in) :: phase integer, intent(in) :: phase
@ -225,16 +208,16 @@ module subroutine phase_damage_getRateAndItsTangents(phiDot, dPhiDot_dPhi, phi,
do so = 1, phase_Nsources(ph) do so = 1, phase_Nsources(ph)
select case(phase_source(so,ph)) select case(phase_source(so,ph))
case (DAMAGE_ISOBRITTLE_ID) case (DAMAGE_ISOBRITTLE_ID)
call source_damage_isobrittle_getRateAndItsTangent (localphiDot, dLocalphiDot_dPhi, phi, ph, me) call isobrittle_getRateAndItsTangent (localphiDot, dLocalphiDot_dPhi, phi, ph, me)
case (DAMAGE_ISODUCTILE_ID) case (DAMAGE_ISODUCTILE_ID)
call source_damage_isoductile_getRateAndItsTangent (localphiDot, dLocalphiDot_dPhi, phi, ph, me) call isoductile_getRateAndItsTangent (localphiDot, dLocalphiDot_dPhi, phi, ph, me)
case (DAMAGE_ANISOBRITTLE_ID) case (DAMAGE_ANISOBRITTLE_ID)
call source_damage_anisobrittle_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, ph, me) call anisobrittle_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, ph, me)
case (DAMAGE_ANISODUCTILE_ID) case (DAMAGE_ANISODUCTILE_ID)
call source_damage_anisoductile_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, ph, me) call anisoductile_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, ph, me)
case default case default
localphiDot = 0.0_pReal localphiDot = 0.0_pReal
@ -281,7 +264,7 @@ module function integrateDamageState(dt,co,ip,el) result(broken)
me = material_phaseMemberAt(co,ip,el) me = material_phaseMemberAt(co,ip,el)
converged_ = .true. converged_ = .true.
broken = phase_damage_collectDotState(co,ip,el,ph,me) broken = phase_damage_collectDotState(ph,me)
if(broken) return if(broken) return
size_so = damageState(ph)%sizeDotState size_so = damageState(ph)%sizeDotState
@ -294,7 +277,7 @@ module function integrateDamageState(dt,co,ip,el) result(broken)
if(nIterationState > 1) source_dotState(1:size_so,2) = source_dotState(1:size_so,1) if(nIterationState > 1) source_dotState(1:size_so,2) = source_dotState(1:size_so,1)
source_dotState(1:size_so,1) = damageState(ph)%dotState(:,me) source_dotState(1:size_so,1) = damageState(ph)%dotState(:,me)
broken = phase_damage_collectDotState(co,ip,el,ph,me) broken = phase_damage_collectDotState(ph,me)
if(broken) exit iteration if(broken) exit iteration
@ -384,39 +367,34 @@ end subroutine damage_results
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief contains the constitutive equation for calculating the rate of change of microstructure !> @brief contains the constitutive equation for calculating the rate of change of microstructure
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
function phase_damage_collectDotState(co,ip,el,ph,me) result(broken) function phase_damage_collectDotState(ph,me) result(broken)
integer, intent(in) :: & integer, intent(in) :: &
co, & !< component-ID me integration point
ip, & !< integration point
el, & !< element
ph, & ph, &
me me !< counter in source loop
integer :: &
so !< counter in source loop
logical :: broken logical :: broken
broken = .false. broken = .false.
SourceLoop: do so = 1, phase_Nsources(ph) if (phase_Nsources(ph)==1) then
sourceType: select case (phase_source(so,ph)) sourceType: select case (phase_source(1,ph))
case (DAMAGE_ISODUCTILE_ID) sourceType case (DAMAGE_ISODUCTILE_ID) sourceType
call isoductile_dotState(co, ip, el) call isoductile_dotState(ph,me)
case (DAMAGE_ANISODUCTILE_ID) sourceType case (DAMAGE_ANISODUCTILE_ID) sourceType
call anisoductile_dotState(co, ip, el) call anisoductile_dotState(ph,me)
case (DAMAGE_ANISOBRITTLE_ID) sourceType case (DAMAGE_ANISOBRITTLE_ID) sourceType
call anisobrittle_dotState(mechanical_S(ph,me),co, ip, el) ! correct stress? call anisobrittle_dotState(mechanical_S(ph,me), ph,me) ! correct stress?
end select sourceType end select sourceType
broken = broken .or. any(IEEE_is_NaN(damageState(ph)%dotState(:,me))) broken = broken .or. any(IEEE_is_NaN(damageState(ph)%dotState(:,me)))
enddo SourceLoop endif
end function phase_damage_collectDotState end function phase_damage_collectDotState
@ -447,7 +425,7 @@ function phase_damage_deltaState(Fe, ph, me) result(broken)
sourceType: select case (phase_source(1,ph)) sourceType: select case (phase_source(1,ph))
case (DAMAGE_ISOBRITTLE_ID) sourceType case (DAMAGE_ISOBRITTLE_ID) sourceType
call source_damage_isoBrittle_deltaState(phase_homogenizedC(ph,me), Fe, ph,me) call isobrittle_deltaState(phase_homogenizedC(ph,me), Fe, ph,me)
broken = any(IEEE_is_NaN(damageState(ph)%deltaState(:,me))) broken = any(IEEE_is_NaN(damageState(ph)%deltaState(:,me)))
if(.not. broken) then if(.not. broken) then
myOffset = damageState(ph)%offsetDeltaState myOffset = damageState(ph)%offsetDeltaState

View File

@ -113,18 +113,14 @@ end function anisobrittle_init
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief calculates derived quantities from state !> @brief calculates derived quantities from state
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
module subroutine anisobrittle_dotState(S, co, ip, el) module subroutine anisobrittle_dotState(S, ph,me)
integer, intent(in) :: & integer, intent(in) :: &
co, & !< component-ID of integration point ph,me
ip, & !< integration point
el !< element
real(pReal), intent(in), dimension(3,3) :: & real(pReal), intent(in), dimension(3,3) :: &
S S
integer :: & integer :: &
ph, &
me, &
sourceOffset, & sourceOffset, &
damageOffset, & damageOffset, &
homog, & homog, &
@ -133,10 +129,6 @@ module subroutine anisobrittle_dotState(S, co, ip, el)
traction_d, traction_t, traction_n, traction_crit traction_d, traction_t, traction_n, traction_crit
ph = material_phaseAt(co,el)
me = material_phasememberAt(co,ip,el)
associate(prm => param(ph)) associate(prm => param(ph))
damageState(ph)%dotState(1,me) = 0.0_pReal damageState(ph)%dotState(1,me) = 0.0_pReal
do i = 1, prm%sum_N_cl do i = 1, prm%sum_N_cl
@ -144,7 +136,7 @@ module subroutine anisobrittle_dotState(S, co, ip, el)
traction_t = math_tensordot(S,prm%cleavage_systems(1:3,1:3,2,i)) traction_t = math_tensordot(S,prm%cleavage_systems(1:3,1:3,2,i))
traction_n = math_tensordot(S,prm%cleavage_systems(1:3,1:3,3,i)) traction_n = math_tensordot(S,prm%cleavage_systems(1:3,1:3,3,i))
traction_crit = prm%g_crit(i)*phase_damage_get_phi(co,ip,el)**2.0_pReal traction_crit = prm%g_crit(i)*damage_phi(ph,me)**2.0_pReal
damageState(ph)%dotState(1,me) = damageState(ph)%dotState(1,me) & damageState(ph)%dotState(1,me) = damageState(ph)%dotState(1,me) &
+ prm%dot_o / prm%s_crit(i) & + prm%dot_o / prm%s_crit(i) &
@ -160,11 +152,11 @@ end subroutine anisobrittle_dotState
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief returns local part of nonlocal damage driving force !> @brief returns local part of nonlocal damage driving force
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
module subroutine source_damage_anisobrittle_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, phase, constituent) module subroutine anisobrittle_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, ph, me)
integer, intent(in) :: & integer, intent(in) :: &
phase, & ph, &
constituent me
real(pReal), intent(in) :: & real(pReal), intent(in) :: &
phi phi
real(pReal), intent(out) :: & real(pReal), intent(out) :: &
@ -172,12 +164,12 @@ module subroutine source_damage_anisobrittle_getRateAndItsTangent(localphiDot, d
dLocalphiDot_dPhi dLocalphiDot_dPhi
dLocalphiDot_dPhi = -damageState(phase)%state(1,constituent) dLocalphiDot_dPhi = -damageState(ph)%state(1,me)
localphiDot = 1.0_pReal & localphiDot = 1.0_pReal &
+ dLocalphiDot_dPhi*phi + dLocalphiDot_dPhi*phi
end subroutine source_damage_anisoBrittle_getRateAndItsTangent end subroutine anisobrittle_getRateAndItsTangent
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
@ -190,6 +182,7 @@ module subroutine anisobrittle_results(phase,group)
integer :: o integer :: o
associate(prm => param(phase), stt => damageState(phase)%state) associate(prm => param(phase), stt => damageState(phase)%state)
outputsLoop: do o = 1,size(prm%output) outputsLoop: do o = 1,size(prm%output)
select case(trim(prm%output(o))) select case(trim(prm%output(o)))

View File

@ -99,20 +99,12 @@ end function anisoductile_init
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief calculates derived quantities from state !> @brief calculates derived quantities from state
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
module subroutine anisoductile_dotState(co, ip, el) module subroutine anisoductile_dotState(ph,me)
integer, intent(in) :: & integer, intent(in) :: &
co, & !< component-ID of integration point
ip, & !< integration point
el !< element
integer :: &
ph, & ph, &
me me
ph = material_phaseAt(co,el)
me = material_phasememberAt(co,ip,el)
associate(prm => param(ph)) associate(prm => param(ph))
damageState(ph)%dotState(1,me) = sum(plasticState(ph)%slipRate(:,me)/(damage_phi(ph,me)**prm%q)/prm%gamma_crit) damageState(ph)%dotState(1,me) = sum(plasticState(ph)%slipRate(:,me)/(damage_phi(ph,me)**prm%q)/prm%gamma_crit)
@ -124,11 +116,11 @@ end subroutine anisoductile_dotState
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief returns local part of nonlocal damage driving force !> @brief returns local part of nonlocal damage driving force
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
module subroutine source_damage_anisoDuctile_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, phase, constituent) module subroutine anisoductile_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, ph,me)
integer, intent(in) :: & integer, intent(in) :: &
phase, & ph, &
constituent me
real(pReal), intent(in) :: & real(pReal), intent(in) :: &
phi phi
real(pReal), intent(out) :: & real(pReal), intent(out) :: &
@ -136,12 +128,12 @@ module subroutine source_damage_anisoDuctile_getRateAndItsTangent(localphiDot, d
dLocalphiDot_dPhi dLocalphiDot_dPhi
dLocalphiDot_dPhi = -damageState(phase)%state(1,constituent) dLocalphiDot_dPhi = -damageState(ph)%state(1,me)
localphiDot = 1.0_pReal & localphiDot = 1.0_pReal &
+ dLocalphiDot_dPhi*phi + dLocalphiDot_dPhi*phi
end subroutine source_damage_anisoDuctile_getRateAndItsTangent end subroutine anisoductile_getRateAndItsTangent
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
@ -154,8 +146,8 @@ module subroutine anisoductile_results(phase,group)
integer :: o integer :: o
associate(prm => param(phase), &
stt => damageState(phase)%state) associate(prm => param(phase), stt => damageState(phase)%state)
outputsLoop: do o = 1,size(prm%output) outputsLoop: do o = 1,size(prm%output)
select case(trim(prm%output(o))) select case(trim(prm%output(o)))
case ('f_phi') case ('f_phi')

View File

@ -86,7 +86,7 @@ end function isobrittle_init
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief calculates derived quantities from state !> @brief calculates derived quantities from state
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
module subroutine source_damage_isoBrittle_deltaState(C, Fe, ph,me) module subroutine isobrittle_deltaState(C, Fe, ph,me)
integer, intent(in) :: ph,me integer, intent(in) :: ph,me
real(pReal), intent(in), dimension(3,3) :: & real(pReal), intent(in), dimension(3,3) :: &
@ -106,24 +106,21 @@ module subroutine source_damage_isoBrittle_deltaState(C, Fe, ph,me)
strainenergy = 2.0_pReal*sum(strain*matmul(C,strain))/prm%W_crit strainenergy = 2.0_pReal*sum(strain*matmul(C,strain))/prm%W_crit
! ToDo: check strainenergy = 2.0_pReal*dot_product(strain,matmul(C,strain))/prm%W_crit ! ToDo: check strainenergy = 2.0_pReal*dot_product(strain,matmul(C,strain))/prm%W_crit
if (strainenergy > damageState(ph)%subState0(1,me)) then damageState(ph)%deltaState(1,me) = merge(strainenergy - damageState(ph)%state(1,me), &
damageState(ph)%deltaState(1,me) = strainenergy - damageState(ph)%state(1,me) damageState(ph)%subState0(1,me) - damageState(ph)%state(1,me), &
else strainenergy > damageState(ph)%subState0(1,me))
damageState(ph)%deltaState(1,me) = damageState(ph)%subState0(1,me) - damageState(ph)%state(1,me)
endif
end associate end associate
end subroutine source_damage_isoBrittle_deltaState end subroutine isobrittle_deltaState
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief returns local part of nonlocal damage driving force !> @brief returns local part of nonlocal damage driving force
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
module subroutine source_damage_isoBrittle_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, phase, constituent) module subroutine isobrittle_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, ph, me)
integer, intent(in) :: & integer, intent(in) :: &
phase, & ph, me
constituent
real(pReal), intent(in) :: & real(pReal), intent(in) :: &
phi phi
real(pReal), intent(out) :: & real(pReal), intent(out) :: &
@ -131,13 +128,13 @@ module subroutine source_damage_isoBrittle_getRateAndItsTangent(localphiDot, dLo
dLocalphiDot_dPhi dLocalphiDot_dPhi
associate(prm => param(phase)) associate(prm => param(ph))
localphiDot = 1.0_pReal & localphiDot = 1.0_pReal &
- phi*damageState(phase)%state(1,constituent) - phi*damageState(ph)%state(1,me)
dLocalphiDot_dPhi = - damageState(phase)%state(1,constituent) dLocalphiDot_dPhi = - damageState(ph)%state(1,me)
end associate end associate
end subroutine source_damage_isoBrittle_getRateAndItsTangent end subroutine isobrittle_getRateAndItsTangent
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------

View File

@ -90,23 +90,16 @@ end function isoductile_init
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief calculates derived quantities from state !> @brief calculates derived quantities from state
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
module subroutine isoductile_dotState(co, ip, el) module subroutine isoductile_dotState(ph, me)
integer, intent(in) :: & integer, intent(in) :: &
co, & !< component-ID of integration point
ip, & !< integration point
el !< element
integer :: &
ph, & ph, &
me me
ph = material_phaseAt(co,el)
me = material_phasememberAt(co,ip,el)
associate(prm => param(ph)) associate(prm => param(ph))
damageState(ph)%dotState(1,me) = sum(plasticState(ph)%slipRate(:,me))/(damage_phi(ph,me)**prm%q)/prm%gamma_crit damageState(ph)%dotState(1,me) = sum(plasticState(ph)%slipRate(:,me)) &
/ (prm%gamma_crit*damage_phi(ph,me)**prm%q)
end associate end associate
end subroutine isoductile_dotState end subroutine isoductile_dotState
@ -115,11 +108,11 @@ end subroutine isoductile_dotState
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
!> @brief returns local part of nonlocal damage driving force !> @brief returns local part of nonlocal damage driving force
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
module subroutine source_damage_isoDuctile_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, phase, constituent) module subroutine isoductile_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, ph, me)
integer, intent(in) :: & integer, intent(in) :: &
phase, & ph, &
constituent me
real(pReal), intent(in) :: & real(pReal), intent(in) :: &
phi phi
real(pReal), intent(out) :: & real(pReal), intent(out) :: &
@ -127,12 +120,12 @@ module subroutine source_damage_isoDuctile_getRateAndItsTangent(localphiDot, dLo
dLocalphiDot_dPhi dLocalphiDot_dPhi
dLocalphiDot_dPhi = -damageState(phase)%state(1,constituent) dLocalphiDot_dPhi = -damageState(ph)%state(1,me)
localphiDot = 1.0_pReal & localphiDot = 1.0_pReal &
+ dLocalphiDot_dPhi*phi + dLocalphiDot_dPhi*phi
end subroutine source_damage_isoDuctile_getRateAndItsTangent end subroutine isoductile_getRateAndItsTangent
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------