From 38ec462f4a1a73608ccbba21b3a359f792dd0cec Mon Sep 17 00:00:00 2001 From: Christoph Kords Date: Thu, 22 Nov 2012 13:04:19 +0000 Subject: [PATCH] order of integration now same for all state integrators: first calculate the dotState and then the deltaState --- code/crystallite.f90 | 275 ++++++++++++++++++++++++------------------- 1 file changed, 154 insertions(+), 121 deletions(-) diff --git a/code/crystallite.f90 b/code/crystallite.f90 index 7b6f54224..3766b6e2e 100644 --- a/code/crystallite.f90 +++ b/code/crystallite.f90 @@ -814,6 +814,7 @@ if(updateJaco) then crystallite_Temperature = crystallite_subTemperature0 crystallite_Fp = crystallite_subFp0 crystallite_Fe = crystallite_subFe0 + crystallite_Lp = crystallite_subLp0 crystallite_Tstar_v = crystallite_subTstar0_v end select @@ -1140,6 +1141,22 @@ do n = 1_pInt,4_pInt !$OMP ENDDO + ! --- state jump --- + + !$OMP DO + do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains + if (crystallite_todo(g,i,e)) then + crystallite_todo(g,i,e) = crystallite_stateJump(g,i,e) + if (.not. crystallite_todo(g,i,e) .and. .not. crystallite_localPlasticity(g,i,e)) then ! if broken non-local... + !$OMP CRITICAL (checkTodo) + crystallite_todo = crystallite_todo .and. crystallite_localPlasticity ! ...all non-locals skipped + !$OMP END CRITICAL (checkTodo) + endif + endif + enddo; enddo; enddo + !$OMP ENDDO + + ! --- update dependent states --- !$OMP DO @@ -1203,22 +1220,6 @@ do n = 1_pInt,4_pInt enddo -! --- STATE JUMP --- - -!$OMP DO - do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains - if (crystallite_todo(g,i,e)) then - crystallite_todo(g,i,e) = crystallite_stateJump(g,i,e) - if (.not. crystallite_todo(g,i,e) .and. .not. crystallite_localPlasticity(g,i,e)) then ! if broken non-local... - !$OMP CRITICAL (checkTodo) - crystallite_todo = crystallite_todo .and. crystallite_localPlasticity ! ...all non-locals skipped - !$OMP END CRITICAL (checkTodo) - endif - endif - enddo; enddo; enddo -!$OMP ENDDO - - ! --- SET CONVERGENCE FLAG --- !$OMP DO @@ -1482,6 +1483,22 @@ do n = 1_pInt,5_pInt enddo; enddo; enddo !$OMP ENDDO + + ! --- state jump --- + + !$OMP DO + do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains + if (crystallite_todo(g,i,e)) then + crystallite_todo(g,i,e) = crystallite_stateJump(g,i,e) + if (.not. crystallite_todo(g,i,e) .and. .not. crystallite_localPlasticity(g,i,e)) then ! if broken non-local... + !$OMP CRITICAL (checkTodo) + crystallite_todo = crystallite_todo .and. crystallite_localPlasticity ! ...all non-locals skipped + !$OMP END CRITICAL (checkTodo) + endif + endif + enddo; enddo; enddo + !$OMP ENDDO + ! --- update dependent states --- @@ -1660,6 +1677,22 @@ relTemperatureResiduum = 0.0_pReal !$OMP ENDDO +! --- STATE JUMP --- + +!$OMP DO + do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains + if (crystallite_todo(g,i,e)) then + crystallite_todo(g,i,e) = crystallite_stateJump(g,i,e) + if (.not. crystallite_todo(g,i,e) .and. .not. crystallite_localPlasticity(g,i,e)) then ! if broken non-local... + !$OMP CRITICAL (checkTodo) + crystallite_todo = crystallite_todo .and. crystallite_localPlasticity ! ...all non-locals skipped + !$OMP END CRITICAL (checkTodo) + endif + endif + enddo; enddo; enddo +!$OMP ENDDO + + ! --- UPDATE DEPENDENT STATES IF RESIDUUM BELOW TOLERANCE --- !$OMP DO @@ -1688,22 +1721,6 @@ relTemperatureResiduum = 0.0_pReal !$OMP ENDDO -! --- STATE JUMP --- - -!$OMP DO - do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains - if (crystallite_todo(g,i,e)) then - crystallite_todo(g,i,e) = crystallite_stateJump(g,i,e) - if (.not. crystallite_todo(g,i,e) .and. .not. crystallite_localPlasticity(g,i,e)) then ! if broken non-local... - !$OMP CRITICAL (checkTodo) - crystallite_todo = crystallite_todo .and. crystallite_localPlasticity ! ...all non-locals skipped - !$OMP END CRITICAL (checkTodo) - endif - endif - enddo; enddo; enddo -!$OMP ENDDO - - ! --- SET CONVERGENCE FLAG --- !$OMP DO @@ -1822,23 +1839,7 @@ endif !$OMP PARALLEL PRIVATE(mySizeDotState) -if (numerics_integrationMode < 2) then - - ! --- STATE JUMP --- - - !$OMP DO - do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains - if (crystallite_todo(g,i,e)) then - crystallite_todo(g,i,e) = crystallite_stateJump(g,i,e) - if (.not. crystallite_todo(g,i,e) .and. .not. crystallite_localPlasticity(g,i,e)) then ! if broken non-local... - !$OMP CRITICAL (checkTodo) - crystallite_todo = crystallite_todo .and. crystallite_localPlasticity ! ...all non-locals skipped - !$OMP END CRITICAL (checkTodo) - endif - endif - enddo; enddo; enddo - !$OMP ENDDO - +if (numerics_integrationMode == 1_pInt) then ! --- DOT STATE AND TEMPERATURE (EULER INTEGRATION) --- @@ -1888,6 +1889,22 @@ if (numerics_integrationMode < 2) then !$OMP ENDDO + ! --- STATE JUMP --- + + !$OMP DO + do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains + if (crystallite_todo(g,i,e)) then + crystallite_todo(g,i,e) = crystallite_stateJump(g,i,e) + if (.not. crystallite_todo(g,i,e) .and. .not. crystallite_localPlasticity(g,i,e)) then ! if broken non-local... + !$OMP CRITICAL (checkTodo) + crystallite_todo = crystallite_todo .and. crystallite_localPlasticity ! ...all non-locals skipped + !$OMP END CRITICAL (checkTodo) + endif + endif + enddo; enddo; enddo + !$OMP ENDDO + + ! --- UPDATE DEPENDENT STATES (EULER INTEGRATION) --- !$OMP DO @@ -1901,7 +1918,6 @@ if (numerics_integrationMode < 2) then endif - ! --- STRESS INTEGRATION (EULER INTEGRATION) --- !$OMP DO @@ -1918,7 +1934,7 @@ endif !$OMP ENDDO -if (numerics_integrationMode < 2) then +if (numerics_integrationMode == 1_pInt) then ! --- DOT STATE AND TEMPERATURE (HEUN METHOD) --- @@ -1949,61 +1965,77 @@ if (numerics_integrationMode < 2) then enddo; enddo; enddo !$OMP ENDDO -endif + ! --- ERROR ESTIMATE FOR STATE AND TEMPERATURE (HEUN METHOD) --- -! --- ERROR ESTIMATE FOR STATE AND TEMPERATURE (HEUN METHOD) --- + relStateResiduum = 0.0_pReal + relTemperatureResiduum = 0.0_pReal + !$OMP DO + do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains + if (crystallite_todo(g,i,e)) then + mySizeDotState = constitutive_sizeDotState(g,i,e) + + + ! --- contribution of heun step to absolute residui --- + + stateResiduum(1:mySizeDotState,g,i,e) = stateResiduum(1:mySizeDotState,g,i,e) & + + 0.5_pReal * constitutive_dotState(g,i,e)%p * crystallite_subdt(g,i,e) ! contribution to absolute residuum in state and temperature + temperatureResiduum(g,i,e) = temperatureResiduum(g,i,e) & + + 0.5_pReal * crystallite_dotTemperature(g,i,e) * crystallite_subdt(g,i,e) + !$OMP FLUSH(stateResiduum,temperatureResiduum) -relStateResiduum = 0.0_pReal -relTemperatureResiduum = 0.0_pReal -!$OMP DO - do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains - if (crystallite_todo(g,i,e)) then - mySizeDotState = constitutive_sizeDotState(g,i,e) - - - ! --- contribution of heun step to absolute residui --- - - stateResiduum(1:mySizeDotState,g,i,e) = stateResiduum(1:mySizeDotState,g,i,e) & - + 0.5_pReal * constitutive_dotState(g,i,e)%p * crystallite_subdt(g,i,e) ! contribution to absolute residuum in state and temperature - temperatureResiduum(g,i,e) = temperatureResiduum(g,i,e) & - + 0.5_pReal * crystallite_dotTemperature(g,i,e) * crystallite_subdt(g,i,e) - !$OMP FLUSH(stateResiduum,temperatureResiduum) + ! --- relative residui --- - ! --- relative residui --- - - forall (s = 1_pInt:mySizeDotState, abs(constitutive_state(g,i,e)%p(s)) > 0.0_pReal) & - relStateResiduum(s,g,i,e) = stateResiduum(s,g,i,e) / constitutive_state(g,i,e)%p(s) - if (crystallite_Temperature(g,i,e) > 0_pInt) & - relTemperatureResiduum(g,i,e) = temperatureResiduum(g,i,e) / crystallite_Temperature(g,i,e) - !$OMP FLUSH(relStateResiduum,relTemperatureResiduum) + forall (s = 1_pInt:mySizeDotState, abs(constitutive_state(g,i,e)%p(s)) > 0.0_pReal) & + relStateResiduum(s,g,i,e) = stateResiduum(s,g,i,e) / constitutive_state(g,i,e)%p(s) + if (crystallite_Temperature(g,i,e) > 0_pInt) & + relTemperatureResiduum(g,i,e) = temperatureResiduum(g,i,e) / crystallite_Temperature(g,i,e) + !$OMP FLUSH(relStateResiduum,relTemperatureResiduum) #ifndef _OPENMP - if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt & - .and. ((e == debug_e .and. i == debug_i .and. g == debug_g)& - .or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0_pInt)) then - write(6,'(a,i8,1x,i2,1x,i3)') '<< CRYST >> updateState at el ip g ',e,i,g - write(6,*) - write(6,'(a,/,(12x,12(f12.1,1x)))') '<< CRYST >> absolute residuum tolerance', & - stateResiduum(1:mySizeDotState,g,i,e) / constitutive_aTolState(g,i,e)%p(1:mySizeDotState) - write(6,*) - write(6,'(a,/,(12x,12(f12.1,1x)))') '<< CRYST >> relative residuum tolerance', & - relStateResiduum(1:mySizeDotState,g,i,e) / rTol_crystalliteState - write(6,*) - write(6,'(a,/,(12x,12(e12.5,1x)))') '<< CRYST >> dotState', constitutive_dotState(g,i,e)%p(1:mySizeDotState) & - - 2.0_pReal * stateResiduum(1:mySizeDotState,g,i,e) / crystallite_subdt(g,i,e) ! calculate former dotstate from higher order solution and state residuum - write(6,*) - write(6,'(a,/,(12x,12(e12.5,1x)))') '<< CRYST >> new state', constitutive_state(g,i,e)%p(1:mySizeDotState) - write(6,*) - endif + if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt & + .and. ((e == debug_e .and. i == debug_i .and. g == debug_g)& + .or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0_pInt)) then + write(6,'(a,i8,1x,i2,1x,i3)') '<< CRYST >> updateState at el ip g ',e,i,g + write(6,*) + write(6,'(a,/,(12x,12(f12.1,1x)))') '<< CRYST >> absolute residuum tolerance', & + stateResiduum(1:mySizeDotState,g,i,e) / constitutive_aTolState(g,i,e)%p(1:mySizeDotState) + write(6,*) + write(6,'(a,/,(12x,12(f12.1,1x)))') '<< CRYST >> relative residuum tolerance', & + relStateResiduum(1:mySizeDotState,g,i,e) / rTol_crystalliteState + write(6,*) + write(6,'(a,/,(12x,12(e12.5,1x)))') '<< CRYST >> dotState', constitutive_dotState(g,i,e)%p(1:mySizeDotState) & + - 2.0_pReal * stateResiduum(1:mySizeDotState,g,i,e) / crystallite_subdt(g,i,e) ! calculate former dotstate from higher order solution and state residuum + write(6,*) + write(6,'(a,/,(12x,12(e12.5,1x)))') '<< CRYST >> new state', constitutive_state(g,i,e)%p(1:mySizeDotState) + write(6,*) + endif #endif - - ! --- converged ? --- + + ! --- converged ? --- - if ( all( abs(relStateResiduum(:,g,i,e)) < rTol_crystalliteState & - .or. abs(stateResiduum(1:mySizeDotState,g,i,e)) < constitutive_aTolState(g,i,e)%p(1:mySizeDotState)) & - .and. abs(relTemperatureResiduum(g,i,e)) < rTol_crystalliteTemperature ) then - crystallite_converged(g,i,e) = .true. ! ... converged per definitionem + if ( all( abs(relStateResiduum(:,g,i,e)) < rTol_crystalliteState & + .or. abs(stateResiduum(1:mySizeDotState,g,i,e)) < constitutive_aTolState(g,i,e)%p(1:mySizeDotState)) & + .and. abs(relTemperatureResiduum(g,i,e)) < rTol_crystalliteTemperature ) then + crystallite_converged(g,i,e) = .true. ! ... converged per definitionem + if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0_pInt) then + !$OMP CRITICAL (distributionState) + debug_StateLoopDistribution(2,numerics_integrationMode) = & + debug_StateLoopDistribution(2,numerics_integrationMode) + 1_pInt + !$OMP END CRITICAL (distributionState) + endif + endif + + endif + enddo; enddo; enddo + !$OMP ENDDO + +elseif (numerics_integrationMode > 1) then ! stiffness calculation + + !$OMP DO + do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains + if (crystallite_todo(g,i,e)) then + crystallite_converged(g,i,e) = .true. ! ... converged per definitionem if (iand(debug_level(debug_crystallite), debug_levelBasic) /= 0_pInt) then !$OMP CRITICAL (distributionState) debug_StateLoopDistribution(2,numerics_integrationMode) = & @@ -2011,10 +2043,11 @@ relTemperatureResiduum = 0.0_pReal !$OMP END CRITICAL (distributionState) endif endif + enddo; enddo; enddo + !$OMP ENDDO + +endif - endif - enddo; enddo; enddo -!$OMP ENDDO !$OMP END PARALLEL @@ -2099,23 +2132,7 @@ endif !$OMP PARALLEL -if (numerics_integrationMode < 2) then ! in stiffness calculation mode we do not need to do the state integration again, since this is not influenced by a small perturbation in F - - ! --- STATE JUMP --- - - !$OMP DO - do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains - if (crystallite_todo(g,i,e)) then - crystallite_todo(g,i,e) = crystallite_stateJump(g,i,e) - if (.not. crystallite_todo(g,i,e) .and. .not. crystallite_localPlasticity(g,i,e)) then ! if broken non-local... - !$OMP CRITICAL (checkTodo) - crystallite_todo = crystallite_todo .and. crystallite_localPlasticity ! ...all non-locals skipped - !$OMP END CRITICAL (checkTodo) - endif - endif - enddo; enddo; enddo - !$OMP ENDDO - +if (numerics_integrationMode == 1_pInt) then ! --- DOT STATE AND TEMPERATURE --- @@ -2157,7 +2174,7 @@ if (numerics_integrationMode < 2) then + constitutive_dotState(g,i,e)%p(1:mySizeDotState) * crystallite_subdt(g,i,e) crystallite_Temperature(g,i,e) = crystallite_subTemperature0(g,i,e) & + crystallite_dotTemperature(g,i,e) * crystallite_subdt(g,i,e) -#ifndef _OPENMP + #ifndef _OPENMP if (iand(debug_level(debug_crystallite), debug_levelExtensive) /= 0_pInt & .and. ((e == debug_e .and. i == debug_i .and. g == debug_g) & .or. .not. iand(debug_level(debug_crystallite), debug_levelSelective) /= 0_pInt)) then @@ -2168,19 +2185,35 @@ if (numerics_integrationMode < 2) then write(6,'(a,/,(12x,12(e12.6,1x)))') '<< CRYST >> new state', constitutive_state(g,i,e)%p(1:mySizeDotState) write(6,*) endif -#endif + #endif endif enddo; enddo; enddo !$OMP ENDDO + ! --- STATE JUMP --- + + !$OMP DO + do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains + if (crystallite_todo(g,i,e)) then + crystallite_todo(g,i,e) = crystallite_stateJump(g,i,e) + if (.not. crystallite_todo(g,i,e) .and. .not. crystallite_localPlasticity(g,i,e)) then ! if broken non-local... + !$OMP CRITICAL (checkTodo) + crystallite_todo = crystallite_todo .and. crystallite_localPlasticity ! ...all non-locals skipped + !$OMP END CRITICAL (checkTodo) + endif + endif + enddo; enddo; enddo + !$OMP ENDDO + + ! --- UPDATE DEPENDENT STATES --- !$OMP DO do e = eIter(1),eIter(2); do i = iIter(1,e),iIter(2,e); do g = gIter(1,e),gIter(2,e) ! iterate over elements, ips and grains if (crystallite_todo(g,i,e)) then call constitutive_microstructure(crystallite_Temperature(g,i,e), crystallite_Fe(1:3,1:3,g,i,e), & - crystallite_Fp(1:3,1:3,g,i,e), g, i, e) ! update dependent state variables to be consistent with basic states + crystallite_Fp(1:3,1:3,g,i,e), g, i, e) ! update dependent state variables to be consistent with basic states endif enddo; enddo; enddo !$OMP ENDDO