diff --git a/src/constitutive.f90 b/src/constitutive.f90 index c85ae0553..1a03f3c50 100644 --- a/src/constitutive.f90 +++ b/src/constitutive.f90 @@ -49,8 +49,7 @@ module constitutive real(pReal), dimension(:,:,:), allocatable, public :: & crystallite_dt !< requested time increment of each grain real(pReal), dimension(:,:,:), allocatable :: & - crystallite_subdt, & !< substepped time increment of each grain - crystallite_subStep !< size of next integration step + crystallite_subdt !< substepped time increment of each grain type(rotation), dimension(:,:,:), allocatable :: & crystallite_orientation !< current orientation real(pReal), dimension(:,:,:,:,:), allocatable :: & @@ -882,7 +881,7 @@ subroutine crystallite_init source = crystallite_partitionedF) allocate(crystallite_dt(cMax,iMax,eMax),source=0.0_pReal) - allocate(crystallite_subdt,crystallite_subStep, & + allocate(crystallite_subdt, & source = crystallite_dt) allocate(crystallite_orientation(cMax,iMax,eMax)) @@ -1015,16 +1014,12 @@ function crystallite_stress(co,ip,el) NiterationCrystallite, & ! number of iterations in crystallite loop s, ph, me logical :: todo - real(pReal) :: subFrac !ToDo: need to set some values to false for different Ngrains + real(pReal) :: subFrac,subStep real(pReal), dimension(3,3) :: & subLp0, & !< plastic velocity grad at start of crystallite inc subLi0 !< intermediate velocity grad at start of crystallite inc -!-------------------------------------------------------------------------------------------------- -! initialize to starting condition - crystallite_subStep(co,ip,el) = 0.0_pReal - ph = material_phaseAt(co,el) me = material_phaseMemberAt(co,ip,el) subLi0 = constitutive_mech_partionedLi0(ph)%data(1:3,1:3,me) @@ -1040,7 +1035,7 @@ function crystallite_stress(co,ip,el) crystallite_subFi0(1:3,1:3,co,ip,el) = constitutive_mech_partionedFi0(ph)%data(1:3,1:3,me) crystallite_subF0(1:3,1:3,co,ip,el) = crystallite_partitionedF0(1:3,1:3,co,ip,el) subFrac = 0.0_pReal - crystallite_subStep(co,ip,el) = 1.0_pReal/num%subStepSizeCryst + subStep = 1.0_pReal/num%subStepSizeCryst todo = .true. crystallite_converged(co,ip,el) = .false. ! pretend failed step of 1/subStepSizeCryst @@ -1052,12 +1047,11 @@ function crystallite_stress(co,ip,el) !-------------------------------------------------------------------------------------------------- ! wind forward if (crystallite_converged(co,ip,el)) then - formerSubStep = crystallite_subStep(co,ip,el) - subFrac = subFrac + crystallite_subStep(co,ip,el) - crystallite_subStep(co,ip,el) = min(1.0_pReal - subFrac, & - num%stepIncreaseCryst * crystallite_subStep(co,ip,el)) + formerSubStep = subStep + subFrac = subFrac + subStep + subStep = min(1.0_pReal - subFrac, num%stepIncreaseCryst * subStep) - todo = crystallite_subStep(co,ip,el) > 0.0_pReal ! still time left to integrate on? + todo = subStep > 0.0_pReal ! still time left to integrate on? if (todo) then crystallite_subF0 (1:3,1:3,co,ip,el) = crystallite_subF(1:3,1:3,co,ip,el) @@ -1076,11 +1070,11 @@ function crystallite_stress(co,ip,el) !-------------------------------------------------------------------------------------------------- ! cut back (reduced time and restore) else - crystallite_subStep(co,ip,el) = num%subStepSizeCryst * crystallite_subStep(co,ip,el) + subStep = num%subStepSizeCryst * subStep constitutive_mech_Fp(ph)%data(1:3,1:3,me) = crystallite_subFp0(1:3,1:3,co,ip,el) constitutive_mech_Fi(ph)%data(1:3,1:3,me) = crystallite_subFi0(1:3,1:3,co,ip,el) crystallite_S (1:3,1:3,co,ip,el) = crystallite_S0 (1:3,1:3,co,ip,el) - if (crystallite_subStep(co,ip,el) < 1.0_pReal) then ! actual (not initial) cutback + if (subStep < 1.0_pReal) then ! actual (not initial) cutback crystallite_Lp (1:3,1:3,co,ip,el) = subLp0 constitutive_mech_Li(ph)%data(1:3,1:3,me) = subLi0 endif @@ -1092,25 +1086,25 @@ function crystallite_stress(co,ip,el) enddo ! cant restore dotState here, since not yet calculated in first cutback after initialization - todo = crystallite_subStep(co,ip,el) > num%subStepMinCryst ! still on track or already done (beyond repair) + todo = subStep > num%subStepMinCryst ! still on track or already done (beyond repair) endif !-------------------------------------------------------------------------------------------------- ! prepare for integration if (todo) then crystallite_subF(1:3,1:3,co,ip,el) = crystallite_subF0(1:3,1:3,co,ip,el) & - + crystallite_subStep(co,ip,el) *( crystallite_partitionedF (1:3,1:3,co,ip,el) & - -crystallite_partitionedF0(1:3,1:3,co,ip,el)) + + subStep *( crystallite_partitionedF (1:3,1:3,co,ip,el) & + -crystallite_partitionedF0(1:3,1:3,co,ip,el)) crystallite_Fe(1:3,1:3,co,ip,el) = matmul(crystallite_subF(1:3,1:3,co,ip,el), & math_inv33(matmul(constitutive_mech_Fi(ph)%data(1:3,1:3,me), & constitutive_mech_Fp(ph)%data(1:3,1:3,me)))) - crystallite_subdt(co,ip,el) = crystallite_subStep(co,ip,el) * crystallite_dt(co,ip,el) + crystallite_subdt(co,ip,el) = subStep * crystallite_dt(co,ip,el) crystallite_converged(co,ip,el) = .false. call integrateState(co,ip,el) call integrateSourceState(co,ip,el) endif - if (.not. crystallite_converged(co,ip,el) .and. crystallite_subStep(co,ip,el) > num%subStepMinCryst) & ! do not try non-converged but fully cutbacked any further + if (.not. crystallite_converged(co,ip,el) .and. subStep > num%subStepMinCryst) & ! do not try non-converged but fully cutbacked any further todo = .true. enddo cutbackLooping