Merge branch 'MiscImprovements' of magit1.mpie.de:/damask/DAMASK into MiscImprovements

This commit is contained in:
Martin Diehl 2019-06-14 17:24:44 +02:00
commit 72298205a7
26 changed files with 210 additions and 224 deletions

View File

@ -262,7 +262,7 @@ subroutine CPFEM_general(mode, parallelExecution, ffn, ffn1, temperature_inp, dt
if (debug_e <= discretization_nElem .and. debug_i <=discretization_nIP) then if (debug_e <= discretization_nElem .and. debug_i <=discretization_nIP) then
write(6,'(a,1x,i8,1x,i2,1x,i4,/,(12x,6(e20.8,1x)),/)') & write(6,'(a,1x,i8,1x,i2,1x,i4,/,(12x,6(e20.8,1x)),/)') &
'<< CPFEM >> aged state of elFE ip grain',debug_e, debug_i, 1, & '<< CPFEM >> aged state of elFE ip grain',debug_e, debug_i, 1, &
plasticState(phaseAt(1,debug_i,debug_e))%state(:,phasememberAt(1,debug_i,debug_e)) plasticState(material_phaseAt(1,debug_e))%state(:,material_phasememberAt(1,debug_i,debug_e))
endif endif
endif endif

View File

@ -89,7 +89,6 @@ subroutine CPFEM_init
fileHandle = HDF5_openFile(trim(getSolverJobName())//trim(rankStr)//'.hdf5') fileHandle = HDF5_openFile(trim(getSolverJobName())//trim(rankStr)//'.hdf5')
call HDF5_read(fileHandle,material_phase, 'recordedPhase')
call HDF5_read(fileHandle,crystallite_F0, 'convergedF') call HDF5_read(fileHandle,crystallite_F0, 'convergedF')
call HDF5_read(fileHandle,crystallite_Fp0, 'convergedFp') call HDF5_read(fileHandle,crystallite_Fp0, 'convergedFp')
call HDF5_read(fileHandle,crystallite_Fi0, 'convergedFi') call HDF5_read(fileHandle,crystallite_Fi0, 'convergedFi')
@ -158,7 +157,6 @@ subroutine CPFEM_age
write(rankStr,'(a1,i0)')'_',worldrank write(rankStr,'(a1,i0)')'_',worldrank
fileHandle = HDF5_openFile(trim(getSolverJobName())//trim(rankStr)//'.hdf5','a') fileHandle = HDF5_openFile(trim(getSolverJobName())//trim(rankStr)//'.hdf5','a')
call HDF5_write(fileHandle,material_phase, 'recordedPhase')
call HDF5_write(fileHandle,crystallite_F0, 'convergedF') call HDF5_write(fileHandle,crystallite_F0, 'convergedF')
call HDF5_write(fileHandle,crystallite_Fp0, 'convergedFp') call HDF5_write(fileHandle,crystallite_Fp0, 'convergedFp')
call HDF5_write(fileHandle,crystallite_Fi0, 'convergedFi') call HDF5_write(fileHandle,crystallite_Fi0, 'convergedFi')

View File

@ -251,7 +251,7 @@ function constitutive_homogenizedC(ipc,ip,el)
ip, & !< integration point ip, & !< integration point
el !< element el !< element
plasticityType: select case (phase_plasticity(material_phase(ipc,ip,el))) plasticityType: select case (phase_plasticity(material_phaseAt(ipc,el)))
case (PLASTICITY_DISLOTWIN_ID) plasticityType case (PLASTICITY_DISLOTWIN_ID) plasticityType
constitutive_homogenizedC = plastic_dislotwin_homogenizedC(ipc,ip,el) constitutive_homogenizedC = plastic_dislotwin_homogenizedC(ipc,ip,el)
case default plasticityType case default plasticityType
@ -280,14 +280,14 @@ subroutine constitutive_microstructure(Fe, Fp, ipc, ip, el)
ho = material_homogenizationAt(el) ho = material_homogenizationAt(el)
tme = thermalMapping(ho)%p(ip,el) tme = thermalMapping(ho)%p(ip,el)
plasticityType: select case (phase_plasticity(material_phase(ipc,ip,el))) plasticityType: select case (phase_plasticity(material_phaseAt(ipc,el)))
case (PLASTICITY_DISLOTWIN_ID) plasticityType case (PLASTICITY_DISLOTWIN_ID) plasticityType
of = phasememberAt(ipc,ip,el) of = material_phasememberAt(ipc,ip,el)
instance = phase_plasticityInstance(material_phase(ipc,ip,el)) instance = phase_plasticityInstance(material_phaseAt(ipc,el))
call plastic_dislotwin_dependentState(temperature(ho)%p(tme),instance,of) call plastic_dislotwin_dependentState(temperature(ho)%p(tme),instance,of)
case (PLASTICITY_DISLOUCLA_ID) plasticityType case (PLASTICITY_DISLOUCLA_ID) plasticityType
of = phasememberAt(ipc,ip,el) of = material_phasememberAt(ipc,ip,el)
instance = phase_plasticityInstance(material_phase(ipc,ip,el)) instance = phase_plasticityInstance(material_phaseAt(ipc,el))
call plastic_disloUCLA_dependentState(instance,of) call plastic_disloUCLA_dependentState(instance,of)
case (PLASTICITY_NONLOCAL_ID) plasticityType case (PLASTICITY_NONLOCAL_ID) plasticityType
call plastic_nonlocal_dependentState (Fe,Fp,ip,el) call plastic_nonlocal_dependentState (Fe,Fp,ip,el)
@ -331,25 +331,25 @@ subroutine constitutive_LpAndItsTangents(Lp, dLp_dS, dLp_dFi, &
Mp = matmul(matmul(transpose(Fi),Fi),S) Mp = matmul(matmul(transpose(Fi),Fi),S)
plasticityType: select case (phase_plasticity(material_phase(ipc,ip,el))) plasticityType: select case (phase_plasticity(material_phaseAt(ipc,el)))
case (PLASTICITY_NONE_ID) plasticityType case (PLASTICITY_NONE_ID) plasticityType
Lp = 0.0_pReal Lp = 0.0_pReal
dLp_dMp = 0.0_pReal dLp_dMp = 0.0_pReal
case (PLASTICITY_ISOTROPIC_ID) plasticityType case (PLASTICITY_ISOTROPIC_ID) plasticityType
of = phasememberAt(ipc,ip,el) of = material_phasememberAt(ipc,ip,el)
instance = phase_plasticityInstance(material_phase(ipc,ip,el)) instance = phase_plasticityInstance(material_phaseAt(ipc,el))
call plastic_isotropic_LpAndItsTangent (Lp,dLp_dMp,Mp,instance,of) call plastic_isotropic_LpAndItsTangent (Lp,dLp_dMp,Mp,instance,of)
case (PLASTICITY_PHENOPOWERLAW_ID) plasticityType case (PLASTICITY_PHENOPOWERLAW_ID) plasticityType
of = phasememberAt(ipc,ip,el) of = material_phasememberAt(ipc,ip,el)
instance = phase_plasticityInstance(material_phase(ipc,ip,el)) instance = phase_plasticityInstance(material_phaseAt(ipc,el))
call plastic_phenopowerlaw_LpAndItsTangent (Lp,dLp_dMp,Mp,instance,of) call plastic_phenopowerlaw_LpAndItsTangent (Lp,dLp_dMp,Mp,instance,of)
case (PLASTICITY_KINEHARDENING_ID) plasticityType case (PLASTICITY_KINEHARDENING_ID) plasticityType
of = phasememberAt(ipc,ip,el) of = material_phasememberAt(ipc,ip,el)
instance = phase_plasticityInstance(material_phase(ipc,ip,el)) instance = phase_plasticityInstance(material_phaseAt(ipc,el))
call plastic_kinehardening_LpAndItsTangent (Lp,dLp_dMp, Mp,instance,of) call plastic_kinehardening_LpAndItsTangent (Lp,dLp_dMp, Mp,instance,of)
case (PLASTICITY_NONLOCAL_ID) plasticityType case (PLASTICITY_NONLOCAL_ID) plasticityType
@ -357,13 +357,13 @@ subroutine constitutive_LpAndItsTangents(Lp, dLp_dS, dLp_dFi, &
temperature(ho)%p(tme),geometry_plastic_nonlocal_IPvolume0(ip,el),ip,el) temperature(ho)%p(tme),geometry_plastic_nonlocal_IPvolume0(ip,el),ip,el)
case (PLASTICITY_DISLOTWIN_ID) plasticityType case (PLASTICITY_DISLOTWIN_ID) plasticityType
of = phasememberAt(ipc,ip,el) of = material_phasememberAt(ipc,ip,el)
instance = phase_plasticityInstance(material_phase(ipc,ip,el)) instance = phase_plasticityInstance(material_phaseAt(ipc,el))
call plastic_dislotwin_LpAndItsTangent (Lp,dLp_dMp,Mp,temperature(ho)%p(tme),instance,of) call plastic_dislotwin_LpAndItsTangent (Lp,dLp_dMp,Mp,temperature(ho)%p(tme),instance,of)
case (PLASTICITY_DISLOUCLA_ID) plasticityType case (PLASTICITY_DISLOUCLA_ID) plasticityType
of = phasememberAt(ipc,ip,el) of = material_phasememberAt(ipc,ip,el)
instance = phase_plasticityInstance(material_phase(ipc,ip,el)) instance = phase_plasticityInstance(material_phaseAt(ipc,el))
call plastic_disloucla_LpAndItsTangent (Lp,dLp_dMp,Mp,temperature(ho)%p(tme),instance,of) call plastic_disloucla_LpAndItsTangent (Lp,dLp_dMp,Mp,temperature(ho)%p(tme),instance,of)
end select plasticityType end select plasticityType
@ -414,10 +414,10 @@ subroutine constitutive_LiAndItsTangents(Li, dLi_dS, dLi_dFi, &
dLi_dS = 0.0_pReal dLi_dS = 0.0_pReal
dLi_dFi = 0.0_pReal dLi_dFi = 0.0_pReal
plasticityType: select case (phase_plasticity(material_phase(ipc,ip,el))) plasticityType: select case (phase_plasticity(material_phaseAt(ipc,el)))
case (PLASTICITY_isotropic_ID) plasticityType case (PLASTICITY_isotropic_ID) plasticityType
of = phasememberAt(ipc,ip,el) of = material_phasememberAt(ipc,ip,el)
instance = phase_plasticityInstance(material_phase(ipc,ip,el)) instance = phase_plasticityInstance(material_phaseAt(ipc,el))
call plastic_isotropic_LiAndItsTangent(my_Li, my_dLi_dS, S ,instance,of) call plastic_isotropic_LiAndItsTangent(my_Li, my_dLi_dS, S ,instance,of)
case default plasticityType case default plasticityType
my_Li = 0.0_pReal my_Li = 0.0_pReal
@ -427,8 +427,8 @@ subroutine constitutive_LiAndItsTangents(Li, dLi_dS, dLi_dFi, &
Li = Li + my_Li Li = Li + my_Li
dLi_dS = dLi_dS + my_dLi_dS dLi_dS = dLi_dS + my_dLi_dS
KinematicsLoop: do k = 1, phase_Nkinematics(material_phase(ipc,ip,el)) KinematicsLoop: do k = 1, phase_Nkinematics(material_phaseAt(ipc,el))
kinematicsType: select case (phase_kinematics(k,material_phase(ipc,ip,el))) kinematicsType: select case (phase_kinematics(k,material_phaseAt(ipc,el)))
case (KINEMATICS_cleavage_opening_ID) kinematicsType case (KINEMATICS_cleavage_opening_ID) kinematicsType
call kinematics_cleavage_opening_LiAndItsTangent(my_Li, my_dLi_dS, S, ipc, ip, el) call kinematics_cleavage_opening_LiAndItsTangent(my_Li, my_dLi_dS, S, ipc, ip, el)
case (KINEMATICS_slipplane_opening_ID) kinematicsType case (KINEMATICS_slipplane_opening_ID) kinematicsType
@ -475,7 +475,7 @@ pure function constitutive_initialFi(ipc, ip, el)
homog, offset homog, offset
constitutive_initialFi = math_I3 constitutive_initialFi = math_I3
phase = material_phase(ipc,ip,el) phase = material_phaseAt(ipc,el)
KinematicsLoop: do k = 1, phase_Nkinematics(phase) !< Warning: small initial strain assumption KinematicsLoop: do k = 1, phase_Nkinematics(phase) !< Warning: small initial strain assumption
kinematicsType: select case (phase_kinematics(k,phase)) kinematicsType: select case (phase_kinematics(k,phase))
@ -546,8 +546,8 @@ subroutine constitutive_hooke_SandItsTangents(S, dS_dFe, dS_dFi, &
ho = material_homogenizationAt(el) ho = material_homogenizationAt(el)
C = math_66toSym3333(constitutive_homogenizedC(ipc,ip,el)) C = math_66toSym3333(constitutive_homogenizedC(ipc,ip,el))
DegradationLoop: do d = 1, phase_NstiffnessDegradations(material_phase(ipc,ip,el)) DegradationLoop: do d = 1, phase_NstiffnessDegradations(material_phaseAt(ipc,el))
degradationType: select case(phase_stiffnessDegradation(d,material_phase(ipc,ip,el))) degradationType: select case(phase_stiffnessDegradation(d,material_phaseAt(ipc,el)))
case (STIFFNESS_DEGRADATION_damage_ID) degradationType case (STIFFNESS_DEGRADATION_damage_ID) degradationType
C = C * damage(ho)%p(damageMapping(ho)%p(ip,el))**2 C = C * damage(ho)%p(damageMapping(ho)%p(ip,el))**2
end select degradationType end select degradationType
@ -597,31 +597,31 @@ subroutine constitutive_collectDotState(S, FeArray, Fi, FpArray, subdt, ipc, ip,
Mp = matmul(matmul(transpose(Fi),Fi),S) Mp = matmul(matmul(transpose(Fi),Fi),S)
plasticityType: select case (phase_plasticity(material_phase(ipc,ip,el))) plasticityType: select case (phase_plasticity(material_phaseAt(ipc,el)))
case (PLASTICITY_ISOTROPIC_ID) plasticityType case (PLASTICITY_ISOTROPIC_ID) plasticityType
of = phasememberAt(ipc,ip,el) of = material_phasememberAt(ipc,ip,el)
instance = phase_plasticityInstance(material_phase(ipc,ip,el)) instance = phase_plasticityInstance(material_phaseAt(ipc,el))
call plastic_isotropic_dotState (Mp,instance,of) call plastic_isotropic_dotState (Mp,instance,of)
case (PLASTICITY_PHENOPOWERLAW_ID) plasticityType case (PLASTICITY_PHENOPOWERLAW_ID) plasticityType
of = phasememberAt(ipc,ip,el) of = material_phasememberAt(ipc,ip,el)
instance = phase_plasticityInstance(material_phase(ipc,ip,el)) instance = phase_plasticityInstance(material_phaseAt(ipc,el))
call plastic_phenopowerlaw_dotState(Mp,instance,of) call plastic_phenopowerlaw_dotState(Mp,instance,of)
case (PLASTICITY_KINEHARDENING_ID) plasticityType case (PLASTICITY_KINEHARDENING_ID) plasticityType
of = phasememberAt(ipc,ip,el) of = material_phasememberAt(ipc,ip,el)
instance = phase_plasticityInstance(material_phase(ipc,ip,el)) instance = phase_plasticityInstance(material_phaseAt(ipc,el))
call plastic_kinehardening_dotState(Mp,instance,of) call plastic_kinehardening_dotState(Mp,instance,of)
case (PLASTICITY_DISLOTWIN_ID) plasticityType case (PLASTICITY_DISLOTWIN_ID) plasticityType
of = phasememberAt(ipc,ip,el) of = material_phasememberAt(ipc,ip,el)
instance = phase_plasticityInstance(material_phase(ipc,ip,el)) instance = phase_plasticityInstance(material_phaseAt(ipc,el))
call plastic_dislotwin_dotState (Mp,temperature(ho)%p(tme),instance,of) call plastic_dislotwin_dotState (Mp,temperature(ho)%p(tme),instance,of)
case (PLASTICITY_DISLOUCLA_ID) plasticityType case (PLASTICITY_DISLOUCLA_ID) plasticityType
of = phasememberAt(ipc,ip,el) of = material_phasememberAt(ipc,ip,el)
instance = phase_plasticityInstance(material_phase(ipc,ip,el)) instance = phase_plasticityInstance(material_phaseAt(ipc,el))
call plastic_disloucla_dotState (Mp,temperature(ho)%p(tme),instance,of) call plastic_disloucla_dotState (Mp,temperature(ho)%p(tme),instance,of)
case (PLASTICITY_NONLOCAL_ID) plasticityType case (PLASTICITY_NONLOCAL_ID) plasticityType
@ -629,9 +629,9 @@ subroutine constitutive_collectDotState(S, FeArray, Fi, FpArray, subdt, ipc, ip,
subdt,ip,el) subdt,ip,el)
end select plasticityType end select plasticityType
SourceLoop: do i = 1, phase_Nsources(material_phase(ipc,ip,el)) SourceLoop: do i = 1, phase_Nsources(material_phaseAt(ipc,el))
sourceType: select case (phase_source(i,material_phase(ipc,ip,el))) sourceType: select case (phase_source(i,material_phaseAt(ipc,el)))
case (SOURCE_damage_anisoBrittle_ID) sourceType case (SOURCE_damage_anisoBrittle_ID) sourceType
call source_damage_anisoBrittle_dotState (S, ipc, ip, el) !< correct stress? call source_damage_anisoBrittle_dotState (S, ipc, ip, el) !< correct stress?
@ -643,8 +643,8 @@ subroutine constitutive_collectDotState(S, FeArray, Fi, FpArray, subdt, ipc, ip,
call source_damage_anisoDuctile_dotState ( ipc, ip, el) call source_damage_anisoDuctile_dotState ( ipc, ip, el)
case (SOURCE_thermal_externalheat_ID) sourceType case (SOURCE_thermal_externalheat_ID) sourceType
of = phasememberAt(ipc,ip,el) of = material_phasememberAt(ipc,ip,el)
call source_thermal_externalheat_dotState(material_phase(ipc,ip,el),of) call source_thermal_externalheat_dotState(material_phaseAt(ipc,el),of)
end select sourceType end select sourceType
@ -674,11 +674,11 @@ subroutine constitutive_collectDeltaState(S, Fe, Fi, ipc, ip, el)
Mp = matmul(matmul(transpose(Fi),Fi),S) Mp = matmul(matmul(transpose(Fi),Fi),S)
plasticityType: select case (phase_plasticity(material_phase(ipc,ip,el))) plasticityType: select case (phase_plasticity(material_phaseAt(ipc,el)))
case (PLASTICITY_KINEHARDENING_ID) plasticityType case (PLASTICITY_KINEHARDENING_ID) plasticityType
of = phasememberAt(ipc,ip,el) of = material_phasememberAt(ipc,ip,el)
instance = phase_plasticityInstance(material_phase(ipc,ip,el)) instance = phase_plasticityInstance(material_phaseAt(ipc,el))
call plastic_kinehardening_deltaState(Mp,instance,of) call plastic_kinehardening_deltaState(Mp,instance,of)
case (PLASTICITY_NONLOCAL_ID) plasticityType case (PLASTICITY_NONLOCAL_ID) plasticityType
@ -686,9 +686,9 @@ subroutine constitutive_collectDeltaState(S, Fe, Fi, ipc, ip, el)
end select plasticityType end select plasticityType
sourceLoop: do i = 1, phase_Nsources(material_phase(ipc,ip,el)) sourceLoop: do i = 1, phase_Nsources(material_phaseAt(ipc,el))
sourceType: select case (phase_source(i,material_phase(ipc,ip,el))) sourceType: select case (phase_source(i,material_phaseAt(ipc,el)))
case (SOURCE_damage_isoBrittle_ID) sourceType case (SOURCE_damage_isoBrittle_ID) sourceType
call source_damage_isoBrittle_deltaState (constitutive_homogenizedC(ipc,ip,el), Fe, & call source_damage_isoBrittle_deltaState (constitutive_homogenizedC(ipc,ip,el), Fe, &
@ -710,8 +710,8 @@ function constitutive_postResults(S, Fi, ipc, ip, el)
ipc, & !< component-ID of integration point ipc, & !< component-ID of integration point
ip, & !< integration point ip, & !< integration point
el !< element el !< element
real(pReal), dimension(plasticState(material_phase(ipc,ip,el))%sizePostResults + & real(pReal), dimension(plasticState(material_phaseAt(ipc,el))%sizePostResults + &
sum(sourceState(material_phase(ipc,ip,el))%p(:)%sizePostResults)) :: & sum(sourceState(material_phaseAt(ipc,el))%p(:)%sizePostResults)) :: &
constitutive_postResults constitutive_postResults
real(pReal), intent(in), dimension(3,3) :: & real(pReal), intent(in), dimension(3,3) :: &
Fi !< intermediate deformation gradient Fi !< intermediate deformation gradient
@ -734,12 +734,12 @@ function constitutive_postResults(S, Fi, ipc, ip, el)
tme = thermalMapping(ho)%p(ip,el) tme = thermalMapping(ho)%p(ip,el)
startPos = 1 startPos = 1
endPos = plasticState(material_phase(ipc,ip,el))%sizePostResults endPos = plasticState(material_phaseAt(ipc,el))%sizePostResults
of = phasememberAt(ipc,ip,el) of = material_phasememberAt(ipc,ip,el)
instance = phase_plasticityInstance(material_phase(ipc,ip,el)) instance = phase_plasticityInstance(material_phaseAt(ipc,el))
plasticityType: select case (phase_plasticity(material_phase(ipc,ip,el))) plasticityType: select case (phase_plasticity(material_phaseAt(ipc,el)))
case (PLASTICITY_ISOTROPIC_ID) plasticityType case (PLASTICITY_ISOTROPIC_ID) plasticityType
constitutive_postResults(startPos:endPos) = & constitutive_postResults(startPos:endPos) = &
plastic_isotropic_postResults(Mp,instance,of) plastic_isotropic_postResults(Mp,instance,of)
@ -762,23 +762,23 @@ function constitutive_postResults(S, Fi, ipc, ip, el)
case (PLASTICITY_NONLOCAL_ID) plasticityType case (PLASTICITY_NONLOCAL_ID) plasticityType
constitutive_postResults(startPos:endPos) = & constitutive_postResults(startPos:endPos) = &
plastic_nonlocal_postResults (material_phase(ipc,ip,el),instance,of) plastic_nonlocal_postResults (material_phaseAt(ipc,el),instance,of)
end select plasticityType end select plasticityType
SourceLoop: do i = 1, phase_Nsources(material_phase(ipc,ip,el)) SourceLoop: do i = 1, phase_Nsources(material_phaseAt(ipc,el))
startPos = endPos + 1 startPos = endPos + 1
endPos = endPos + sourceState(material_phase(ipc,ip,el))%p(i)%sizePostResults endPos = endPos + sourceState(material_phaseAt(ipc,el))%p(i)%sizePostResults
of = phasememberAt(ipc,ip,el) of = material_phasememberAt(ipc,ip,el)
sourceType: select case (phase_source(i,material_phase(ipc,ip,el))) sourceType: select case (phase_source(i,material_phaseAt(ipc,el)))
case (SOURCE_damage_isoBrittle_ID) sourceType case (SOURCE_damage_isoBrittle_ID) sourceType
constitutive_postResults(startPos:endPos) = source_damage_isoBrittle_postResults(material_phase(ipc,ip,el),of) constitutive_postResults(startPos:endPos) = source_damage_isoBrittle_postResults(material_phaseAt(ipc,el),of)
case (SOURCE_damage_isoDuctile_ID) sourceType case (SOURCE_damage_isoDuctile_ID) sourceType
constitutive_postResults(startPos:endPos) = source_damage_isoDuctile_postResults(material_phase(ipc,ip,el),of) constitutive_postResults(startPos:endPos) = source_damage_isoDuctile_postResults(material_phaseAt(ipc,el),of)
case (SOURCE_damage_anisoBrittle_ID) sourceType case (SOURCE_damage_anisoBrittle_ID) sourceType
constitutive_postResults(startPos:endPos) = source_damage_anisoBrittle_postResults(material_phase(ipc,ip,el),of) constitutive_postResults(startPos:endPos) = source_damage_anisoBrittle_postResults(material_phaseAt(ipc,el),of)
case (SOURCE_damage_anisoDuctile_ID) sourceType case (SOURCE_damage_anisoDuctile_ID) sourceType
constitutive_postResults(startPos:endPos) = source_damage_anisoDuctile_postResults(material_phase(ipc,ip,el),of) constitutive_postResults(startPos:endPos) = source_damage_anisoDuctile_postResults(material_phaseAt(ipc,el),of)
end select sourceType end select sourceType
enddo SourceLoop enddo SourceLoop

View File

@ -386,7 +386,7 @@ subroutine crystallite_init
crystallite_Fp0(1:3,1:3,c,i,e) = math_EulerToR(material_EulerAngles(1:3,c,i,e)) ! plastic def gradient reflects init orientation crystallite_Fp0(1:3,1:3,c,i,e) = math_EulerToR(material_EulerAngles(1:3,c,i,e)) ! plastic def gradient reflects init orientation
crystallite_Fi0(1:3,1:3,c,i,e) = constitutive_initialFi(c,i,e) crystallite_Fi0(1:3,1:3,c,i,e) = constitutive_initialFi(c,i,e)
crystallite_F0(1:3,1:3,c,i,e) = math_I3 crystallite_F0(1:3,1:3,c,i,e) = math_I3
crystallite_localPlasticity(c,i,e) = phase_localPlasticity(material_phase(c,i,e)) crystallite_localPlasticity(c,i,e) = phase_localPlasticity(material_phaseAt(c,e))
crystallite_Fe(1:3,1:3,c,i,e) = math_inv33(matmul(crystallite_Fi0(1:3,1:3,c,i,e), & crystallite_Fe(1:3,1:3,c,i,e) = math_inv33(matmul(crystallite_Fi0(1:3,1:3,c,i,e), &
crystallite_Fp0(1:3,1:3,c,i,e))) ! assuming that euler angles are given in internal strain free configuration crystallite_Fp0(1:3,1:3,c,i,e))) ! assuming that euler angles are given in internal strain free configuration
crystallite_Fp(1:3,1:3,c,i,e) = crystallite_Fp0(1:3,1:3,c,i,e) crystallite_Fp(1:3,1:3,c,i,e) = crystallite_Fp0(1:3,1:3,c,i,e)
@ -483,12 +483,12 @@ function crystallite_stress(dummyArgumentToPreventInternalCompilerErrorWithGCC)
elementLooping1: do e = FEsolving_execElem(1),FEsolving_execElem(2) elementLooping1: do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e); do c = 1,homogenization_Ngrains(material_homogenizationAt(e)) do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e); do c = 1,homogenization_Ngrains(material_homogenizationAt(e))
homogenizationRequestsCalculation: if (crystallite_requested(c,i,e)) then homogenizationRequestsCalculation: if (crystallite_requested(c,i,e)) then
plasticState (phaseAt(c,i,e))%subState0( :,phasememberAt(c,i,e)) = & plasticState (material_phaseAt(c,e))%subState0( :,material_phaseMemberAt(c,i,e)) = &
plasticState (phaseAt(c,i,e))%partionedState0(:,phasememberAt(c,i,e)) plasticState (material_phaseAt(c,e))%partionedState0(:,material_phaseMemberAt(c,i,e))
do s = 1, phase_Nsources(phaseAt(c,i,e)) do s = 1, phase_Nsources(material_phaseAt(c,e))
sourceState(phaseAt(c,i,e))%p(s)%subState0( :,phasememberAt(c,i,e)) = & sourceState(material_phaseAt(c,e))%p(s)%subState0( :,material_phaseMemberAt(c,i,e)) = &
sourceState(phaseAt(c,i,e))%p(s)%partionedState0(:,phasememberAt(c,i,e)) sourceState(material_phaseAt(c,e))%p(s)%partionedState0(:,material_phaseMemberAt(c,i,e))
enddo enddo
crystallite_subFp0(1:3,1:3,c,i,e) = crystallite_partionedFp0(1:3,1:3,c,i,e) crystallite_subFp0(1:3,1:3,c,i,e) = crystallite_partionedFp0(1:3,1:3,c,i,e)
crystallite_subLp0(1:3,1:3,c,i,e) = crystallite_partionedLp0(1:3,1:3,c,i,e) crystallite_subLp0(1:3,1:3,c,i,e) = crystallite_partionedLp0(1:3,1:3,c,i,e)
@ -543,11 +543,11 @@ function crystallite_stress(dummyArgumentToPreventInternalCompilerErrorWithGCC)
crystallite_subFi0(1:3,1:3,c,i,e) = crystallite_Fi (1:3,1:3,c,i,e) crystallite_subFi0(1:3,1:3,c,i,e) = crystallite_Fi (1:3,1:3,c,i,e)
crystallite_subS0 (1:3,1:3,c,i,e) = crystallite_S (1:3,1:3,c,i,e) crystallite_subS0 (1:3,1:3,c,i,e) = crystallite_S (1:3,1:3,c,i,e)
!if abbrevation, make c and p private in omp !if abbrevation, make c and p private in omp
plasticState( phaseAt(c,i,e))%subState0(:,phasememberAt(c,i,e)) & plasticState( material_phaseAt(c,e))%subState0(:,material_phaseMemberAt(c,i,e)) &
= plasticState(phaseAt(c,i,e))%state( :,phasememberAt(c,i,e)) = plasticState(material_phaseAt(c,e))%state( :,material_phaseMemberAt(c,i,e))
do s = 1, phase_Nsources(phaseAt(c,i,e)) do s = 1, phase_Nsources(material_phaseAt(c,e))
sourceState( phaseAt(c,i,e))%p(s)%subState0(:,phasememberAt(c,i,e)) & sourceState( material_phaseAt(c,e))%p(s)%subState0(:,material_phaseMemberAt(c,i,e)) &
= sourceState(phaseAt(c,i,e))%p(s)%state( :,phasememberAt(c,i,e)) = sourceState(material_phaseAt(c,e))%p(s)%state( :,material_phaseMemberAt(c,i,e))
enddo enddo
#ifdef DEBUG #ifdef DEBUG
if (iand(debug_level(debug_crystallite),debug_levelBasic) /= 0 & if (iand(debug_level(debug_crystallite),debug_levelBasic) /= 0 &
@ -572,11 +572,11 @@ function crystallite_stress(dummyArgumentToPreventInternalCompilerErrorWithGCC)
crystallite_Lp (1:3,1:3,c,i,e) = crystallite_subLp0(1:3,1:3,c,i,e) crystallite_Lp (1:3,1:3,c,i,e) = crystallite_subLp0(1:3,1:3,c,i,e)
crystallite_Li (1:3,1:3,c,i,e) = crystallite_subLi0(1:3,1:3,c,i,e) crystallite_Li (1:3,1:3,c,i,e) = crystallite_subLi0(1:3,1:3,c,i,e)
endif endif
plasticState (phaseAt(c,i,e))%state( :,phasememberAt(c,i,e)) & plasticState (material_phaseAt(c,e))%state( :,material_phaseMemberAt(c,i,e)) &
= plasticState(phaseAt(c,i,e))%subState0(:,phasememberAt(c,i,e)) = plasticState(material_phaseAt(c,e))%subState0(:,material_phaseMemberAt(c,i,e))
do s = 1, phase_Nsources(phaseAt(c,i,e)) do s = 1, phase_Nsources(material_phaseAt(c,e))
sourceState( phaseAt(c,i,e))%p(s)%state( :,phasememberAt(c,i,e)) & sourceState( material_phaseAt(c,e))%p(s)%state( :,material_phaseMemberAt(c,i,e)) &
= sourceState(phaseAt(c,i,e))%p(s)%subState0(:,phasememberAt(c,i,e)) = sourceState(material_phaseAt(c,e))%p(s)%subState0(:,material_phaseMemberAt(c,i,e))
enddo enddo
! cant restore dotState here, since not yet calculated in first cutback after initialization ! cant restore dotState here, since not yet calculated in first cutback after initialization
@ -839,7 +839,7 @@ subroutine crystallite_orientations
!$OMP PARALLEL DO !$OMP PARALLEL DO
do e = FEsolving_execElem(1),FEsolving_execElem(2) do e = FEsolving_execElem(1),FEsolving_execElem(2)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
if (plasticState(material_phase(1,i,e))%nonLocal) & ! if nonlocal model if (plasticState(material_phaseAt(1,e))%nonLocal) & ! if nonlocal model
call plastic_nonlocal_updateCompatibility(crystallite_orientation,i,e) call plastic_nonlocal_updateCompatibility(crystallite_orientation,i,e)
enddo; enddo enddo; enddo
!$OMP END PARALLEL DO !$OMP END PARALLEL DO
@ -879,8 +879,8 @@ function crystallite_postResults(ipc, ip, el)
ipc !< grain index ipc !< grain index
real(pReal), dimension(1+crystallite_sizePostResults(microstructure_crystallite(discretization_microstructureAt(el))) + & real(pReal), dimension(1+crystallite_sizePostResults(microstructure_crystallite(discretization_microstructureAt(el))) + &
1+plasticState(material_phase(ipc,ip,el))%sizePostResults + & 1+plasticState(material_phaseAt(ipc,el))%sizePostResults + &
sum(sourceState(material_phase(ipc,ip,el))%p(:)%sizePostResults)) :: & sum(sourceState(material_phaseAt(ipc,el))%p(:)%sizePostResults)) :: &
crystallite_postResults crystallite_postResults
integer :: & integer :: &
o, & o, &
@ -893,7 +893,7 @@ function crystallite_postResults(ipc, ip, el)
crystID = microstructure_crystallite(discretization_microstructureAt(el)) crystID = microstructure_crystallite(discretization_microstructureAt(el))
crystallite_postResults = 0.0_pReal crystallite_postResults = 0.0_pReal
crystallite_postResults(1) = real(crystallite_sizePostResults(crystID),pReal) ! header-like information (length) crystallite_postResults(1) = real(crystallite_sizePostResults(crystID),pReal) ! header-like information (length)
c = 1 c = 1
do o = 1,crystallite_Noutput(crystID) do o = 1,crystallite_Noutput(crystID)
@ -901,7 +901,7 @@ function crystallite_postResults(ipc, ip, el)
select case(crystallite_outputID(o,crystID)) select case(crystallite_outputID(o,crystID))
case (phase_ID) case (phase_ID)
mySize = 1 mySize = 1
crystallite_postResults(c+1) = real(material_phase(ipc,ip,el),pReal) ! phaseID of grain crystallite_postResults(c+1) = real(material_phaseAt(ipc,el),pReal) ! phaseID of grain
case (texture_ID) case (texture_ID)
mySize = 1 mySize = 1
crystallite_postResults(c+1) = real(material_texture(ipc,ip,el),pReal) ! textureID of grain crystallite_postResults(c+1) = real(material_texture(ipc,ip,el),pReal) ! textureID of grain
@ -967,7 +967,7 @@ function crystallite_postResults(ipc, ip, el)
c = c + mySize c = c + mySize
enddo enddo
crystallite_postResults(c+1) = real(plasticState(material_phase(ipc,ip,el))%sizePostResults,pReal) ! size of constitutive results crystallite_postResults(c+1) = real(plasticState(material_phaseAt(ipc,el))%sizePostResults,pReal) ! size of constitutive results
c = c + 1 c = c + 1
if (size(crystallite_postResults)-c > 0) & if (size(crystallite_postResults)-c > 0) &
crystallite_postResults(c+1:size(crystallite_postResults)) = & crystallite_postResults(c+1:size(crystallite_postResults)) = &
@ -1555,7 +1555,7 @@ subroutine integrateStateFPI
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
do g = 1,homogenization_Ngrains(material_homogenizationAt(e)) do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) then if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) then
p = phaseAt(g,i,e); c = phasememberAt(g,i,e) p = material_phaseAt(g,e); c = material_phaseMemberAt(g,i,e)
plasticState(p)%previousDotState2(:,c) = merge(plasticState(p)%previousDotState(:,c),& plasticState(p)%previousDotState2(:,c) = merge(plasticState(p)%previousDotState(:,c),&
0.0_pReal,& 0.0_pReal,&
@ -1583,7 +1583,7 @@ subroutine integrateStateFPI
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
do g = 1,homogenization_Ngrains(material_homogenizationAt(e)) do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) then if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) then
p = phaseAt(g,i,e); c = phasememberAt(g,i,e) p = material_phaseAt(g,e); c = material_phaseMemberAt(g,i,e)
sizeDotState = plasticState(p)%sizeDotState sizeDotState = plasticState(p)%sizeDotState
zeta = damper(plasticState(p)%dotState (:,c), & zeta = damper(plasticState(p)%dotState (:,c), &
@ -1746,7 +1746,7 @@ subroutine integrateStateAdaptiveEuler
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
do g = 1,homogenization_Ngrains(material_homogenizationAt(e)) do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
if (crystallite_todo(g,i,e)) then if (crystallite_todo(g,i,e)) then
p = phaseAt(g,i,e); c = phasememberAt(g,i,e) p = material_phaseAt(g,e); c = material_phaseMemberAt(g,i,e)
sizeDotState = plasticState(p)%sizeDotState sizeDotState = plasticState(p)%sizeDotState
residuum_plastic(1:sizeDotState,g,i,e) = plasticState(p)%dotstate(1:sizeDotState,c) & residuum_plastic(1:sizeDotState,g,i,e) = plasticState(p)%dotstate(1:sizeDotState,c) &
@ -1775,7 +1775,7 @@ subroutine integrateStateAdaptiveEuler
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
do g = 1,homogenization_Ngrains(material_homogenizationAt(e)) do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
if (crystallite_todo(g,i,e)) then if (crystallite_todo(g,i,e)) then
p = phaseAt(g,i,e); c = phasememberAt(g,i,e) p = material_phaseAt(g,e); c = material_phaseMemberAt(g,i,e)
sizeDotState = plasticState(p)%sizeDotState sizeDotState = plasticState(p)%sizeDotState
residuum_plastic(1:sizeDotState,g,i,e) = residuum_plastic(1:sizeDotState,g,i,e) & residuum_plastic(1:sizeDotState,g,i,e) = residuum_plastic(1:sizeDotState,g,i,e) &
@ -1835,7 +1835,7 @@ subroutine integrateStateRK4
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
do g = 1,homogenization_Ngrains(material_homogenizationAt(e)) do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
if (crystallite_todo(g,i,e)) then if (crystallite_todo(g,i,e)) then
p = phaseAt(g,i,e); c = phasememberAt(g,i,e) p = material_phaseAt(g,e); c = material_phaseMemberAt(g,i,e)
plasticState(p)%RK4dotState(:,c) = WEIGHT(n)*plasticState(p)%dotState(:,c) & plasticState(p)%RK4dotState(:,c) = WEIGHT(n)*plasticState(p)%dotState(:,c) &
+ merge(plasticState(p)%RK4dotState(:,c),0.0_pReal,n>1) + merge(plasticState(p)%RK4dotState(:,c),0.0_pReal,n>1)
@ -1926,7 +1926,7 @@ subroutine integrateStateRKCK45
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
do g = 1,homogenization_Ngrains(material_homogenizationAt(e)) do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
if (crystallite_todo(g,i,e)) then if (crystallite_todo(g,i,e)) then
p = phaseAt(g,i,e); cc = phasememberAt(g,i,e) p = material_phaseAt(g,e); cc = material_phaseMemberAt(g,i,e)
plasticState(p)%RKCK45dotState(stage,:,cc) = plasticState(p)%dotState(:,cc) plasticState(p)%RKCK45dotState(stage,:,cc) = plasticState(p)%dotState(:,cc)
plasticState(p)%dotState(:,cc) = A(1,stage) * plasticState(p)%RKCK45dotState(1,:,cc) plasticState(p)%dotState(:,cc) = A(1,stage) * plasticState(p)%RKCK45dotState(1,:,cc)
@ -1966,7 +1966,7 @@ subroutine integrateStateRKCK45
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
do g = 1,homogenization_Ngrains(material_homogenizationAt(e)) do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
if (crystallite_todo(g,i,e)) then if (crystallite_todo(g,i,e)) then
p = phaseAt(g,i,e); cc = phasememberAt(g,i,e) p = material_phaseAt(g,e); cc = material_phaseMemberAt(g,i,e)
sizeDotState = plasticState(p)%sizeDotState sizeDotState = plasticState(p)%sizeDotState
@ -2005,7 +2005,7 @@ subroutine integrateStateRKCK45
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
do g = 1,homogenization_Ngrains(material_homogenizationAt(e)) do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
if (crystallite_todo(g,i,e)) then if (crystallite_todo(g,i,e)) then
p = phaseAt(g,i,e); cc = phasememberAt(g,i,e) p = material_phaseAt(g,e); cc = material_phaseMemberAt(g,i,e)
sizeDotState = plasticState(p)%sizeDotState sizeDotState = plasticState(p)%sizeDotState
@ -2163,7 +2163,7 @@ subroutine update_state(timeFraction)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e) do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
do g = 1,homogenization_Ngrains(material_homogenizationAt(e)) do g = 1,homogenization_Ngrains(material_homogenizationAt(e))
if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) then if (crystallite_todo(g,i,e) .and. .not. crystallite_converged(g,i,e)) then
p = phaseAt(g,i,e); c = phasememberAt(g,i,e) p = material_phaseAt(g,e); c = material_phaseMemberAt(g,i,e)
mySize = plasticState(p)%sizeDotState mySize = plasticState(p)%sizeDotState
plasticState(p)%state(1:mySize,c) = plasticState(p)%subState0(1:mySize,c) & plasticState(p)%state(1:mySize,c) = plasticState(p)%subState0(1:mySize,c) &
@ -2214,7 +2214,7 @@ subroutine update_dotState(timeFraction)
crystallite_Fi(1:3,1:3,g,i,e), & crystallite_Fi(1:3,1:3,g,i,e), &
crystallite_Fp, & crystallite_Fp, &
crystallite_subdt(g,i,e)*timeFraction, g,i,e) crystallite_subdt(g,i,e)*timeFraction, g,i,e)
p = phaseAt(g,i,e); c = phasememberAt(g,i,e) p = material_phaseAt(g,e); c = material_phaseMemberAt(g,i,e)
NaN = any(IEEE_is_NaN(plasticState(p)%dotState(:,c))) NaN = any(IEEE_is_NaN(plasticState(p)%dotState(:,c)))
do s = 1, phase_Nsources(p) do s = 1, phase_Nsources(p)
NaN = NaN .or. any(IEEE_is_NaN(sourceState(p)%p(s)%dotState(:,c))) NaN = NaN .or. any(IEEE_is_NaN(sourceState(p)%p(s)%dotState(:,c)))
@ -2259,7 +2259,7 @@ subroutine update_deltaState
crystallite_Fe(1:3,1:3,g,i,e), & crystallite_Fe(1:3,1:3,g,i,e), &
crystallite_Fi(1:3,1:3,g,i,e), & crystallite_Fi(1:3,1:3,g,i,e), &
g,i,e) g,i,e)
p = phaseAt(g,i,e); c = phasememberAt(g,i,e) p = material_phaseAt(g,e); c = material_phaseMemberAt(g,i,e)
myOffset = plasticState(p)%offsetDeltaState myOffset = plasticState(p)%offsetDeltaState
mySize = plasticState(p)%sizeDeltaState mySize = plasticState(p)%sizeDeltaState
NaN = any(IEEE_is_NaN(plasticState(p)%deltaState(1:mySize,c))) NaN = any(IEEE_is_NaN(plasticState(p)%deltaState(1:mySize,c)))
@ -2311,8 +2311,8 @@ logical function stateJump(ipc,ip,el)
myOffset, & myOffset, &
mySize mySize
c = phasememberAt(ipc,ip,el) c = material_phaseMemberAt(ipc,ip,el)
p = phaseAt(ipc,ip,el) p = material_phaseAt(ipc,el)
call constitutive_collectDeltaState(crystallite_S(1:3,1:3,ipc,ip,el), & call constitutive_collectDeltaState(crystallite_S(1:3,1:3,ipc,ip,el), &
crystallite_Fe(1:3,1:3,ipc,ip,el), & crystallite_Fe(1:3,1:3,ipc,ip,el), &

View File

@ -178,8 +178,8 @@ subroutine damage_local_getSourceAndItsTangent(phiDot, dPhiDot_dPhi, phi, ip, el
phiDot = 0.0_pReal phiDot = 0.0_pReal
dPhiDot_dPhi = 0.0_pReal dPhiDot_dPhi = 0.0_pReal
do grain = 1, homogenization_Ngrains(material_homogenizationAt(el)) do grain = 1, homogenization_Ngrains(material_homogenizationAt(el))
phase = phaseAt(grain,ip,el) phase = material_phaseAt(grain,el)
constituent = phasememberAt(grain,ip,el) constituent = material_phasememberAt(grain,ip,el)
do source = 1, phase_Nsources(phase) do source = 1, phase_Nsources(phase)
select case(phase_source(source,phase)) select case(phase_source(source,phase))
case (SOURCE_damage_isoBrittle_ID) case (SOURCE_damage_isoBrittle_ID)

View File

@ -144,8 +144,8 @@ subroutine damage_nonlocal_getSourceAndItsTangent(phiDot, dPhiDot_dPhi, phi, ip,
phiDot = 0.0_pReal phiDot = 0.0_pReal
dPhiDot_dPhi = 0.0_pReal dPhiDot_dPhi = 0.0_pReal
do grain = 1, homogenization_Ngrains(material_homogenizationAt(el)) do grain = 1, homogenization_Ngrains(material_homogenizationAt(el))
phase = phaseAt(grain,ip,el) phase = material_phaseAt(grain,el)
constituent = phasememberAt(grain,ip,el) constituent = material_phasememberAt(grain,ip,el)
do source = 1, phase_Nsources(phase) do source = 1, phase_Nsources(phase)
select case(phase_source(source,phase)) select case(phase_source(source,phase))
case (SOURCE_damage_isoBrittle_ID) case (SOURCE_damage_isoBrittle_ID)

View File

@ -326,11 +326,11 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt)
do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e); do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e);
do g = 1,myNgrains do g = 1,myNgrains
plasticState (phaseAt(g,i,e))%partionedState0(:,phasememberAt(g,i,e)) = & plasticState (material_phaseAt(g,e))%partionedState0(:,material_phasememberAt(g,i,e)) = &
plasticState (phaseAt(g,i,e))%state0( :,phasememberAt(g,i,e)) plasticState (material_phaseAt(g,e))%state0( :,material_phasememberAt(g,i,e))
do mySource = 1, phase_Nsources(phaseAt(g,i,e)) do mySource = 1, phase_Nsources(material_phaseAt(g,e))
sourceState(phaseAt(g,i,e))%p(mySource)%partionedState0(:,phasememberAt(g,i,e)) = & sourceState(material_phaseAt(g,e))%p(mySource)%partionedState0(:,material_phasememberAt(g,i,e)) = &
sourceState(phaseAt(g,i,e))%p(mySource)%state0( :,phasememberAt(g,i,e)) sourceState(material_phaseAt(g,e))%p(mySource)%state0( :,material_phasememberAt(g,i,e))
enddo enddo
crystallite_partionedFp0(1:3,1:3,g,i,e) = crystallite_Fp0(1:3,1:3,g,i,e) crystallite_partionedFp0(1:3,1:3,g,i,e) = crystallite_Fp0(1:3,1:3,g,i,e)
@ -412,11 +412,11 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt)
crystallite_S (1:3,1:3,1:myNgrains,i,e) crystallite_S (1:3,1:3,1:myNgrains,i,e)
do g = 1,myNgrains do g = 1,myNgrains
plasticState (phaseAt(g,i,e))%partionedState0(:,phasememberAt(g,i,e)) = & plasticState (material_phaseAt(g,e))%partionedState0(:,material_phasememberAt(g,i,e)) = &
plasticState (phaseAt(g,i,e))%state (:,phasememberAt(g,i,e)) plasticState (material_phaseAt(g,e))%state (:,material_phasememberAt(g,i,e))
do mySource = 1, phase_Nsources(phaseAt(g,i,e)) do mySource = 1, phase_Nsources(material_phaseAt(g,e))
sourceState(phaseAt(g,i,e))%p(mySource)%partionedState0(:,phasememberAt(g,i,e)) = & sourceState(material_phaseAt(g,e))%p(mySource)%partionedState0(:,material_phasememberAt(g,i,e)) = &
sourceState(phaseAt(g,i,e))%p(mySource)%state (:,phasememberAt(g,i,e)) sourceState(material_phaseAt(g,e))%p(mySource)%state (:,material_phasememberAt(g,i,e))
enddo enddo
enddo enddo
@ -475,11 +475,11 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt)
crystallite_S(1:3,1:3,1:myNgrains,i,e) = & crystallite_S(1:3,1:3,1:myNgrains,i,e) = &
crystallite_partionedS0(1:3,1:3,1:myNgrains,i,e) crystallite_partionedS0(1:3,1:3,1:myNgrains,i,e)
do g = 1, myNgrains do g = 1, myNgrains
plasticState (phaseAt(g,i,e))%state( :,phasememberAt(g,i,e)) = & plasticState (material_phaseAt(g,e))%state( :,material_phasememberAt(g,i,e)) = &
plasticState (phaseAt(g,i,e))%partionedState0(:,phasememberAt(g,i,e)) plasticState (material_phaseAt(g,e))%partionedState0(:,material_phasememberAt(g,i,e))
do mySource = 1, phase_Nsources(phaseAt(g,i,e)) do mySource = 1, phase_Nsources(material_phaseAt(g,e))
sourceState(phaseAt(g,i,e))%p(mySource)%state( :,phasememberAt(g,i,e)) = & sourceState(material_phaseAt(g,e))%p(mySource)%state( :,material_phasememberAt(g,i,e)) = &
sourceState(phaseAt(g,i,e))%p(mySource)%partionedState0(:,phasememberAt(g,i,e)) sourceState(material_phaseAt(g,e))%p(mySource)%partionedState0(:,material_phasememberAt(g,i,e))
enddo enddo
enddo enddo
if(homogState(material_homogenizationAt(e))%sizeState > 0) & if(homogState(material_homogenizationAt(e))%sizeState > 0) &

View File

@ -138,7 +138,7 @@ subroutine kinematics_cleavage_opening_LiAndItsTangent(Ld, dLd_dTstar, S, ipc, i
traction_d, traction_t, traction_n, traction_crit, & traction_d, traction_t, traction_n, traction_crit, &
udotd, dudotd_dt, udott, dudott_dt, udotn, dudotn_dt udotd, dudotd_dt, udott, dudott_dt, udotn, dudotn_dt
phase = material_phase(ipc,ip,el) phase = material_phaseAt(ipc,el)
instance = kinematics_cleavage_opening_instance(phase) instance = kinematics_cleavage_opening_instance(phase)
homog = material_homogenizationAt(el) homog = material_homogenizationAt(el)
damageOffset = damageMapping(homog)%p(ip,el) damageOffset = damageMapping(homog)%p(ip,el)

View File

@ -124,7 +124,7 @@ subroutine kinematics_slipplane_opening_LiAndItsTangent(Ld, dLd_dTstar, S, ipc,
traction_d, traction_t, traction_n, traction_crit, & traction_d, traction_t, traction_n, traction_crit, &
udotd, dudotd_dt, udott, dudott_dt, udotn, dudotn_dt udotd, dudotd_dt, udott, dudott_dt, udotn, dudotn_dt
phase = material_phase(ipc,ip,el) phase = material_phaseAt(ipc,el)
instance = kinematics_slipplane_opening_instance(phase) instance = kinematics_slipplane_opening_instance(phase)
homog = material_homogenizationAt(el) homog = material_homogenizationAt(el)
damageOffset = damageMapping(homog)%p(ip,el) damageOffset = damageMapping(homog)%p(ip,el)

View File

@ -112,7 +112,7 @@ subroutine kinematics_thermal_expansion_LiAndItsTangent(Li, dLi_dTstar, ipc, ip,
real(pReal) :: & real(pReal) :: &
T, TRef, TDot T, TRef, TDot
phase = material_phase(ipc,ip,el) phase = material_phaseAt(ipc,el)
homog = material_homogenizationAt(el) homog = material_homogenizationAt(el)
offset = thermalMapping(homog)%p(ip,el) offset = thermalMapping(homog)%p(ip,el)
T = temperature(homog)%p(offset) T = temperature(homog)%p(offset)

View File

@ -180,9 +180,6 @@ module material
homogenization_active homogenization_active
! BEGIN DEPRECATED ! BEGIN DEPRECATED
integer, dimension(:,:,:), allocatable, public :: phaseAt !< phase ID of every material point (ipc,ip,el)
integer, dimension(:,:,:), allocatable, public :: phasememberAt !< memberID of given phase at every material point (ipc,ip,el)
integer, dimension(:,:,:), allocatable, public, target :: mappingHomogenization !< mapping from material points to offset in heterogenous state/field integer, dimension(:,:,:), allocatable, public, target :: mappingHomogenization !< mapping from material points to offset in heterogenous state/field
integer, dimension(:,:), allocatable, private, target :: mappingHomogenizationConst !< mapping from material points to offset in constant state/field integer, dimension(:,:), allocatable, private, target :: mappingHomogenizationConst !< mapping from material points to offset in constant state/field
! END DEPRECATED ! END DEPRECATED
@ -247,11 +244,7 @@ contains
subroutine material_init subroutine material_init
integer, parameter :: FILEUNIT = 210 integer, parameter :: FILEUNIT = 210
integer :: m,c,h, myDebug, myPhase, myHomog integer :: i,e,m,c,h, myDebug, myPhase, myHomog
integer :: &
g, & !< grain number
i, & !< integration point number
e !< element number
integer, dimension(:), allocatable :: & integer, dimension(:), allocatable :: &
CounterPhase, & CounterPhase, &
CounterHomogenization CounterHomogenization
@ -375,26 +368,15 @@ subroutine material_init
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! BEGIN DEPRECATED ! BEGIN DEPRECATED
allocate(phaseAt ( homogenization_maxNgrains,discretization_nIP,discretization_nElem),source=0)
allocate(phasememberAt ( homogenization_maxNgrains,discretization_nIP,discretization_nElem),source=0)
allocate(mappingHomogenization (2, discretization_nIP,discretization_nElem),source=0) allocate(mappingHomogenization (2, discretization_nIP,discretization_nElem),source=0)
allocate(mappingHomogenizationConst( discretization_nIP,discretization_nElem),source=1) allocate(mappingHomogenizationConst( discretization_nIP,discretization_nElem),source=1)
CounterHomogenization=0 CounterHomogenization=0
CounterPhase =0
do e = 1,discretization_nElem do e = 1,discretization_nElem
myHomog = discretization_homogenizationAt(e) myHomog = discretization_homogenizationAt(e)
do i = 1, discretization_nIP do i = 1, discretization_nIP
CounterHomogenization(myHomog) = CounterHomogenization(myHomog) + 1 CounterHomogenization(myHomog) = CounterHomogenization(myHomog) + 1
mappingHomogenization(1:2,i,e) = [CounterHomogenization(myHomog),huge(1)] mappingHomogenization(1:2,i,e) = [CounterHomogenization(myHomog),huge(1)]
do g = 1,homogenization_Ngrains(myHomog)
myPhase = material_phase(g,i,e)
CounterPhase(myPhase) = CounterPhase(myPhase)+1 ! not distinguishing between instances of same phase
phaseAt(g,i,e) = myPhase
phasememberAt(g,i,e) = CounterPhase(myPhase)
enddo
enddo enddo
enddo enddo
! END DEPRECATED ! END DEPRECATED

View File

@ -13,6 +13,7 @@ module plastic_disloUCLA
use material use material
use config use config
use lattice use lattice
use discretization
use results use results
implicit none implicit none
@ -295,7 +296,7 @@ subroutine plastic_disloUCLA_init()
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! allocate state arrays ! allocate state arrays
NipcMyPhase = count(material_phase == p) NipcMyPhase = count(material_phaseAt == p) * discretization_nIP
sizeDotState = size(['rho_mob ','rho_dip ','gamma_sl']) * prm%sum_N_sl sizeDotState = size(['rho_mob ','rho_dip ','gamma_sl']) * prm%sum_N_sl
sizeState = sizeDotState sizeState = sizeDotState

View File

@ -15,6 +15,7 @@ module plastic_dislotwin
use material use material
use config use config
use lattice use lattice
use discretization
#if defined(PETSc) || defined(DAMASK_HDF5) #if defined(PETSc) || defined(DAMASK_HDF5)
use results use results
#endif #endif
@ -494,7 +495,7 @@ subroutine plastic_dislotwin_init
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! allocate state arrays ! allocate state arrays
NipcMyPhase = count(material_phase == p) NipcMyPhase = count(material_phaseAt == p) * discretization_nIP
sizeDotState = size(['rho_mob ','rho_dip ','gamma_sl']) * prm%sum_N_sl & sizeDotState = size(['rho_mob ','rho_dip ','gamma_sl']) * prm%sum_N_sl &
+ size(['f_tw']) * prm%sum_N_tw & + size(['f_tw']) * prm%sum_N_tw &
+ size(['f_tr']) * prm%sum_N_tr + size(['f_tr']) * prm%sum_N_tr
@ -581,9 +582,9 @@ function plastic_dislotwin_homogenizedC(ipc,ip,el) result(homogenizedC)
of of
real(pReal) :: f_unrotated real(pReal) :: f_unrotated
of = phasememberAt(ipc,ip,el) of = material_phasememberAt(ipc,ip,el)
associate(prm => param(phase_plasticityInstance(material_phase(ipc,ip,el))),& associate(prm => param(phase_plasticityInstance(material_phaseAt(ipc,el))),&
stt => state(phase_plasticityInstance(material_phase(ipc,ip,el)))) stt => state(phase_plasticityInstance(material_phaseAT(ipc,el))))
f_unrotated = 1.0_pReal & f_unrotated = 1.0_pReal &
- sum(stt%f_tw(1:prm%sum_N_tw,of)) & - sum(stt%f_tw(1:prm%sum_N_tw,of)) &

View File

@ -14,9 +14,8 @@ module plastic_isotropic
use IO use IO
use material use material
use config use config
#if defined(PETSc) || defined(DAMASK_HDF5) use discretization
use results use results
#endif
implicit none implicit none
private private
@ -127,8 +126,8 @@ subroutine plastic_isotropic_init
config => config_phase(p)) config => config_phase(p))
#ifdef DEBUG #ifdef DEBUG
if (p==material_phase(debug_g,debug_i,debug_e)) & if (p==material_phaseAt(debug_g,debug_e)) &
prm%of_debug = phasememberAt(debug_g,debug_i,debug_e) prm%of_debug = material_phasememberAt(debug_g,debug_i,debug_e)
#endif #endif
prm%xi_0 = config%getFloat('tau0') prm%xi_0 = config%getFloat('tau0')
@ -190,7 +189,7 @@ subroutine plastic_isotropic_init
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! allocate state arrays ! allocate state arrays
NipcMyPhase = count(material_phase == p) NipcMyPhase = count(material_phaseAt == p) * discretization_nIP
sizeDotState = size(['xi ','accumulated_shear']) sizeDotState = size(['xi ','accumulated_shear'])
sizeState = sizeDotState sizeState = sizeDotState

View File

@ -13,9 +13,8 @@ module plastic_kinehardening
use material use material
use config use config
use lattice use lattice
#if defined(PETSc) || defined(DAMASK_HDF5) use discretization
use results use results
#endif
implicit none implicit none
private private
@ -146,8 +145,8 @@ subroutine plastic_kinehardening_init
config => config_phase(p)) config => config_phase(p))
#ifdef DEBUG #ifdef DEBUG
if (p==material_phase(debug_g,debug_i,debug_e)) then if (p==material_phaseAt(debug_g,debug_e)) then
prm%of_debug = phasememberAt(debug_g,debug_i,debug_e) prm%of_debug = material_phasememberAt(debug_g,debug_i,debug_e)
endif endif
#endif #endif
@ -257,7 +256,7 @@ subroutine plastic_kinehardening_init
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! allocate state arrays ! allocate state arrays
NipcMyPhase = count(material_phase == p) NipcMyPhase = count(material_phaseAt == p) * discretization_nIP
sizeDotState = size(['crss ','crss_back', 'accshear ']) * prm%totalNslip sizeDotState = size(['crss ','crss_back', 'accshear ']) * prm%totalNslip
sizeDeltaState = size(['sense ', 'chi0 ', 'gamma0' ]) * prm%totalNslip sizeDeltaState = size(['sense ', 'chi0 ', 'gamma0' ]) * prm%totalNslip
sizeState = sizeDotState + sizeDeltaState sizeState = sizeDotState + sizeDeltaState

View File

@ -6,6 +6,7 @@
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
module plastic_none module plastic_none
use material use material
use discretization
use debug use debug
implicit none implicit none
@ -36,7 +37,7 @@ subroutine plastic_none_init
do p = 1, size(phase_plasticity) do p = 1, size(phase_plasticity)
if (phase_plasticity(p) /= PLASTICITY_NONE_ID) cycle if (phase_plasticity(p) /= PLASTICITY_NONE_ID) cycle
NipcMyPhase = count(material_phase == p) NipcMyPhase = count(material_phaseAt == p) * discretization_nIP
call material_allocatePlasticState(p,NipcMyPhase,0,0,0, & call material_allocatePlasticState(p,NipcMyPhase,0,0,0, &
0,0,0) 0,0,0)
plasticState(p)%sizePostResults = 0 plasticState(p)%sizePostResults = 0

View File

@ -556,7 +556,7 @@ subroutine plastic_nonlocal_init
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! allocate state arrays ! allocate state arrays
NofMyPhase=count(material_phase==p) NofMyPhase = count(material_phaseAt==p) * discretization_nIP
sizeDotState = size([ 'rhoSglEdgePosMobile ','rhoSglEdgeNegMobile ', & sizeDotState = size([ 'rhoSglEdgePosMobile ','rhoSglEdgeNegMobile ', &
'rhoSglScrewPosMobile ','rhoSglScrewNegMobile ', & 'rhoSglScrewPosMobile ','rhoSglScrewNegMobile ', &
'rhoSglEdgePosImmobile ','rhoSglEdgeNegImmobile ', & 'rhoSglEdgePosImmobile ','rhoSglEdgeNegImmobile ', &
@ -677,7 +677,7 @@ subroutine plastic_nonlocal_init
allocate(iD(maxval(totalNslip),2,maxNinstances), source=0) allocate(iD(maxval(totalNslip),2,maxNinstances), source=0)
initializeInstances: do p = 1, size(phase_plasticity) initializeInstances: do p = 1, size(phase_plasticity)
NofMyPhase=count(material_phase==p) NofMyPhase = count(material_phaseAt==p) * discretization_nIP
myPhase2: if (phase_plasticity(p) == PLASTICITY_NONLOCAL_ID) then myPhase2: if (phase_plasticity(p) == PLASTICITY_NONLOCAL_ID) then
!*** determine indices to state array !*** determine indices to state array
@ -766,7 +766,7 @@ subroutine plastic_nonlocal_init
! get the total volume of the instance ! get the total volume of the instance
do e = 1,discretization_nElem do e = 1,discretization_nElem
do i = 1,discretization_nIP do i = 1,discretization_nIP
if (material_phase(1,i,e) == phase) volume(phasememberAt(1,i,e)) = IPvolume(i,e) if (material_phaseAt(1,e) == phase) volume(material_phasememberAt(1,i,e)) = IPvolume(i,e)
enddo enddo
enddo enddo
totalVolume = sum(volume) totalVolume = sum(volume)
@ -854,29 +854,29 @@ subroutine plastic_nonlocal_dependentState(Fe, Fp, ip, el)
invConnections invConnections
real(pReal), dimension(3,nIPneighbors) :: & real(pReal), dimension(3,nIPneighbors) :: &
connection_latticeConf connection_latticeConf
real(pReal), dimension(2,totalNslip(phase_plasticityInstance(material_phase(1,ip,el)))) :: & real(pReal), dimension(2,totalNslip(phase_plasticityInstance(material_phaseAt(1,el)))) :: &
rhoExcess rhoExcess
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el)))) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,el)))) :: &
rho_edg_delta, & rho_edg_delta, &
rho_scr_delta rho_scr_delta
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el))),10) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,el))),10) :: &
rho, & rho, &
rho_neighbor rho_neighbor
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el))), & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,el))), &
totalNslip(phase_plasticityInstance(material_phase(1,ip,el)))) :: & totalNslip(phase_plasticityInstance(material_phaseAt(1,el)))) :: &
myInteractionMatrix ! corrected slip interaction matrix myInteractionMatrix ! corrected slip interaction matrix
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el))),nIPneighbors) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,el))),nIPneighbors) :: &
rho_edg_delta_neighbor, & rho_edg_delta_neighbor, &
rho_scr_delta_neighbor rho_scr_delta_neighbor
real(pReal), dimension(2,maxval(totalNslip),nIPneighbors) :: & real(pReal), dimension(2,maxval(totalNslip),nIPneighbors) :: &
neighbor_rhoExcess, & ! excess density at neighboring material point neighbor_rhoExcess, & ! excess density at neighboring material point
neighbor_rhoTotal ! total density at neighboring material point neighbor_rhoTotal ! total density at neighboring material point
real(pReal), dimension(3,totalNslip(phase_plasticityInstance(material_phase(1,ip,el))),2) :: & real(pReal), dimension(3,totalNslip(phase_plasticityInstance(material_phaseAt(1,el))),2) :: &
m ! direction of dislocation motion m ! direction of dislocation motion
ph = phaseAt(1,ip,el) ph = material_phaseAt(1,el)
of = phasememberAt(1,ip,el) of = material_phasememberAt(1,ip,el)
instance = phase_plasticityInstance(ph) instance = phase_plasticityInstance(ph)
associate(prm => param(instance),dst => microstructure(instance), stt => state(instance)) associate(prm => param(instance),dst => microstructure(instance), stt => state(instance))
@ -935,9 +935,9 @@ subroutine plastic_nonlocal_dependentState(Fe, Fp, ip, el)
do n = 1,nIPneighbors do n = 1,nIPneighbors
neighbor_el = IPneighborhood(1,n,ip,el) neighbor_el = IPneighborhood(1,n,ip,el)
neighbor_ip = IPneighborhood(2,n,ip,el) neighbor_ip = IPneighborhood(2,n,ip,el)
no = phasememberAt(1,neighbor_ip,neighbor_el) no = material_phasememberAt(1,neighbor_ip,neighbor_el)
if (neighbor_el > 0 .and. neighbor_ip > 0) then if (neighbor_el > 0 .and. neighbor_ip > 0) then
neighbor_instance = phase_plasticityInstance(material_phase(1,neighbor_ip,neighbor_el)) neighbor_instance = phase_plasticityInstance(material_phaseAt(1,neighbor_el))
if (neighbor_instance == instance) then if (neighbor_instance == instance) then
nRealNeighbors = nRealNeighbors + 1.0_pReal nRealNeighbors = nRealNeighbors + 1.0_pReal
@ -1202,22 +1202,22 @@ subroutine plastic_nonlocal_LpAndItsTangent(Lp, dLp_dMp, &
of, & !offset of, & !offset
t, & !< dislocation type t, & !< dislocation type
s !< index of my current slip system s !< index of my current slip system
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el))),8) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,el))),8) :: &
rhoSgl !< single dislocation densities (including blocked) rhoSgl !< single dislocation densities (including blocked)
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el))),10) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,el))),10) :: &
rho rho
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el))),4) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,el))),4) :: &
v, & !< velocity v, & !< velocity
tauNS, & !< resolved shear stress including non Schmid and backstress terms tauNS, & !< resolved shear stress including non Schmid and backstress terms
dv_dtau, & !< velocity derivative with respect to the shear stress dv_dtau, & !< velocity derivative with respect to the shear stress
dv_dtauNS !< velocity derivative with respect to the shear stress dv_dtauNS !< velocity derivative with respect to the shear stress
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el)))) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,el)))) :: &
tau, & !< resolved shear stress including backstress terms tau, & !< resolved shear stress including backstress terms
gdotTotal !< shear rate gdotTotal !< shear rate
!*** shortcut for mapping !*** shortcut for mapping
ph = phaseAt(1,ip,el) ph = material_phaseAt(1,el)
of = phasememberAt(1,ip,el) of = material_phasememberAt(1,ip,el)
instance = phase_plasticityInstance(ph) instance = phase_plasticityInstance(ph)
associate(prm => param(instance),dst=>microstructure(instance)) associate(prm => param(instance),dst=>microstructure(instance))
@ -1323,23 +1323,23 @@ subroutine plastic_nonlocal_deltaState(Mp,ip,el)
c, & ! character of dislocation c, & ! character of dislocation
t, & ! type of dislocation t, & ! type of dislocation
s ! index of my current slip system s ! index of my current slip system
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el))),10) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,el))),10) :: &
deltaRhoRemobilization, & ! density increment by remobilization deltaRhoRemobilization, & ! density increment by remobilization
deltaRhoDipole2SingleStress ! density increment by dipole dissociation (by stress change) deltaRhoDipole2SingleStress ! density increment by dipole dissociation (by stress change)
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el))),10) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,el))),10) :: &
rho ! current dislocation densities rho ! current dislocation densities
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el))),4) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,el))),4) :: &
v ! dislocation glide velocity v ! dislocation glide velocity
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el)))) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,el)))) :: &
tau ! current resolved shear stress tau ! current resolved shear stress
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el))),2) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,el))),2) :: &
rhoDip, & ! current dipole dislocation densities (screw and edge dipoles) rhoDip, & ! current dipole dislocation densities (screw and edge dipoles)
dUpper, & ! current maximum stable dipole distance for edges and screws dUpper, & ! current maximum stable dipole distance for edges and screws
dUpperOld, & ! old maximum stable dipole distance for edges and screws dUpperOld, & ! old maximum stable dipole distance for edges and screws
deltaDUpper ! change in maximum stable dipole distance for edges and screws deltaDUpper ! change in maximum stable dipole distance for edges and screws
ph = phaseAt(1,ip,el) ph = material_phaseAt(1,el)
of = phasememberAt(1,ip,el) of = material_phasememberAt(1,ip,el)
instance = phase_plasticityInstance(ph) instance = phase_plasticityInstance(ph)
associate(prm => param(instance),dst => microstructure(instance),del => deltaState(instance)) associate(prm => param(instance),dst => microstructure(instance),del => deltaState(instance))
@ -1459,7 +1459,7 @@ subroutine plastic_nonlocal_dotState(Mp, Fe, Fp, Temperature, &
np,& !< neighbour phase shortcut np,& !< neighbour phase shortcut
topp, & !< type of dislocation with opposite sign to t topp, & !< type of dislocation with opposite sign to t
s !< index of my current slip system s !< index of my current slip system
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el))),10) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,el))),10) :: &
rho, & rho, &
rhoDot, & !< density evolution rhoDot, & !< density evolution
rhoDotMultiplication, & !< density evolution by multiplication rhoDotMultiplication, & !< density evolution by multiplication
@ -1467,24 +1467,24 @@ subroutine plastic_nonlocal_dotState(Mp, Fe, Fp, Temperature, &
rhoDotSingle2DipoleGlide, & !< density evolution by dipole formation (by glide) rhoDotSingle2DipoleGlide, & !< density evolution by dipole formation (by glide)
rhoDotAthermalAnnihilation, & !< density evolution by athermal annihilation rhoDotAthermalAnnihilation, & !< density evolution by athermal annihilation
rhoDotThermalAnnihilation !< density evolution by thermal annihilation rhoDotThermalAnnihilation !< density evolution by thermal annihilation
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el))),8) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,el))),8) :: &
rhoSgl, & !< current single dislocation densities (positive/negative screw and edge without dipoles) rhoSgl, & !< current single dislocation densities (positive/negative screw and edge without dipoles)
neighbor_rhoSgl, & !< current single dislocation densities of neighboring ip (positive/negative screw and edge without dipoles) neighbor_rhoSgl, & !< current single dislocation densities of neighboring ip (positive/negative screw and edge without dipoles)
my_rhoSgl !< single dislocation densities of central ip (positive/negative screw and edge without dipoles) my_rhoSgl !< single dislocation densities of central ip (positive/negative screw and edge without dipoles)
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el))),4) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,el))),4) :: &
v, & !< current dislocation glide velocity v, & !< current dislocation glide velocity
my_v, & !< dislocation glide velocity of central ip my_v, & !< dislocation glide velocity of central ip
neighbor_v, & !< dislocation glide velocity of enighboring ip neighbor_v, & !< dislocation glide velocity of enighboring ip
gdot !< shear rates gdot !< shear rates
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el)))) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,el)))) :: &
tau, & !< current resolved shear stress tau, & !< current resolved shear stress
vClimb !< climb velocity of edge dipoles vClimb !< climb velocity of edge dipoles
real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phase(1,ip,el))),2) :: & real(pReal), dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,el))),2) :: &
rhoDip, & !< current dipole dislocation densities (screw and edge dipoles) rhoDip, & !< current dipole dislocation densities (screw and edge dipoles)
dLower, & !< minimum stable dipole distance for edges and screws dLower, & !< minimum stable dipole distance for edges and screws
dUpper !< current maximum stable dipole distance for edges and screws dUpper !< current maximum stable dipole distance for edges and screws
real(pReal), dimension(3,totalNslip(phase_plasticityInstance(material_phase(1,ip,el))),4) :: & real(pReal), dimension(3,totalNslip(phase_plasticityInstance(material_phaseAt(1,el))),4) :: &
m !< direction of dislocation motion m !< direction of dislocation motion
real(pReal), dimension(3,3) :: & real(pReal), dimension(3,3) :: &
my_F, & !< my total deformation gradient my_F, & !< my total deformation gradient
@ -1507,15 +1507,15 @@ subroutine plastic_nonlocal_dotState(Mp, Fe, Fp, Temperature, &
considerEnteringFlux, & considerEnteringFlux, &
considerLeavingFlux considerLeavingFlux
p = phaseAt(1,ip,el) p = material_phaseAt(1,el)
o = phasememberAt(1,ip,el) o = material_phasememberAt(1,ip,el)
if (timestep <= 0.0_pReal) then if (timestep <= 0.0_pReal) then
plasticState(p)%dotState = 0.0_pReal plasticState(p)%dotState = 0.0_pReal
return return
endif endif
ph = material_phase(1,ip,el) ph = material_phaseAt(1,el)
instance = phase_plasticityInstance(ph) instance = phase_plasticityInstance(ph)
associate(prm => param(instance),dst => microstructure(instance),dot => dotState(instance),stt => state(instance)) associate(prm => param(instance),dst => microstructure(instance),dot => dotState(instance),stt => state(instance))
ns = totalNslip(instance) ns = totalNslip(instance)
@ -1592,7 +1592,7 @@ subroutine plastic_nonlocal_dotState(Mp, Fe, Fp, Temperature, &
!**************************************************************************** !****************************************************************************
!*** calculate dislocation fluxes (only for nonlocal plasticity) !*** calculate dislocation fluxes (only for nonlocal plasticity)
rhoDotFlux = 0.0_pReal rhoDotFlux = 0.0_pReal
if (.not. phase_localPlasticity(material_phase(1,ip,el))) then if (.not. phase_localPlasticity(material_phaseAt(1,el))) then
!*** check CFL (Courant-Friedrichs-Lewy) condition for flux !*** check CFL (Courant-Friedrichs-Lewy) condition for flux
if (any( abs(gdot) > 0.0_pReal & ! any active slip system ... if (any( abs(gdot) > 0.0_pReal & ! any active slip system ...
@ -1630,8 +1630,8 @@ subroutine plastic_nonlocal_dotState(Mp, Fe, Fp, Temperature, &
neighbor_el = IPneighborhood(1,n,ip,el) neighbor_el = IPneighborhood(1,n,ip,el)
neighbor_ip = IPneighborhood(2,n,ip,el) neighbor_ip = IPneighborhood(2,n,ip,el)
neighbor_n = IPneighborhood(3,n,ip,el) neighbor_n = IPneighborhood(3,n,ip,el)
np = phaseAt(1,neighbor_ip,neighbor_el) np = material_phaseAt(1,neighbor_el)
no = phasememberAt(1,neighbor_ip,neighbor_el) no = material_phasememberAt(1,neighbor_ip,neighbor_el)
opposite_neighbor = n + mod(n,2) - mod(n+1,2) opposite_neighbor = n + mod(n,2) - mod(n+1,2)
opposite_el = IPneighborhood(1,opposite_neighbor,ip,el) opposite_el = IPneighborhood(1,opposite_neighbor,ip,el)
@ -1639,7 +1639,7 @@ subroutine plastic_nonlocal_dotState(Mp, Fe, Fp, Temperature, &
opposite_n = IPneighborhood(3,opposite_neighbor,ip,el) opposite_n = IPneighborhood(3,opposite_neighbor,ip,el)
if (neighbor_n > 0) then ! if neighbor exists, average deformation gradient if (neighbor_n > 0) then ! if neighbor exists, average deformation gradient
neighbor_instance = phase_plasticityInstance(material_phase(1,neighbor_ip,neighbor_el)) neighbor_instance = phase_plasticityInstance(material_phaseAt(1,neighbor_el))
neighbor_Fe = Fe(1:3,1:3,1,neighbor_ip,neighbor_el) neighbor_Fe = Fe(1:3,1:3,1,neighbor_ip,neighbor_el)
neighbor_F = matmul(neighbor_Fe, Fp(1:3,1:3,1,neighbor_ip,neighbor_el)) neighbor_F = matmul(neighbor_Fe, Fp(1:3,1:3,1,neighbor_ip,neighbor_el))
Favg = 0.5_pReal * (my_F + neighbor_F) Favg = 0.5_pReal * (my_F + neighbor_F)
@ -1661,7 +1661,7 @@ subroutine plastic_nonlocal_dotState(Mp, Fe, Fp, Temperature, &
neighbor_v = 0.0_pReal ! needed for check of sign change in flux density below neighbor_v = 0.0_pReal ! needed for check of sign change in flux density below
neighbor_rhoSgl = 0.0_pReal neighbor_rhoSgl = 0.0_pReal
if (neighbor_n > 0) then if (neighbor_n > 0) then
if (phase_plasticity(material_phase(1,neighbor_ip,neighbor_el)) == PLASTICITY_NONLOCAL_ID & if (phase_plasticity(material_phaseAt(1,neighbor_el)) == PLASTICITY_NONLOCAL_ID &
.and. any(compatibility(:,:,:,n,ip,el) > 0.0_pReal)) & .and. any(compatibility(:,:,:,n,ip,el) > 0.0_pReal)) &
considerEnteringFlux = .true. considerEnteringFlux = .true.
endif endif
@ -1714,7 +1714,7 @@ subroutine plastic_nonlocal_dotState(Mp, Fe, Fp, Temperature, &
considerLeavingFlux = .true. considerLeavingFlux = .true.
if (opposite_n > 0) then if (opposite_n > 0) then
if (phase_plasticity(material_phase(1,opposite_ip,opposite_el)) /= PLASTICITY_NONLOCAL_ID) & if (phase_plasticity(material_phaseAt(1,opposite_el)) /= PLASTICITY_NONLOCAL_ID) &
considerLeavingFlux = .false. considerLeavingFlux = .false.
endif endif
@ -1905,20 +1905,20 @@ subroutine plastic_nonlocal_updateCompatibility(orientation,i,e)
s2 ! slip system index (my neighbor) s2 ! slip system index (my neighbor)
real(pReal), dimension(4) :: & real(pReal), dimension(4) :: &
absoluteMisorientation ! absolute misorientation (without symmetry) between me and my neighbor absoluteMisorientation ! absolute misorientation (without symmetry) between me and my neighbor
real(pReal), dimension(2,totalNslip(phase_plasticityInstance(material_phase(1,i,e))),& real(pReal), dimension(2,totalNslip(phase_plasticityInstance(material_phaseAt(1,e))),&
totalNslip(phase_plasticityInstance(material_phase(1,i,e))),& totalNslip(phase_plasticityInstance(material_phaseAt(1,e))),&
nIPneighbors) :: & nIPneighbors) :: &
my_compatibility ! my_compatibility for current element and ip my_compatibility ! my_compatibility for current element and ip
real(pReal) :: & real(pReal) :: &
my_compatibilitySum, & my_compatibilitySum, &
thresholdValue, & thresholdValue, &
nThresholdValues nThresholdValues
logical, dimension(totalNslip(phase_plasticityInstance(material_phase(1,i,e)))) :: & logical, dimension(totalNslip(phase_plasticityInstance(material_phaseAt(1,e)))) :: &
belowThreshold belowThreshold
type(rotation) :: rot type(rotation) :: rot
Nneighbors = nIPneighbors Nneighbors = nIPneighbors
ph = material_phase(1,i,e) ph = material_phaseAt(1,e)
textureID = material_texture(1,i,e) textureID = material_texture(1,i,e)
instance = phase_plasticityInstance(ph) instance = phase_plasticityInstance(ph)
ns = totalNslip(instance) ns = totalNslip(instance)
@ -1950,7 +1950,7 @@ subroutine plastic_nonlocal_updateCompatibility(orientation,i,e)
!* we consider this to be a real "physical" phase boundary, so completely incompatible. !* we consider this to be a real "physical" phase boundary, so completely incompatible.
!* If one of the two phases has a local plasticity law, !* If one of the two phases has a local plasticity law,
!* we do not consider this to be a phase boundary, so completely compatible. !* we do not consider this to be a phase boundary, so completely compatible.
neighbor_phase = material_phase(1,neighbor_i,neighbor_e) neighbor_phase = material_phaseAt(1,neighbor_e)
if (neighbor_phase /= ph) then if (neighbor_phase /= ph) then
if (.not. phase_localPlasticity(neighbor_phase) .and. .not. phase_localPlasticity(ph))& if (.not. phase_localPlasticity(neighbor_phase) .and. .not. phase_localPlasticity(ph))&
forall(s1 = 1:ns) my_compatibility(1:2,s1,s1,n) = 0.0_pReal forall(s1 = 1:ns) my_compatibility(1:2,s1,s1,n) = 0.0_pReal

View File

@ -12,9 +12,8 @@ module plastic_phenopowerlaw
use material use material
use config use config
use lattice use lattice
#if defined(PETSc) || defined(DAMASK_HDF5) use discretization
use results use results
#endif
implicit none implicit none
private private
@ -314,7 +313,7 @@ subroutine plastic_phenopowerlaw_init
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
! allocate state arrays ! allocate state arrays
NipcMyPhase = count(material_phase == p) NipcMyPhase = count(material_phaseAt == p) * discretization_nIP
sizeDotState = size(['tau_slip ','gamma_slip']) * prm%totalNslip & sizeDotState = size(['tau_slip ','gamma_slip']) * prm%totalNslip &
+ size(['tau_twin ','gamma_twin']) * prm%totalNtwin + size(['tau_twin ','gamma_twin']) * prm%totalNtwin
sizeState = sizeDotState sizeState = sizeDotState

View File

@ -10,6 +10,7 @@ module source_damage_anisoBrittle
use IO use IO
use math use math
use material use material
use discretization
use config use config
use lattice use lattice
@ -164,7 +165,7 @@ subroutine source_damage_anisoBrittle_init
end associate end associate
phase = p phase = p
NofMyPhase=count(material_phase==phase) NofMyPhase=count(material_phaseAt==phase) * discretization_nIP
instance = source_damage_anisoBrittle_instance(phase) instance = source_damage_anisoBrittle_instance(phase)
sourceOffset = source_damage_anisoBrittle_offset(phase) sourceOffset = source_damage_anisoBrittle_offset(phase)
@ -202,8 +203,8 @@ subroutine source_damage_anisoBrittle_dotState(S, ipc, ip, el)
real(pReal) :: & real(pReal) :: &
traction_d, traction_t, traction_n, traction_crit traction_d, traction_t, traction_n, traction_crit
phase = phaseAt(ipc,ip,el) phase = material_phaseAt(ipc,el)
constituent = phasememberAt(ipc,ip,el) constituent = material_phasememberAt(ipc,ip,el)
instance = source_damage_anisoBrittle_instance(phase) instance = source_damage_anisoBrittle_instance(phase)
sourceOffset = source_damage_anisoBrittle_offset(phase) sourceOffset = source_damage_anisoBrittle_offset(phase)
homog = material_homogenizationAt(el) homog = material_homogenizationAt(el)

View File

@ -9,6 +9,7 @@ module source_damage_anisoDuctile
use debug use debug
use IO use IO
use math use math
use discretization
use material use material
use config use config
@ -150,7 +151,7 @@ subroutine source_damage_anisoDuctile_init
phase = p phase = p
NofMyPhase=count(material_phase==phase) NofMyPhase=count(material_phaseAt==phase) * discretization_nIP
instance = source_damage_anisoDuctile_instance(phase) instance = source_damage_anisoDuctile_instance(phase)
sourceOffset = source_damage_anisoDuctile_offset(phase) sourceOffset = source_damage_anisoDuctile_offset(phase)
@ -180,8 +181,8 @@ subroutine source_damage_anisoDuctile_dotState(ipc, ip, el)
instance, & instance, &
i i
phase = phaseAt(ipc,ip,el) phase = material_phaseAt(ipc,el)
constituent = phasememberAt(ipc,ip,el) constituent = material_phasememberAt(ipc,ip,el)
instance = source_damage_anisoDuctile_instance(phase) instance = source_damage_anisoDuctile_instance(phase)
sourceOffset = source_damage_anisoDuctile_offset(phase) sourceOffset = source_damage_anisoDuctile_offset(phase)
homog = material_homogenizationAt(el) homog = material_homogenizationAt(el)

View File

@ -9,6 +9,7 @@ module source_damage_isoBrittle
use debug use debug
use IO use IO
use math use math
use discretization
use material use material
use config use config
@ -133,7 +134,7 @@ subroutine source_damage_isoBrittle_init
phase = p phase = p
NofMyPhase=count(material_phase==phase) NofMyPhase = count(material_phaseAt==phase) * discretization_nIP
instance = source_damage_isoBrittle_instance(phase) instance = source_damage_isoBrittle_instance(phase)
sourceOffset = source_damage_isoBrittle_offset(phase) sourceOffset = source_damage_isoBrittle_offset(phase)
@ -164,8 +165,8 @@ subroutine source_damage_isoBrittle_deltaState(C, Fe, ipc, ip, el)
strain(6), & strain(6), &
strainenergy strainenergy
phase = phaseAt(ipc,ip,el) !< phase ID at ipc,ip,el phase = material_phaseAt(ipc,el) !< phase ID at ipc,ip,el
constituent = phasememberAt(ipc,ip,el) !< state array offset for phase ID at ipc,ip,el constituent = material_phasememberAt(ipc,ip,el) !< state array offset for phase ID at ipc,ip,el
! ToDo: capability for multiple instances of SAME source within given phase. Needs Ninstance loop from here on! ! ToDo: capability for multiple instances of SAME source within given phase. Needs Ninstance loop from here on!
instance = source_damage_isoBrittle_instance(phase) !< instance of damage_isoBrittle source instance = source_damage_isoBrittle_instance(phase) !< instance of damage_isoBrittle source
sourceOffset = source_damage_isoBrittle_offset(phase) sourceOffset = source_damage_isoBrittle_offset(phase)

View File

@ -8,6 +8,7 @@ module source_damage_isoDuctile
use prec use prec
use debug use debug
use IO use IO
use discretization
use material use material
use config use config
@ -132,7 +133,7 @@ subroutine source_damage_isoDuctile_init
end associate end associate
phase = p phase = p
NofMyPhase=count(material_phase==phase) NofMyPhase=count(material_phaseAt==phase) * discretization_nIP
instance = source_damage_isoDuctile_instance(phase) instance = source_damage_isoDuctile_instance(phase)
sourceOffset = source_damage_isoDuctile_offset(phase) sourceOffset = source_damage_isoDuctile_offset(phase)
@ -157,8 +158,8 @@ subroutine source_damage_isoDuctile_dotState(ipc, ip, el)
integer :: & integer :: &
phase, constituent, instance, homog, sourceOffset, damageOffset phase, constituent, instance, homog, sourceOffset, damageOffset
phase = phaseAt(ipc,ip,el) phase = material_phaseAt(ipc,el)
constituent = phasememberAt(ipc,ip,el) constituent = material_phasememberAt(ipc,ip,el)
instance = source_damage_isoDuctile_instance(phase) instance = source_damage_isoDuctile_instance(phase)
sourceOffset = source_damage_isoDuctile_offset(phase) sourceOffset = source_damage_isoDuctile_offset(phase)
homog = material_homogenizationAt(el) homog = material_homogenizationAt(el)

View File

@ -7,6 +7,7 @@
module source_thermal_dissipation module source_thermal_dissipation
use prec use prec
use debug use debug
use discretization
use material use material
use config use config
@ -75,7 +76,7 @@ subroutine source_thermal_dissipation_init
if (all(phase_source(:,p) /= SOURCE_THERMAL_DISSIPATION_ID)) cycle if (all(phase_source(:,p) /= SOURCE_THERMAL_DISSIPATION_ID)) cycle
instance = source_thermal_dissipation_instance(p) instance = source_thermal_dissipation_instance(p)
param(instance)%kappa = config_phase(p)%getFloat('dissipation_coldworkcoeff') param(instance)%kappa = config_phase(p)%getFloat('dissipation_coldworkcoeff')
NofMyPhase=count(material_phase==p) NofMyPhase = count(material_phaseAt==p) * discretization_nIP
sourceOffset = source_thermal_dissipation_offset(p) sourceOffset = source_thermal_dissipation_offset(p)
call material_allocateSourceState(p,sourceOffset,NofMyPhase,0,0,0) call material_allocateSourceState(p,sourceOffset,NofMyPhase,0,0,0)

View File

@ -7,6 +7,7 @@
module source_thermal_externalheat module source_thermal_externalheat
use prec use prec
use debug use debug
use discretization
use material use material
use config use config
@ -83,7 +84,7 @@ subroutine source_thermal_externalheat_init
if (all(phase_source(:,p) /= SOURCE_thermal_externalheat_ID)) cycle if (all(phase_source(:,p) /= SOURCE_thermal_externalheat_ID)) cycle
instance = source_thermal_externalheat_instance(p) instance = source_thermal_externalheat_instance(p)
sourceOffset = source_thermal_externalheat_offset(p) sourceOffset = source_thermal_externalheat_offset(p)
NofMyPhase=count(material_phase==p) NofMyPhase = count(material_phaseAt==p) * discretization_nIP
param(instance)%time = config_phase(p)%getFloats('externalheat_time') param(instance)%time = config_phase(p)%getFloats('externalheat_time')
param(instance)%nIntervals = size(param(instance)%time) - 1 param(instance)%nIntervals = size(param(instance)%time) - 1

View File

@ -167,8 +167,8 @@ subroutine thermal_adiabatic_getSourceAndItsTangent(Tdot, dTdot_dT, T, ip, el)
Tdot = 0.0_pReal Tdot = 0.0_pReal
dTdot_dT = 0.0_pReal dTdot_dT = 0.0_pReal
do grain = 1, homogenization_Ngrains(homog) do grain = 1, homogenization_Ngrains(homog)
phase = phaseAt(grain,ip,el) phase = material_phaseAt(grain,el)
constituent = phasememberAt(grain,ip,el) constituent = material_phasememberAt(grain,ip,el)
do source = 1, phase_Nsources(phase) do source = 1, phase_Nsources(phase)
select case(phase_source(source,phase)) select case(phase_source(source,phase))
case (SOURCE_thermal_dissipation_ID) case (SOURCE_thermal_dissipation_ID)

View File

@ -132,8 +132,8 @@ subroutine thermal_conduction_getSourceAndItsTangent(Tdot, dTdot_dT, T, ip, el)
Tdot = 0.0_pReal Tdot = 0.0_pReal
dTdot_dT = 0.0_pReal dTdot_dT = 0.0_pReal
do grain = 1, homogenization_Ngrains(homog) do grain = 1, homogenization_Ngrains(homog)
phase = phaseAt(grain,ip,el) phase = material_phaseAt(grain,el)
constituent = phasememberAt(grain,ip,el) constituent = material_phasememberAt(grain,ip,el)
do source = 1, phase_Nsources(phase) do source = 1, phase_Nsources(phase)
select case(phase_source(source,phase)) select case(phase_source(source,phase))
case (SOURCE_thermal_dissipation_ID) case (SOURCE_thermal_dissipation_ID)