223031012
/
DAMASK_EICMD
mirror of https://github.com/achalhp/DAMASK_EICMD.git
1
0
Fork 0
Code Issues Projects Releases Wiki Activity

common functions to be clubbed together

This commit is contained in:
Sharan Roongta 2020-07-10 17:10:23 +02:00
parent 957c51fb07
commit 80fb571fb4
5 changed files with 164 additions and 302 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

268
src/constitutive.f90
View File

@ -35,6 +35,103 @@ module constitutive
module subroutine thermal_init
end subroutine thermal_init
module subroutine plastic_isotropic_dotState(Mp,instance,of)
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
integer, intent(in) :: &
instance, &
of
end subroutine plastic_isotropic_dotState
module subroutine plastic_phenopowerlaw_dotState(Mp,instance,of)
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
integer, intent(in) :: &
instance, &
of
end subroutine plastic_phenopowerlaw_dotState
module subroutine plastic_kinehardening_dotState(Mp,instance,of)
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
integer, intent(in) :: &
instance, &
of
end subroutine plastic_kinehardening_dotState
module subroutine plastic_dislotwin_dotState(Mp,T,instance,of)
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
real(pReal), intent(in) :: &
T
integer, intent(in) :: &
instance, &
of
end subroutine plastic_dislotwin_dotState
module subroutine plastic_disloUCLA_dotState(Mp,T,instance,of)
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
real(pReal), intent(in) :: &
T
integer, intent(in) :: &
instance, &
of
end subroutine plastic_disloUCLA_dotState
module subroutine plastic_nonlocal_dotState(Mp, F, Fp, Temperature,timestep, &
instance,of,ip,el)
real(pReal), dimension(3,3), intent(in) ::&
Mp !< MandelStress
real(pReal), dimension(3,3,homogenization_maxNgrains,discretization_nIP,discretization_nElem), intent(in) :: &
F, & !< deformation gradient
Fp !< plastic deformation gradient
real(pReal), intent(in) :: &
Temperature, & !< temperature
timestep !< substepped crystallite time increment
integer, intent(in) :: &
instance, &
of, &
ip, & !< current integration point
el !< current element number
end subroutine plastic_nonlocal_dotState
module subroutine source_damage_anisoBrittle_dotState(S, ipc, ip, el)
integer, intent(in) :: &
ipc, & !< component-ID of integration point
ip, & !< integration point
el !< element
real(pReal), intent(in), dimension(3,3) :: &
S
end subroutine source_damage_anisoBrittle_dotState
module subroutine source_damage_anisoDuctile_dotState(ipc, ip, el)
integer, intent(in) :: &
ipc, & !< component-ID of integration point
ip, & !< integration point
el !< element
end subroutine source_damage_anisoDuctile_dotState
module subroutine source_damage_isoDuctile_dotState(ipc, ip, el)
integer, intent(in) :: &
ipc, & !< component-ID of integration point
ip, & !< integration point
el !< element
end subroutine source_damage_isoDuctile_dotState
module subroutine source_thermal_externalheat_dotState(phase, of)
integer, intent(in) :: &
phase, &
of
end subroutine source_thermal_externalheat_dotState
pure module function kinematics_thermal_expansion_initialStrain(homog,phase,offset) result(initialStrain)
integer, intent(in) :: &
@ -46,14 +143,14 @@ module constitutive
end function kinematics_thermal_expansion_initialStrain
module function plastic_homogenizedC(ipc,ip,el) result(homogenizedC)
module function constitutive_homogenizedC(ipc,ip,el) result(homogenizedC)
real(pReal), dimension(6,6) :: &
homogenizedC
integer, intent(in) :: &
ipc, & !< component-ID of integration point
ip, & !< integration point
el !< element
end function plastic_homogenizedC
end function constitutive_homogenizedC
module subroutine plastic_isotropic_LiAndItsTangent(Li,dLi_dMi,Mi,instance,of)
@ -109,65 +206,6 @@ module constitutive
dLi_dTstar !< derivative of Li with respect to Tstar (4th-order tensor defined to be zero)
end subroutine kinematics_thermal_expansion_LiAndItsTangent
module function plastic_dotState(S, FArray, Fi, FpArray, subdt, ipc, ip, el,phase,of) result(broken_plastic)
integer, intent(in) :: &
ipc, & !< component-ID of integration point
ip, & !< integration point
el, & !< element
phase, &
of
real(pReal), intent(in) :: &
subdt !< timestep
real(pReal), intent(in), dimension(3,3,homogenization_maxNgrains,discretization_nIP,discretization_nElem) :: &
FArray, & !< elastic deformation gradient
FpArray !< plastic deformation gradient
real(pReal), intent(in), dimension(3,3) :: &
Fi !< intermediate deformation gradient
real(pReal), intent(in), dimension(3,3) :: &
S !< 2nd Piola Kirchhoff stress (vector notation)
logical :: broken_plastic
end function plastic_dotState
module function damage_dotState(S, FArray, Fi, FpArray, subdt, ipc, ip, el,phase,of) result(broken_damage)
integer, intent(in) :: &
ipc, & !< component-ID of integration point
ip, & !< integration point
el, & !< element
phase, &
of
real(pReal), intent(in) :: &
subdt !< timestep
real(pReal), intent(in), dimension(3,3,homogenization_maxNgrains,discretization_nIP,discretization_nElem) :: &
FArray, & !< elastic deformation gradient
FpArray !< plastic deformation gradient
real(pReal), intent(in), dimension(3,3) :: &
Fi !< intermediate deformation gradient
real(pReal), intent(in), dimension(3,3) :: &
S !< 2nd Piola Kirchhoff stress (vector notation)
logical :: broken_damage
end function damage_dotState
module function thermal_dotState(S, FArray, Fi, FpArray, subdt, ipc, ip, el,phase,of) result(broken_thermal)
integer, intent(in) :: &
ipc, & !< component-ID of integration point
ip, & !< integration point
el, & !< element
phase, &
of
real(pReal), intent(in) :: &
subdt !< timestep
real(pReal), intent(in), dimension(3,3,homogenization_maxNgrains,discretization_nIP,discretization_nElem) :: &
FArray, & !< elastic deformation gradient
FpArray !< plastic deformation gradient
real(pReal), intent(in), dimension(3,3) :: &
Fi !< intermediate deformation gradient
real(pReal), intent(in), dimension(3,3) :: &
S !< 2nd Piola Kirchhoff stress (vector notation)
logical :: broken_thermal
end function thermal_dotState
module subroutine plastic_kinehardening_deltaState(Mp,instance,of)
real(pReal), dimension(3,3), intent(in) :: &
@ -261,7 +299,7 @@ module constitutive
Fp !< plastic deformation gradient
end subroutine plastic_dependentState
module subroutine plastic_LpAndItsTangents(Lp, dLp_dS, dLp_dFi, &
module subroutine constitutive_plastic_LpAndItsTangents(Lp, dLp_dS, dLp_dFi, &
S, Fi, ipc, ip, el)
integer, intent(in) :: &
ipc, & !< component-ID of integration point
@ -276,7 +314,7 @@ module constitutive
dLp_dS, &
dLp_dFi !< derivative of Lp with respect to Fi
end subroutine plastic_LpAndItsTangents
end subroutine constitutive_plastic_LpAndItsTangents
module subroutine plastic_nonlocal_updateCompatibility(orientation,instance,i,e)
integer, intent(in) :: &
@ -331,7 +369,7 @@ module constitutive
public :: &
constitutive_init, &
constitutive_homogenizedC, &
constitutive_LpAndItsTangents, &
constitutive_plastic_LpAndItsTangents, &
constitutive_dependentState, &
constitutive_LiAndItsTangents, &
constitutive_initialFi, &
@ -393,25 +431,6 @@ subroutine constitutive_init
end subroutine constitutive_init
subroutine constitutive_LpAndItsTangents(Lp, dLp_dS, dLp_dFi, &
S, Fi, ipc, ip, el)
integer, intent(in) :: &
ipc, & !< component-ID of integration point
ip, & !< integration point
el !< element
real(pReal), intent(in), dimension(3,3) :: &
S, & !< 2nd Piola-Kirchhoff stress
Fi !< intermediate deformation gradient
real(pReal), intent(out), dimension(3,3) :: &
Lp !< plastic velocity gradient
real(pReal), intent(out), dimension(3,3,3,3) :: &
dLp_dS, &
dLp_dFi !< derivative of Lp with respect to Fi
call plastic_LpAndItsTangents(Lp, dLp_dS, dLp_dFi, S, Fi, ipc, ip, el)
end subroutine constitutive_LpAndItsTangents
!--------------------------------------------------------------------------------------------------
!> @brief calls microstructure function of the different constitutive models
!--------------------------------------------------------------------------------------------------
@ -430,23 +449,6 @@ subroutine constitutive_dependentState(F, Fp, ipc, ip, el)
end subroutine constitutive_dependentState
!--------------------------------------------------------------------------------------------------
!> @brief returns the homogenize elasticity matrix
!> ToDo: homogenizedC66 would be more consistent
!--------------------------------------------------------------------------------------------------
function constitutive_homogenizedC(ipc,ip,el)
real(pReal), dimension(6,6) :: constitutive_homogenizedC
integer, intent(in) :: &
ipc, & !< component-ID of integration point
ip, & !< integration point
el !< element
constitutive_homogenizedC = plastic_homogenizedC(ipc,ip,el)
end function constitutive_homogenizedC
!--------------------------------------------------------------------------------------------------
!> @brief contains the constitutive equation for calculating the velocity gradient
! ToDo: MD: S is Mi?
@ -654,19 +656,65 @@ function constitutive_collectDotState(S, FArray, Fi, FpArray, subdt, ipc, ip, el
Fi !< intermediate deformation gradient
real(pReal), intent(in), dimension(3,3) :: &
S !< 2nd Piola Kirchhoff stress (vector notation)
real(pReal), dimension(3,3) :: &
Mp
integer :: &
ho, &
tme, &
i, &
instance
logical :: broken
logical :: &
broken_plastic, &
broken_thermal, &
broken_damage
ho = material_homogenizationAt(el)
tme = thermalMapping(ho)%p(ip,el)
instance = phase_plasticityInstance(phase)
broken_plastic = plastic_dotState(S, FArray, Fi, FpArray, subdt, ipc, ip, el,phase,of)
broken = broken_plastic
broken_damage = damage_dotState(S, FArray, Fi, FpArray, subdt, ipc, ip, el,phase,of)
broken = broken .or. broken_damage
broken_thermal = thermal_dotState(S, FArray, Fi, FpArray, subdt, ipc, ip, el,phase,of)
broken = broken .or. broken_thermal
Mp = matmul(matmul(transpose(Fi),Fi),S)
plasticityType: select case (phase_plasticity(phase))
case (PLASTICITY_ISOTROPIC_ID) plasticityType
call plastic_isotropic_dotState (Mp,instance,of)
case (PLASTICITY_PHENOPOWERLAW_ID) plasticityType
call plastic_phenopowerlaw_dotState(Mp,instance,of)
case (PLASTICITY_KINEHARDENING_ID) plasticityType
call plastic_kinehardening_dotState(Mp,instance,of)
case (PLASTICITY_DISLOTWIN_ID) plasticityType
call plastic_dislotwin_dotState (Mp,temperature(ho)%p(tme),instance,of)
case (PLASTICITY_DISLOUCLA_ID) plasticityType
call plastic_disloucla_dotState (Mp,temperature(ho)%p(tme),instance,of)
case (PLASTICITY_NONLOCAL_ID) plasticityType
call plastic_nonlocal_dotState (Mp,FArray,FpArray,temperature(ho)%p(tme),subdt, &
instance,of,ip,el)
end select plasticityType
broken = any(IEEE_is_NaN(plasticState(phase)%dotState(:,of)))
SourceLoop: do i = 1, phase_Nsources(phase)
sourceType: select case (phase_source(i,phase))
case (SOURCE_damage_anisoBrittle_ID) sourceType
call source_damage_anisoBrittle_dotState (S, ipc, ip, el) !< correct stress?
case (SOURCE_damage_isoDuctile_ID) sourceType
call source_damage_isoDuctile_dotState ( ipc, ip, el)
case (SOURCE_damage_anisoDuctile_ID) sourceType
call source_damage_anisoDuctile_dotState ( ipc, ip, el)
case (SOURCE_thermal_externalheat_ID) sourceType
call source_thermal_externalheat_dotState(phase,of)
end select sourceType
broken = broken .or. any(IEEE_is_NaN(sourceState(phase)%p(i)%dotState(:,of)))
enddo SourceLoop
end function constitutive_collectDotState

54
src/constitutive_damage.f90
View File

@ -22,32 +22,6 @@ submodule(constitutive) constitutive_damage
module subroutine kinematics_slipplane_opening_init
end subroutine kinematics_slipplane_opening_init
module subroutine source_damage_anisoBrittle_dotState(S, ipc, ip, el)
integer, intent(in) :: &
ipc, & !< component-ID of integration point
ip, & !< integration point
el !< element
real(pReal), intent(in), dimension(3,3) :: &
S
end subroutine source_damage_anisoBrittle_dotState
module subroutine source_damage_anisoDuctile_dotState(ipc, ip, el)
integer, intent(in) :: &
ipc, & !< component-ID of integration point
ip, & !< integration point
el !< element
end subroutine source_damage_anisoDuctile_dotState
module subroutine source_damage_isoDuctile_dotState(ipc, ip, el)
integer, intent(in) :: &
ipc, & !< component-ID of integration point
ip, & !< integration point
el !< element
end subroutine source_damage_isoDuctile_dotState
module subroutine source_damage_anisobrittle_getRateAndItsTangent(localphiDot, dLocalphiDot_dPhi, phi, phase, constituent)
@ -134,34 +108,6 @@ module subroutine damage_init
end subroutine damage_init
!--------------------------------------------------------------------------------------------------
!> @brief contains the constitutive equation for calculating the rate of change of microstructure
!--------------------------------------------------------------------------------------------------
module procedure damage_dotState
integer :: i
SourceLoop: do i = 1, phase_Nsources(phase)
sourceType: select case (phase_source(i,phase))
case (SOURCE_damage_anisoBrittle_ID) sourceType
call source_damage_anisoBrittle_dotState (S, ipc, ip, el) !< correct stress?
case (SOURCE_damage_isoDuctile_ID) sourceType
call source_damage_isoDuctile_dotState ( ipc, ip, el)
case (SOURCE_damage_anisoDuctile_ID) sourceType
call source_damage_anisoDuctile_dotState ( ipc, ip, el)
end select sourceType
enddo sourceLoop
broken_damage = any(IEEE_is_NaN(sourceState(phase)%p(i)%dotState(:,of)))
end procedure damage_dotState
module procedure damage_source_getRateAndItsTangents

112
src/constitutive_plastic.f90
View File

@ -25,67 +25,6 @@ submodule(constitutive) constitutive_plastic
module subroutine plastic_nonlocal_init
end subroutine plastic_nonlocal_init
module subroutine plastic_isotropic_dotState(Mp,instance,of)
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
integer, intent(in) :: &
instance, &
of
end subroutine plastic_isotropic_dotState
module subroutine plastic_phenopowerlaw_dotState(Mp,instance,of)
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
integer, intent(in) :: &
instance, &
of
end subroutine plastic_phenopowerlaw_dotState
module subroutine plastic_kinehardening_dotState(Mp,instance,of)
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
integer, intent(in) :: &
instance, &
of
end subroutine plastic_kinehardening_dotState
module subroutine plastic_dislotwin_dotState(Mp,T,instance,of)
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
real(pReal), intent(in) :: &
T
integer, intent(in) :: &
instance, &
of
end subroutine plastic_dislotwin_dotState
module subroutine plastic_disloUCLA_dotState(Mp,T,instance,of)
real(pReal), dimension(3,3), intent(in) :: &
Mp !< Mandel stress
real(pReal), intent(in) :: &
T
integer, intent(in) :: &
instance, &
of
end subroutine plastic_disloUCLA_dotState
module subroutine plastic_nonlocal_dotState(Mp, F, Fp, Temperature,timestep, &
instance,of,ip,el)
real(pReal), dimension(3,3), intent(in) ::&
Mp !< MandelStress
real(pReal), dimension(3,3,homogenization_maxNgrains,discretization_nIP,discretization_nElem), intent(in) :: &
F, & !< deformation gradient
Fp !< plastic deformation gradient
real(pReal), intent(in) :: &
Temperature, & !< temperature
timestep !< substepped crystallite time increment
integer, intent(in) :: &
instance, &
of, &
ip, & !< current integration point
el !< current element number
end subroutine plastic_nonlocal_dotState
module subroutine plastic_isotropic_LpAndItsTangent(Lp,dLp_dMp,Mp,instance,of)
real(pReal), dimension(3,3), intent(out) :: &
@ -236,55 +175,12 @@ module subroutine plastic_init
end subroutine plastic_init
!--------------------------------------------------------------------------------------------------
!> @brief contains the constitutive equation for calculating the rate of change of microstructure
!--------------------------------------------------------------------------------------------------
module procedure plastic_dotState
real(pReal), dimension(3,3) :: &
Mp
integer :: &
ho, & !< homogenization
tme, & !< thermal member position
instance
ho = material_homogenizationAt(el)
tme = thermalMapping(ho)%p(ip,el)
instance = phase_plasticityInstance(phase)
Mp = matmul(matmul(transpose(Fi),Fi),S)
plasticityType: select case (phase_plasticity(phase))
case (PLASTICITY_ISOTROPIC_ID) plasticityType
call plastic_isotropic_dotState (Mp,instance,of)
case (PLASTICITY_PHENOPOWERLAW_ID) plasticityType
call plastic_phenopowerlaw_dotState(Mp,instance,of)
case (PLASTICITY_KINEHARDENING_ID) plasticityType
call plastic_kinehardening_dotState(Mp,instance,of)
case (PLASTICITY_DISLOTWIN_ID) plasticityType
call plastic_dislotwin_dotState (Mp,temperature(ho)%p(tme),instance,of)
case (PLASTICITY_DISLOUCLA_ID) plasticityType
call plastic_disloucla_dotState (Mp,temperature(ho)%p(tme),instance,of)
case (PLASTICITY_NONLOCAL_ID) plasticityType
call plastic_nonlocal_dotState (Mp,FArray,FpArray,temperature(ho)%p(tme),subdt, &
instance,of,ip,el)
end select plasticityType
broken_plastic = any(IEEE_is_NaN(plasticState(phase)%dotState(:,of)))
end procedure plastic_dotState
!--------------------------------------------------------------------------------------------------
!> @brief returns the homogenize elasticity matrix
!> ToDo: homogenizedC66 would be more consistent
!--------------------------------------------------------------------------------------------------
module procedure plastic_homogenizedC
module procedure constitutive_homogenizedC
plasticityType: select case (phase_plasticity(material_phaseAt(ipc,el)))
case (PLASTICITY_DISLOTWIN_ID) plasticityType
@ -293,7 +189,7 @@ module procedure plastic_homogenizedC
homogenizedC = lattice_C66(1:6,1:6,material_phaseAt(ipc,el))
end select plasticityType
end procedure plastic_homogenizedC
end procedure constitutive_homogenizedC
!--------------------------------------------------------------------------------------------------
@ -327,7 +223,7 @@ end procedure plastic_dependentState
! ToDo: Discuss whether it makes sense if crystallite handles the configuration conversion, i.e.
! Mp in, dLp_dMp out
!--------------------------------------------------------------------------------------------------
module procedure plastic_LpAndItsTangents
module procedure constitutive_plastic_LpAndItsTangents
real(pReal), dimension(3,3,3,3) :: &
dLp_dMp !< derivative of Lp with respect to Mandel stress
@ -378,7 +274,7 @@ module procedure plastic_LpAndItsTangents
dLp_dS(i,j,1:3,1:3) = matmul(matmul(transpose(Fi),Fi),dLp_dMp(i,j,1:3,1:3)) ! ToDo: @PS: why not: dLp_dMp:(FiT Fi)
enddo; enddo
end procedure plastic_LpAndItsTangents
end procedure constitutive_plastic_LpAndItsTangents
end submodule constitutive_plastic

28
src/constitutive_thermal.f90
View File

@ -13,12 +13,6 @@ submodule(constitutive) constitutive_thermal
module subroutine kinematics_thermal_expansion_init
end subroutine kinematics_thermal_expansion_init
module subroutine source_thermal_externalheat_dotState(phase, of)
integer, intent(in) :: &
phase, &
of
end subroutine source_thermal_externalheat_dotState
module subroutine source_thermal_dissipation_getRateAndItsTangent(TDot, dTDot_dT, Tstar, Lp, phase)
@ -61,28 +55,6 @@ module subroutine thermal_init
end subroutine thermal_init
!--------------------------------------------------------------------------------------------------
!> @brief contains the constitutive equation for calculating the rate of change of microstructure
!--------------------------------------------------------------------------------------------------
module procedure thermal_dotState
integer :: i
SourceLoop: do i = 1, phase_Nsources(phase)
sourceType: select case (phase_source(i,phase))
case (SOURCE_thermal_externalheat_ID) sourceType
call source_thermal_externalheat_dotState(phase,of)
end select sourceType
enddo sourceLoop
broken_thermal = any(IEEE_is_NaN(sourceState(phase)%p(i)%dotState(:,of)))
end procedure thermal_dotState
module procedure thermal_source_getRateAndItsTangents

4
src/crystallite.f90
View File

@ -554,7 +554,7 @@ subroutine crystallite_stressTangent
dLidS = math_mul3333xx3333(dLidFi,dFidS) + dLidS
endif
call constitutive_LpAndItsTangents(devNull,dLpdS,dLpdFi, &
call constitutive_plastic_LpAndItsTangents(devNull,dLpdS,dLpdFi, &
crystallite_S (1:3,1:3,c,i,e), &
crystallite_Fi(1:3,1:3,c,i,e),c,i,e) ! call constitutive law to calculate Lp tangent in lattice configuration
dLpdS = math_mul3333xx3333(dLpdFi,dFidS) + dLpdS
@ -915,7 +915,7 @@ function integrateStress(ipc,ip,el,timeFraction) result(broken)
call constitutive_SandItsTangents(S, dS_dFe, dS_dFi, &
Fe, Fi_new, ipc, ip, el)
call constitutive_LpAndItsTangents(Lp_constitutive, dLp_dS, dLp_dFi, &
call constitutive_plastic_LpAndItsTangents(Lp_constitutive, dLp_dS, dLp_dFi, &
S, Fi_new, ipc, ip, el)
!* update current residuum and check for convergence of loop