This commit is contained in:
Martin Diehl 2019-04-05 20:45:56 +02:00
parent adebbcf5df
commit 4aa52fa83f
3 changed files with 297 additions and 315 deletions

View File

@ -6,7 +6,6 @@
!--------------------------------------------------------------------------------------------------
module homogenization
use prec, only: &
pInt, &
pReal
!--------------------------------------------------------------------------------------------------
@ -21,7 +20,7 @@ module homogenization
materialpoint_dPdF !< tangent of first P--K stress at IP
real(pReal), dimension(:,:,:), allocatable, public :: &
materialpoint_results !< results array of material point
integer(pInt), public, protected :: &
integer, public, protected :: &
materialpoint_sizeResults, &
homogenization_maxSizePostResults, &
thermal_maxSizePostResults, &
@ -92,10 +91,10 @@ subroutine homogenization_init
worldrank
implicit none
integer(pInt), parameter :: FILEUNIT = 200_pInt
integer(pInt) :: e,i,p
integer(pInt), dimension(:,:), pointer :: thisSize
integer(pInt), dimension(:) , pointer :: thisNoutput
integer, parameter :: FILEUNIT = 200
integer :: e,i,p
integer, dimension(:,:), pointer :: thisSize
integer, dimension(:) , pointer :: thisNoutput
character(len=64), dimension(:,:), pointer :: thisOutput
character(len=32) :: outputName !< name of output, intermediate fix until HDF5 output is ready
logical :: valid
@ -232,9 +231,9 @@ subroutine homogenization_init
!--------------------------------------------------------------------------------------------------
! allocate and initialize global state and postresutls variables
homogenization_maxSizePostResults = 0_pInt
thermal_maxSizePostResults = 0_pInt
damage_maxSizePostResults = 0_pInt
homogenization_maxSizePostResults = 0
thermal_maxSizePostResults = 0
damage_maxSizePostResults = 0
do p = 1,size(config_homogenization)
homogenization_maxSizePostResults = max(homogenization_maxSizePostResults,homogState (p)%sizePostResults)
thermal_maxSizePostResults = max(thermal_maxSizePostResults, thermalState (p)%sizePostResults)
@ -252,7 +251,7 @@ subroutine homogenization_init
write(6,'(/,a)') ' <<<+- homogenization init -+>>>'
if (iand(debug_level(debug_homogenization), debug_levelBasic) /= 0_pInt) then
if (iand(debug_level(debug_homogenization), debug_levelBasic) /= 0) then
#ifdef TODO
write(6,'(a32,1x,7(i8,1x))') 'homogenization_state0: ', shape(homogenization_state0)
write(6,'(a32,1x,7(i8,1x))') 'homogenization_subState0: ', shape(homogenization_subState0)
@ -275,7 +274,7 @@ subroutine homogenization_init
flush(6)
if (debug_g < 1 .or. debug_g > homogenization_Ngrains(mesh_element(3,debug_e))) &
call IO_error(602_pInt,ext_msg='constituent', el=debug_e, g=debug_g)
call IO_error(602,ext_msg='constituent', el=debug_e, g=debug_g)
end subroutine homogenization_init
@ -344,7 +343,7 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt)
implicit none
real(pReal), intent(in) :: dt !< time increment
logical, intent(in) :: updateJaco !< initiating Jacobian update
integer(pInt) :: &
integer :: &
NiterationHomog, &
NiterationMPstate, &
g, & !< grain number
@ -354,7 +353,7 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt)
myNgrains
#ifdef DEBUG
if (iand(debug_level(debug_homogenization), debug_levelBasic) /= 0_pInt) then
if (iand(debug_level(debug_homogenization), debug_levelBasic) /= 0) then
write(6,'(/a,i5,1x,i2)') '<< HOMOG >> Material Point start at el ip ', debug_e, debug_i
write(6,'(a,/,3(12x,3(f14.9,1x)/))') '<< HOMOG >> F0', &
@ -372,7 +371,7 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt)
plasticState (phaseAt(g,i,e))%partionedState0(:,phasememberAt(g,i,e)) = &
plasticState (phaseAt(g,i,e))%state0( :,phasememberAt(g,i,e))
do mySource = 1_pInt, phase_Nsources(phaseAt(g,i,e))
do mySource = 1, phase_Nsources(phaseAt(g,i,e))
sourceState(phaseAt(g,i,e))%p(mySource)%partionedState0(:,phasememberAt(g,i,e)) = &
sourceState(phaseAt(g,i,e))%p(mySource)%state0( :,phasememberAt(g,i,e))
enddo
@ -393,19 +392,19 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt)
materialpoint_requested(i,e) = .true. ! everybody requires calculation
endforall
forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), &
homogState(material_homogenizationAt(e))%sizeState > 0_pInt) &
homogState(material_homogenizationAt(e))%sizeState > 0) &
homogState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = &
homogState(material_homogenizationAt(e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal homogenization state
forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), &
thermalState(material_homogenizationAt(e))%sizeState > 0_pInt) &
thermalState(material_homogenizationAt(e))%sizeState > 0) &
thermalState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = &
thermalState(material_homogenizationAt(e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal thermal state
forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), &
damageState(material_homogenizationAt(e))%sizeState > 0_pInt) &
damageState(material_homogenizationAt(e))%sizeState > 0) &
damageState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = &
damageState(material_homogenizationAt(e))%State0( :,mappingHomogenization(1,i,e)) ! ...internal damage state
enddo
NiterationHomog = 0_pInt
NiterationHomog = 0
cutBackLooping: do while (.not. terminallyIll .and. &
any(materialpoint_subStep(:,FEsolving_execELem(1):FEsolving_execElem(2)) > subStepMinHomog))
@ -417,9 +416,9 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt)
converged: if ( materialpoint_converged(i,e) ) then
#ifdef DEBUG
if (iand(debug_level(debug_homogenization), debug_levelExtensive) /= 0_pInt &
if (iand(debug_level(debug_homogenization), debug_levelExtensive) /= 0 &
.and. ((e == debug_e .and. i == debug_i) &
.or. .not. iand(debug_level(debug_homogenization),debug_levelSelective) /= 0_pInt)) then
.or. .not. iand(debug_level(debug_homogenization),debug_levelSelective) /= 0)) then
write(6,'(a,1x,f12.8,1x,a,1x,f12.8,1x,a,i8,1x,i2/)') '<< HOMOG >> winding forward from', &
materialpoint_subFrac(i,e), 'to current materialpoint_subFrac', &
materialpoint_subFrac(i,e)+materialpoint_subStep(i,e),'in materialpoint_stressAndItsTangent at el ip',e,i
@ -456,29 +455,29 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt)
do g = 1,myNgrains
plasticState (phaseAt(g,i,e))%partionedState0(:,phasememberAt(g,i,e)) = &
plasticState (phaseAt(g,i,e))%state( :,phasememberAt(g,i,e))
do mySource = 1_pInt, phase_Nsources(phaseAt(g,i,e))
do mySource = 1, phase_Nsources(phaseAt(g,i,e))
sourceState(phaseAt(g,i,e))%p(mySource)%partionedState0(:,phasememberAt(g,i,e)) = &
sourceState(phaseAt(g,i,e))%p(mySource)%state( :,phasememberAt(g,i,e))
enddo
enddo
forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), &
homogState(material_homogenizationAt(e))%sizeState > 0_pInt) &
homogState(material_homogenizationAt(e))%sizeState > 0) &
homogState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = &
homogState(material_homogenizationAt(e))%State( :,mappingHomogenization(1,i,e)) ! ...internal homogenization state
forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), &
thermalState(material_homogenizationAt(e))%sizeState > 0_pInt) &
thermalState(material_homogenizationAt(e))%sizeState > 0) &
thermalState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = &
thermalState(material_homogenizationAt(e))%State( :,mappingHomogenization(1,i,e)) ! ...internal thermal state
forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), &
damageState(material_homogenizationAt(e))%sizeState > 0_pInt) &
damageState(material_homogenizationAt(e))%sizeState > 0) &
damageState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) = &
damageState(material_homogenizationAt(e))%State( :,mappingHomogenization(1,i,e)) ! ...internal damage state
materialpoint_subF0(1:3,1:3,i,e) = materialpoint_subF(1:3,1:3,i,e) ! ...def grad
endif steppingNeeded
else converged
if ( (myNgrains == 1_pInt .and. materialpoint_subStep(i,e) <= 1.0 ) .or. & ! single grain already tried internal subStepping in crystallite
if ( (myNgrains == 1 .and. materialpoint_subStep(i,e) <= 1.0 ) .or. & ! single grain already tried internal subStepping in crystallite
subStepSizeHomog * materialpoint_subStep(i,e) <= subStepMinHomog ) then ! would require too small subStep
! cutback makes no sense
!$OMP FLUSH(terminallyIll)
@ -494,9 +493,9 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt)
materialpoint_subStep(i,e) = subStepSizeHomog * materialpoint_subStep(i,e) ! crystallite had severe trouble, so do a significant cutback
#ifdef DEBUG
if (iand(debug_level(debug_homogenization), debug_levelExtensive) /= 0_pInt &
if (iand(debug_level(debug_homogenization), debug_levelExtensive) /= 0 &
.and. ((e == debug_e .and. i == debug_i) &
.or. .not. iand(debug_level(debug_homogenization), debug_levelSelective) /= 0_pInt)) then
.or. .not. iand(debug_level(debug_homogenization), debug_levelSelective) /= 0)) then
write(6,'(a,1x,f12.8,a,i8,1x,i2/)') &
'<< HOMOG >> cutback step in materialpoint_stressAndItsTangent with new materialpoint_subStep:',&
materialpoint_subStep(i,e),' at el ip',e,i
@ -518,21 +517,21 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt)
do g = 1, myNgrains
plasticState (phaseAt(g,i,e))%state( :,phasememberAt(g,i,e)) = &
plasticState (phaseAt(g,i,e))%partionedState0(:,phasememberAt(g,i,e))
do mySource = 1_pInt, phase_Nsources(phaseAt(g,i,e))
do mySource = 1, phase_Nsources(phaseAt(g,i,e))
sourceState(phaseAt(g,i,e))%p(mySource)%state( :,phasememberAt(g,i,e)) = &
sourceState(phaseAt(g,i,e))%p(mySource)%partionedState0(:,phasememberAt(g,i,e))
enddo
enddo
forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), &
homogState(material_homogenizationAt(e))%sizeState > 0_pInt) &
homogState(material_homogenizationAt(e))%sizeState > 0) &
homogState(material_homogenizationAt(e))%State( :,mappingHomogenization(1,i,e)) = &
homogState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) ! ...internal homogenization state
forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), &
thermalState(material_homogenizationAt(e))%sizeState > 0_pInt) &
thermalState(material_homogenizationAt(e))%sizeState > 0) &
thermalState(material_homogenizationAt(e))%State( :,mappingHomogenization(1,i,e)) = &
thermalState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) ! ...internal thermal state
forall(i = FEsolving_execIP(1,e):FEsolving_execIP(2,e), &
damageState(material_homogenizationAt(e))%sizeState > 0_pInt) &
damageState(material_homogenizationAt(e))%sizeState > 0) &
damageState(material_homogenizationAt(e))%State( :,mappingHomogenization(1,i,e)) = &
damageState(material_homogenizationAt(e))%subState0(:,mappingHomogenization(1,i,e)) ! ...internal damage state
endif
@ -550,7 +549,7 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt)
enddo elementLooping1
!$OMP END PARALLEL DO
NiterationMPstate = 0_pInt
NiterationMPstate = 0
convergenceLooping: do while (.not. terminallyIll .and. &
any( materialpoint_requested(:,FEsolving_execELem(1):FEsolving_execElem(2)) &
@ -606,7 +605,7 @@ subroutine materialpoint_stressAndItsTangent(updateJaco,dt)
enddo convergenceLooping
NiterationHomog = NiterationHomog + 1_pInt
NiterationHomog = NiterationHomog + 1
enddo cutBackLooping
@ -652,7 +651,7 @@ subroutine materialpoint_postResults
crystallite_postResults
implicit none
integer(pInt) :: &
integer :: &
thePos, &
theSize, &
myNgrains, &
@ -666,21 +665,21 @@ subroutine materialpoint_postResults
myNgrains = homogenization_Ngrains(mesh_element(3,e))
myCrystallite = microstructure_crystallite(mesh_element(4,e))
IpLooping: do i = FEsolving_execIP(1,e),FEsolving_execIP(2,e)
thePos = 0_pInt
thePos = 0
theSize = homogState (material_homogenizationAt(e))%sizePostResults &
+ thermalState (material_homogenizationAt(e))%sizePostResults &
+ damageState (material_homogenizationAt(e))%sizePostResults
materialpoint_results(thePos+1,i,e) = real(theSize,pReal) ! tell size of homogenization results
thePos = thePos + 1_pInt
thePos = thePos + 1
if (theSize > 0_pInt) then ! any homogenization results to mention?
if (theSize > 0) then ! any homogenization results to mention?
materialpoint_results(thePos+1:thePos+theSize,i,e) = postResults(i,e) ! tell homogenization results
thePos = thePos + theSize
endif
materialpoint_results(thePos+1,i,e) = real(myNgrains,pReal) ! tell number of grains at materialpoint
thePos = thePos + 1_pInt
thePos = thePos + 1
grainLooping :do g = 1,myNgrains
theSize = 1 + crystallite_sizePostResults(myCrystallite) + &
@ -716,7 +715,7 @@ subroutine partitionDeformation(ip,el)
homogenization_RGC_partitionDeformation
implicit none
integer(pInt), intent(in) :: &
integer, intent(in) :: &
ip, & !< integration point
el !< element number
@ -769,7 +768,7 @@ function updateState(ip,el)
damage_local_updateState
implicit none
integer(pInt), intent(in) :: &
integer, intent(in) :: &
ip, & !< integration point
el !< element number
logical, dimension(2) :: updateState
@ -831,7 +830,7 @@ subroutine averageStressAndItsTangent(ip,el)
homogenization_RGC_averageStressAndItsTangent
implicit none
integer(pInt), intent(in) :: &
integer, intent(in) :: &
ip, & !< integration point
el !< element number
@ -900,20 +899,20 @@ function postResults(ip,el)
damage_nonlocal_postResults
implicit none
integer(pInt), intent(in) :: &
integer, intent(in) :: &
ip, & !< integration point
el !< element number
real(pReal), dimension( homogState (material_homogenizationAt(el))%sizePostResults &
+ thermalState (material_homogenizationAt(el))%sizePostResults &
+ damageState (material_homogenizationAt(el))%sizePostResults) :: &
postResults
integer(pInt) :: &
integer :: &
startPos, endPos ,&
of, instance, homog
postResults = 0.0_pReal
startPos = 1_pInt
startPos = 1
endPos = homogState(material_homogenizationAt(el))%sizePostResults
chosenHomogenization: select case (homogenization_type(mesh_element(3,el)))
@ -924,7 +923,7 @@ function postResults(ip,el)
end select chosenHomogenization
startPos = endPos + 1_pInt
startPos = endPos + 1
endPos = endPos + thermalState(material_homogenizationAt(el))%sizePostResults
chosenThermal: select case (thermal_type(mesh_element(3,el)))
@ -939,7 +938,7 @@ function postResults(ip,el)
end select chosenThermal
startPos = endPos + 1_pInt
startPos = endPos + 1
endPos = endPos + damageState(material_homogenizationAt(el))%sizePostResults
chosenDamage: select case (damage_type(mesh_element(3,el)))

File diff suppressed because it is too large Load Diff

View File

@ -17,7 +17,7 @@ module homogenization_isostrain
end enum
type, private :: tParameters !< container type for internal constitutive parameters
integer(pInt) :: &
integer :: &
Nconstituents
integer(kind(average_ID)) :: &
mapping
@ -53,7 +53,7 @@ subroutine homogenization_isostrain_init()
config_homogenization
implicit none
integer(pInt) :: &
integer :: &
Ninstance, &
h, &
NofMyHomog
@ -63,12 +63,12 @@ subroutine homogenization_isostrain_init()
write(6,'(/,a)') ' <<<+- homogenization_'//HOMOGENIZATION_ISOSTRAIN_label//' init -+>>>'
Ninstance = int(count(homogenization_type == HOMOGENIZATION_ISOSTRAIN_ID),pInt)
if (iand(debug_level(debug_HOMOGENIZATION),debug_levelBasic) /= 0_pInt) &
if (iand(debug_level(debug_HOMOGENIZATION),debug_levelBasic) /= 0) &
write(6,'(a16,1x,i5,/)') '# instances:',Ninstance
allocate(param(Ninstance)) ! one container of parameters per instance
do h = 1_pInt, size(homogenization_type)
do h = 1, size(homogenization_type)
if (homogenization_type(h) /= HOMOGENIZATION_ISOSTRAIN_ID) cycle
associate(prm => param(homogenization_typeInstance(h)),&
@ -82,15 +82,15 @@ subroutine homogenization_isostrain_init()
case ('avg')
prm%mapping = average_ID
case default
call IO_error(211_pInt,ext_msg=trim(tag)//' ('//HOMOGENIZATION_isostrain_label//')')
call IO_error(211,ext_msg=trim(tag)//' ('//HOMOGENIZATION_isostrain_label//')')
end select
NofMyHomog = count(material_homogenizationAt == h)
homogState(h)%sizeState = 0_pInt
homogState(h)%sizePostResults = 0_pInt
allocate(homogState(h)%state0 (0_pInt,NofMyHomog))
allocate(homogState(h)%subState0(0_pInt,NofMyHomog))
allocate(homogState(h)%state (0_pInt,NofMyHomog))
homogState(h)%sizeState = 0
homogState(h)%sizePostResults = 0
allocate(homogState(h)%state0 (0,NofMyHomog))
allocate(homogState(h)%subState0(0,NofMyHomog))
allocate(homogState(h)%state (0,NofMyHomog))
end associate
@ -129,7 +129,7 @@ subroutine homogenization_isostrain_averageStressAndItsTangent(avgP,dAvgPdAvgF,P
real(pReal), dimension (:,:,:), intent(in) :: P !< partitioned stresses
real(pReal), dimension (:,:,:,:,:), intent(in) :: dPdF !< partitioned stiffnesses
integer(pInt), intent(in) :: instance
integer, intent(in) :: instance
associate(prm => param(instance))