simplified

This commit is contained in:
Martin Diehl 2021-01-25 13:49:29 +01:00
parent 2f5c988a89
commit 5592f5aa34
1 changed files with 27 additions and 48 deletions

View File

@ -242,7 +242,7 @@ subroutine materialpoint_stressAndItsTangent(dt,FEsolving_execIP,FEsolving_execE
doneAndHappy doneAndHappy
!$OMP PARALLEL !$OMP PARALLEL
!$OMP DO PRIVATE(ce,me,ho,myNgrains,NiterationMPstate,subFrac,converged,subStep,doneAndHappy) !$OMP DO PRIVATE(ce,me,ho,myNgrains,NiterationMPstate,converged,doneAndHappy)
do el = FEsolving_execElem(1),FEsolving_execElem(2) do el = FEsolving_execElem(1),FEsolving_execElem(2)
ho = material_homogenizationAt(el) ho = material_homogenizationAt(el)
myNgrains = homogenization_Nconstituents(ho) myNgrains = homogenization_Nconstituents(ho)
@ -252,42 +252,18 @@ subroutine materialpoint_stressAndItsTangent(dt,FEsolving_execIP,FEsolving_execE
call constitutive_initializeRestorationPoints(ip,el) call constitutive_initializeRestorationPoints(ip,el)
subFrac = 0.0_pReal call constitutive_restore(ce,.false.) ! wrong name (is more a forward function)
converged = .false. ! pretend failed step ...
subStep = 1.0_pReal/num%subStepSizeHomog ! ... larger then the requested calculation
if (homogState(ho)%sizeState > 0) homogState(ho)%subState0(:,me) = homogState(ho)%State0(:,me) if(homogState(ho)%sizeState > 0) homogState(ho)%State(:,me) = homogState(ho)%State0(:,me)
if (damageState_h(ho)%sizeState > 0) damageState_h(ho)%subState0(:,me) = damageState_h(ho)%State0(:,me) if(damageState_h(ho)%sizeState > 0) damageState_h(ho)%State(:,me) = damageState_h(ho)%State0(:,me)
cutBackLooping: do while (.not. terminallyIll .and. subStep > num%subStepMinHomog) doneAndHappy = [.false.,.true.]
if (converged) then
subFrac = subFrac + subStep
subStep = min(1.0_pReal-subFrac,num%stepIncreaseHomog*subStep) ! introduce flexibility for step increase/acceleration
elseif (subStep <= 1.0 ) then
! cutback makes no sense
if (.not. terminallyIll) & ! so first signals terminally ill...
print*, ' Integration point ', ip,' at element ', el, ' terminally ill'
terminallyIll = .true. ! ...and kills all others
else ! cutback makes sense
subStep = num%subStepSizeHomog * subStep ! crystallite had severe trouble, so do a significant cutback
call constitutive_restore(ce,subStep < 1.0_pReal)
if(homogState(ho)%sizeState > 0) homogState(ho)%State(:,me) = homogState(ho)%subState0(:,me)
if(damageState_h(ho)%sizeState > 0) damageState_h(ho)%State(:,me) = damageState_h(ho)%subState0(:,me)
endif
if (subStep > num%subStepMinHomog) doneAndHappy = [.false.,.true.]
NiterationMPstate = 0 NiterationMPstate = 0
convergenceLooping: do while (.not. (terminallyIll .or. doneAndHappy(1)) & convergenceLooping: do while (.not. (terminallyIll .or. doneAndHappy(1)) &
.and. NiterationMPstate < num%nMPstate) .and. NiterationMPstate < num%nMPstate)
NiterationMPstate = NiterationMPstate + 1 NiterationMPstate = NiterationMPstate + 1
!--------------------------------------------------------------------------------------------------
! deformation partitioning
if (.not. doneAndHappy(1)) then if (.not. doneAndHappy(1)) then
call mech_partition(homogenization_F(1:3,1:3,ce),ip,el) call mech_partition(homogenization_F(1:3,1:3,ce),ip,el)
@ -305,7 +281,10 @@ subroutine materialpoint_stressAndItsTangent(dt,FEsolving_execIP,FEsolving_execE
endif endif
enddo convergenceLooping enddo convergenceLooping
enddo cutBackLooping if (.not. converged) then
if (.not. terminallyIll) print*, ' Integration point ', ip,' at element ', el, ' terminally ill'
terminallyIll = .true.
endif
enddo enddo
enddo enddo
!$OMP END DO !$OMP END DO
@ -337,7 +316,7 @@ subroutine materialpoint_stressAndItsTangent(dt,FEsolving_execIP,FEsolving_execE
ho = material_homogenizationAt(el) ho = material_homogenizationAt(el)
do ip = FEsolving_execIP(1),FEsolving_execIP(2) do ip = FEsolving_execIP(1),FEsolving_execIP(2)
ce = (el-1)*discretization_nIPs + ip ce = (el-1)*discretization_nIPs + ip
call damage_partition(ce) ! call damage_partition(ce)
! do co = 1, homogenization_Nconstituents(ho) ! do co = 1, homogenization_Nconstituents(ho)
! ph = material_phaseAt(co,el) ! ph = material_phaseAt(co,el)
! if (.not. thermal_stress(dt,ph,material_phaseMemberAt(co,ip,el))) then ! if (.not. thermal_stress(dt,ph,material_phaseMemberAt(co,ip,el))) then