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

polishing/unifying

This commit is contained in:
Martin Diehl 2020-02-28 10:25:07 +01:00
parent cf0f5f0fee
commit 991d0fe020
9 changed files with 131 additions and 138 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

35
src/kinematics_cleavage_opening.f90
View File

@ -63,40 +63,41 @@ subroutine kinematics_cleavage_opening_init
integer, allocatable, dimension(:) :: tempInt
real(pReal), allocatable, dimension(:) :: tempFloat
integer :: maxNinstance,p,instance
integer :: Ninstance,p,instance
write(6,'(/,a)') ' <<<+- kinematics_'//KINEMATICS_cleavage_opening_LABEL//' init -+>>>'; flush(6)
maxNinstance = count(phase_kinematics == KINEMATICS_cleavage_opening_ID)
if (maxNinstance == 0) return
Ninstance = count(phase_kinematics == KINEMATICS_cleavage_opening_ID)
if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) &
write(6,'(a16,1x,i5,/)') '# instances:',maxNinstance
write(6,'(a16,1x,i5,/)') '# instances:',Ninstance
allocate(kinematics_cleavage_opening_instance(size(config_phase)), source=0)
do p = 1, size(config_phase)
kinematics_cleavage_opening_instance(p) = count(phase_kinematics(:,1:p) == kinematics_cleavage_opening_ID) ! ToDo: count correct?
enddo
allocate(kinematics_cleavage_opening_critDisp(lattice_maxNcleavageFamily,maxNinstance), source=0.0_pReal)
allocate(kinematics_cleavage_opening_critLoad(lattice_maxNcleavageFamily,maxNinstance), source=0.0_pReal)
allocate(kinematics_cleavage_opening_Ncleavage(lattice_maxNcleavageFamily,maxNinstance), source=0)
allocate(kinematics_cleavage_opening_totalNcleavage(maxNinstance), source=0)
allocate(kinematics_cleavage_opening_sdot_0(maxNinstance), source=0.0_pReal)
allocate(kinematics_cleavage_opening_N(maxNinstance), source=0.0_pReal)
allocate(kinematics_cleavage_opening_critDisp(lattice_maxNcleavageFamily,Ninstance), source=0.0_pReal)
allocate(kinematics_cleavage_opening_critLoad(lattice_maxNcleavageFamily,Ninstance), source=0.0_pReal)
allocate(kinematics_cleavage_opening_Ncleavage(lattice_maxNcleavageFamily,Ninstance), source=0)
allocate(kinematics_cleavage_opening_totalNcleavage(Ninstance), source=0)
allocate(kinematics_cleavage_opening_sdot_0(Ninstance), source=0.0_pReal)
allocate(kinematics_cleavage_opening_N(Ninstance), source=0.0_pReal)
do p = 1, size(config_phase)
if (all(phase_kinematics(:,p) /= KINEMATICS_cleavage_opening_ID)) cycle
associate(config => config_phase(p))
instance = kinematics_cleavage_opening_instance(p)
kinematics_cleavage_opening_sdot_0(instance) = config_phase(p)%getFloat('anisobrittle_sdot0')
kinematics_cleavage_opening_N(instance) = config_phase(p)%getFloat('anisobrittle_ratesensitivity')
tempInt = config_phase(p)%getInts('ncleavage')
kinematics_cleavage_opening_sdot_0(instance) = config%getFloat('anisobrittle_sdot0')
kinematics_cleavage_opening_N(instance) = config%getFloat('anisobrittle_ratesensitivity')
tempInt = config%getInts('ncleavage')
kinematics_cleavage_opening_Ncleavage(1:size(tempInt),instance) = tempInt
tempFloat = config_phase(p)%getFloats('anisobrittle_criticaldisplacement',requiredSize=size(tempInt))
tempFloat = config%getFloats('anisobrittle_criticaldisplacement',requiredSize=size(tempInt))
kinematics_cleavage_opening_critDisp(1:size(tempInt),instance) = tempFloat
tempFloat = config_phase(p)%getFloats('anisobrittle_criticalload',requiredSize=size(tempInt))
tempFloat = config%getFloats('anisobrittle_criticalload',requiredSize=size(tempInt))
kinematics_cleavage_opening_critLoad(1:size(tempInt),instance) = tempFloat
kinematics_cleavage_opening_Ncleavage(1:lattice_maxNcleavageFamily,instance) = &
@ -111,6 +112,8 @@ subroutine kinematics_cleavage_opening_init
call IO_error(211,el=instance,ext_msg='critical_load ('//KINEMATICS_cleavage_opening_LABEL//')')
if (kinematics_cleavage_opening_N(instance) <= 0.0_pReal) &
call IO_error(211,el=instance,ext_msg='rate_sensitivity ('//KINEMATICS_cleavage_opening_LABEL//')')
end associate
enddo
end subroutine kinematics_cleavage_opening_init

34
src/kinematics_slipplane_opening.f90
View File

@ -49,37 +49,30 @@ contains
!--------------------------------------------------------------------------------------------------
subroutine kinematics_slipplane_opening_init
integer :: maxNinstance,p,instance
integer :: Ninstance,p,instance
write(6,'(/,a)') ' <<<+- kinematics_'//KINEMATICS_slipplane_opening_LABEL//' init -+>>>'; flush(6)
maxNinstance = count(phase_kinematics == KINEMATICS_slipplane_opening_ID)
if (maxNinstance == 0) return
Ninstance = count(phase_kinematics == KINEMATICS_slipplane_opening_ID)
if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) &
write(6,'(a16,1x,i5,/)') '# instances:',maxNinstance
write(6,'(a16,1x,i5,/)') '# instances:',Ninstance
allocate(kinematics_slipplane_opening_instance(size(config_phase)), source=0)
do p = 1, size(config_phase)
kinematics_slipplane_opening_instance(p) = count(phase_kinematics(:,1:p) == kinematics_slipplane_opening_ID) ! ToDo: count correct?
enddo
allocate(param(maxNinstance))
allocate(param(Ninstance))
do p = 1, size(config_phase)
if (all(phase_kinematics(:,p) /= KINEMATICS_slipplane_opening_ID)) cycle
associate(prm => param(kinematics_slipplane_opening_instance(p)), &
config => config_phase(p))
instance = kinematics_slipplane_opening_instance(p)
prm%sdot0 = config_phase(p)%getFloat('anisoductile_sdot0')
prm%n = config_phase(p)%getFloat('anisoductile_ratesensitivity')
prm%sdot0 = config%getFloat('anisoductile_sdot0')
prm%n = config%getFloat('anisoductile_ratesensitivity')
prm%Nslip = config%getInts('nslip')
prm%critLoad = config_phase(p)%getFloats('anisoductile_criticalload',requiredSize=size(prm%Nslip ))
prm%critLoad = math_expand(prm%critLoad, prm%Nslip)
prm%slip_direction = lattice_slip_direction (prm%Nslip,config%getString('lattice_structure'),&
config%getFloat('c/a',defaultVal=0.0_pReal))
prm%slip_normal = lattice_slip_normal (prm%Nslip,config%getString('lattice_structure'),&
@ -87,6 +80,9 @@ subroutine kinematics_slipplane_opening_init
prm%slip_transverse = lattice_slip_transverse(prm%Nslip,config%getString('lattice_structure'),&
config%getFloat('c/a',defaultVal=0.0_pReal))
prm%critLoad = config%getFloats('anisoductile_criticalload',requiredSize=size(prm%Nslip))
prm%critLoad = math_expand(prm%critLoad, prm%Nslip)
! if (kinematics_slipplane_opening_sdot_0(instance) <= 0.0_pReal) &
! call IO_error(211,el=instance,ext_msg='sdot_0 ('//KINEMATICS_slipplane_opening_LABEL//')')
! if (any(kinematics_slipplane_opening_critPlasticStrain(:,instance) < 0.0_pReal)) &
@ -151,9 +147,21 @@ subroutine kinematics_slipplane_opening_LiAndItsTangent(Ld, dLd_dTstar, S, ipc,
udotn = prm%sdot0* ( max(0.0_pReal,traction_n)/traction_crit &
- max(0.0_pReal,traction_n)/prm%critLoad(i))**prm%n
if (dNeq0(traction_d)) then
dudotd_dt = udotd*prm%n/traction_d
else
dudotd_dt = 0.0_pReal
endif
if (dNeq0(traction_t)) then
dudott_dt = udott*prm%n/traction_t
else
dudott_dt = 0.0_pReal
endif
if (dNeq0(traction_n)) then
dudotn_dt = udotn*prm%n/traction_n
else
dudotn_dt = 0.0_pReal
endif
forall (k=1:3,l=1:3,m=1:3,n=1:3) &
dLd_dTstar(k,l,m,n) = dLd_dTstar(k,l,m,n) + dudotd_dt*projection_d(k,l)*projection_d(m,n) &

16
src/kinematics_thermal_expansion.f90
View File

@ -36,35 +36,34 @@ contains
!--------------------------------------------------------------------------------------------------
subroutine kinematics_thermal_expansion_init
integer :: &
Ninstance, &
p, i
integer :: Ninstance,p,i,instance
real(pReal), dimension(:), allocatable :: &
temp
write(6,'(/,a)') ' <<<+- kinematics_'//KINEMATICS_thermal_expansion_LABEL//' init -+>>>'; flush(6)
Ninstance = count(phase_kinematics == KINEMATICS_thermal_expansion_ID)
if (iand(debug_level(debug_constitutive),debug_levelBasic) /= 0) &
write(6,'(a16,1x,i5,/)') '# instances:',Ninstance
allocate(param(Ninstance))
do p = 1, size(phase_kinematics)
do p = 1, size(config_phase)
if (all(phase_kinematics(:,p) /= KINEMATICS_thermal_expansion_ID)) cycle
associate(config => config_phase(p))
! ToDo: Here we need to decide how to extend the concept of instances to
! kinetics and sources. I would suggest that the same mechanism exists at maximum once per phase
! read up to three parameters (constant, linear, quadratic with T)
temp = config_phase(p)%getFloats('thermal_expansion11')
temp = config%getFloats('thermal_expansion11')
!lattice_thermalExpansion33(1,1,1:size(temp),p) = temp
temp = config_phase(p)%getFloats('thermal_expansion22', &
temp = config%getFloats('thermal_expansion22', &
defaultVal=[(0.0_pReal, i=1,size(temp))],requiredSize=size(temp))
!lattice_thermalExpansion33(2,2,1:size(temp),p) = temp
temp = config_phase(p)%getFloats('thermal_expansion33', &
temp = config%getFloats('thermal_expansion33', &
defaultVal=[(0.0_pReal, i=1,size(temp))],requiredSize=size(temp))
end associate
enddo
end subroutine kinematics_thermal_expansion_init
@ -81,7 +80,6 @@ pure function kinematics_thermal_expansion_initialStrain(homog,phase,offset)
real(pReal), dimension(3,3) :: &
kinematics_thermal_expansion_initialStrain !< initial thermal strain (should be small strain, though)
kinematics_thermal_expansion_initialStrain = &
(temperature(homog)%p(offset) - lattice_referenceTemperature(phase))**1 / 1. * &
lattice_thermalExpansion33(1:3,1:3,1,phase) + & ! constant coefficient

15
src/source_damage_anisoBrittle.f90
View File

@ -95,21 +95,13 @@ subroutine source_damage_anisoBrittle_init
enddo
if (all(phase_source(:,p) /= SOURCE_DAMAGE_ANISOBRITTLE_ID)) cycle
associate(prm => param(source_damage_anisoBrittle_instance(p)), &
config => config_phase(p))
prm%aTol = config%getFloat('anisobrittle_atol',defaultVal = 1.0e-3_pReal)
prm%N = config%getFloat('anisobrittle_ratesensitivity')
prm%sdot_0 = config%getFloat('anisobrittle_sdot0')
! sanity checks
if (prm%aTol < 0.0_pReal) extmsg = trim(extmsg)//' anisobrittle_atol'
if (prm%N <= 0.0_pReal) extmsg = trim(extmsg)//' anisobrittle_ratesensitivity'
if (prm%sdot_0 <= 0.0_pReal) extmsg = trim(extmsg)//' anisobrittle_sdot0'
prm%Ncleavage = config%getInts('ncleavage',defaultVal=emptyIntArray)
prm%critDisp = config%getFloats('anisobrittle_criticaldisplacement',requiredSize=size(prm%Ncleavage))
prm%critLoad = config%getFloats('anisobrittle_criticalload', requiredSize=size(prm%Ncleavage))
@ -120,6 +112,9 @@ subroutine source_damage_anisoBrittle_init
prm%critDisp = math_expand(prm%critDisp, prm%Ncleavage)
prm%critLoad = math_expand(prm%critLoad, prm%Ncleavage)
if (prm%aTol < 0.0_pReal) extmsg = trim(extmsg)//' anisobrittle_atol'
if (prm%N <= 0.0_pReal) extmsg = trim(extmsg)//' anisobrittle_ratesensitivity'
if (prm%sdot_0 <= 0.0_pReal) extmsg = trim(extmsg)//' anisobrittle_sdot0'
if (any(prm%critLoad < 0.0_pReal)) extmsg = trim(extmsg)//' anisobrittle_criticalload'
if (any(prm%critDisp < 0.0_pReal)) extmsg = trim(extmsg)//' anisobrittle_criticaldisplacement'
@ -143,8 +138,6 @@ subroutine source_damage_anisoBrittle_init
enddo
end associate
NofMyPhase=count(material_phaseAt==p) * discretization_nIP
instance = source_damage_anisoBrittle_instance(p)
sourceOffset = source_damage_anisoBrittle_offset(p)
@ -153,6 +146,8 @@ subroutine source_damage_anisoBrittle_init
sourceState(p)%p(sourceOffset)%aTolState=param(instance)%aTol
source_damage_anisoBrittle_Ncleavage(1:size(param(instance)%Ncleavage),instance) = param(instance)%Ncleavage
end associate
enddo
end subroutine source_damage_anisoBrittle_init

20
src/source_damage_anisoDuctile.f90
View File

@ -85,25 +85,20 @@ subroutine source_damage_anisoDuctile_init
enddo
if (all(phase_source(:,p) /= SOURCE_DAMAGE_ANISODUCTILE_ID)) cycle
associate(prm => param(source_damage_anisoDuctile_instance(p)), &
config => config_phase(p))
prm%aTol = config%getFloat('anisoductile_atol',defaultVal = 1.0e-3_pReal)
prm%N = config%getFloat('anisoductile_ratesensitivity')
prm%totalNslip = sum(prm%Nslip)
prm%Nslip = config%getInts('nslip',defaultVal=emptyIntArray)
prm%critPlasticStrain = config%getFloats('anisoductile_criticalplasticstrain',requiredSize=size(prm%Nslip))
prm%totalNslip = sum(prm%Nslip) ! expand: family => system
prm%critPlasticStrain = math_expand(prm%critPlasticStrain, prm%Nslip)
! sanity checks
if (prm%aTol < 0.0_pReal) extmsg = trim(extmsg)//' anisoductile_atol'
if (prm%N <= 0.0_pReal) extmsg = trim(extmsg)//' anisoductile_ratesensitivity'
prm%Nslip = config%getInts('nslip',defaultVal=emptyIntArray)
prm%critPlasticStrain = config%getFloats('anisoductile_criticalplasticstrain',requiredSize=size(prm%Nslip))
! expand: family => system
prm%critPlasticStrain = math_expand(prm%critPlasticStrain, prm%Nslip)
if (any(prm%critPlasticStrain < 0.0_pReal)) extmsg = trim(extmsg)//' anisoductile_criticalplasticstrain'
!--------------------------------------------------------------------------------------------------
@ -125,8 +120,6 @@ subroutine source_damage_anisoDuctile_init
enddo
end associate
NofMyPhase=count(material_phaseAt==p) * discretization_nIP
instance = source_damage_anisoDuctile_instance(p)
sourceOffset = source_damage_anisoDuctile_offset(p)
@ -134,6 +127,7 @@ subroutine source_damage_anisoDuctile_init
call material_allocateSourceState(p,sourceOffset,NofMyPhase,1,1,0)
sourceState(p)%p(sourceOffset)%aTolState=param(instance)%aTol
end associate
enddo
end subroutine source_damage_anisoDuctile_init

6
src/source_damage_isoBrittle.f90
View File

@ -80,18 +80,15 @@ subroutine source_damage_isoBrittle_init
enddo
if (all(phase_source(:,p) /= SOURCE_DAMAGE_ISOBRITTLE_ID)) cycle
associate(prm => param(source_damage_isoBrittle_instance(p)), &
config => config_phase(p))
prm%aTol = config%getFloat('isobrittle_atol',defaultVal = 1.0e-3_pReal)
prm%N = config%getFloat('isobrittle_n')
prm%critStrainEnergy = config%getFloat('isobrittle_criticalstrainenergy')
! sanity checks
if (prm%aTol < 0.0_pReal) extmsg = trim(extmsg)//' isobrittle_atol'
if (prm%N <= 0.0_pReal) extmsg = trim(extmsg)//' isobrittle_n'
if (prm%critStrainEnergy <= 0.0_pReal) extmsg = trim(extmsg)//' isobrittle_criticalstrainenergy'
@ -115,8 +112,6 @@ subroutine source_damage_isoBrittle_init
enddo
end associate
NofMyPhase = count(material_phaseAt==p) * discretization_nIP
instance = source_damage_isoBrittle_instance(p)
sourceOffset = source_damage_isoBrittle_offset(p)
@ -124,6 +119,7 @@ subroutine source_damage_isoBrittle_init
call material_allocateSourceState(p,sourceOffset,NofMyPhase,1,1,1)
sourceState(p)%p(sourceOffset)%aTolState=param(instance)%aTol
end associate
enddo
end subroutine source_damage_isoBrittle_init

7
src/source_damage_isoDuctile.f90
View File

@ -77,18 +77,15 @@ subroutine source_damage_isoDuctile_init
enddo
if (all(phase_source(:,p) /= SOURCE_DAMAGE_ISODUCTILE_ID)) cycle
associate(prm => param(source_damage_isoDuctile_instance(p)), &
config => config_phase(p))
prm%aTol = config%getFloat('isoductile_atol',defaultVal = 1.0e-3_pReal)
prm%N = config%getFloat('isoductile_ratesensitivity')
prm%critPlasticStrain = config%getFloat('isoductile_criticalplasticstrain')
! sanity checks
if (prm%aTol < 0.0_pReal) extmsg = trim(extmsg)//' isoductile_atol'
if (prm%N <= 0.0_pReal) extmsg = trim(extmsg)//' isoductile_ratesensitivity'
if (prm%critPlasticStrain <= 0.0_pReal) extmsg = trim(extmsg)//' isoductile_criticalplasticstrain'
@ -112,9 +109,6 @@ subroutine source_damage_isoDuctile_init
enddo
end associate
NofMyPhase=count(material_phaseAt==p) * discretization_nIP
instance = source_damage_isoDuctile_instance(p)
sourceOffset = source_damage_isoDuctile_offset(p)
@ -122,6 +116,7 @@ subroutine source_damage_isoDuctile_init
call material_allocateSourceState(p,sourceOffset,NofMyPhase,1,1,0)
sourceState(p)%p(sourceOffset)%aTolState=param(instance)%aTol
end associate
enddo
end subroutine source_damage_isoDuctile_init

10
src/source_thermal_dissipation.f90
View File

@ -59,14 +59,16 @@ subroutine source_thermal_dissipation_init
enddo
if (all(phase_source(:,p) /= SOURCE_THERMAL_DISSIPATION_ID)) cycle
associate(config => config_phase(p))
instance = source_thermal_dissipation_instance(p)
param(instance)%kappa = config_phase(p)%getFloat('dissipation_coldworkcoeff')
param(instance)%kappa = config%getFloat('dissipation_coldworkcoeff')
NofMyPhase = count(material_phaseAt==p) * discretization_nIP
sourceOffset = source_thermal_dissipation_offset(p)
call material_allocateSourceState(p,sourceOffset,NofMyPhase,0,0,0)
end associate
enddo
end subroutine source_thermal_dissipation_init
@ -75,7 +77,7 @@ end subroutine source_thermal_dissipation_init
!--------------------------------------------------------------------------------------------------
!> @brief Ninstances dissipation rate
!--------------------------------------------------------------------------------------------------
subroutine source_thermal_dissipation_getRateAndItsTangent(TDot, dTDOT_dT, Tstar, Lp, phase)
subroutine source_thermal_dissipation_getRateAndItsTangent(TDot, dTDot_dT, Tstar, Lp, phase)
integer, intent(in) :: &
phase
@ -85,14 +87,14 @@ subroutine source_thermal_dissipation_getRateAndItsTangent(TDot, dTDOT_dT, Tstar
Lp
real(pReal), intent(out) :: &
TDot, &
dTDOT_dT
dTDot_dT
integer :: &
instance
instance = source_thermal_dissipation_instance(phase)
TDot = param(instance)%kappa*sum(abs(Tstar*Lp))
dTDOT_dT = 0.0_pReal
dTDot_dT = 0.0_pReal
end subroutine source_thermal_dissipation_getRateAndItsTangent

6
src/source_thermal_externalheat.f90
View File

@ -63,19 +63,21 @@ subroutine source_thermal_externalheat_init
enddo
if (all(phase_source(:,p) /= SOURCE_thermal_externalheat_ID)) cycle
associate(config => config_phase(p))
instance = source_thermal_externalheat_instance(p)
sourceOffset = source_thermal_externalheat_offset(p)
NofMyPhase = count(material_phaseAt==p) * discretization_nIP
param(instance)%time = config_phase(p)%getFloats('externalheat_time')
param(instance)%time = config%getFloats('externalheat_time')
param(instance)%nIntervals = size(param(instance)%time) - 1
param(instance)%heat_rate = config_phase(p)%getFloats('externalheat_rate',requiredSize = size(param(instance)%time))
param(instance)%heat_rate = config%getFloats('externalheat_rate',requiredSize = size(param(instance)%time))
call material_allocateSourceState(p,sourceOffset,NofMyPhase,1,1,0)
end associate
enddo
end subroutine source_thermal_externalheat_init