From b64b5b82de263fb3eab2a450ae230e13b247869f Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 5 Aug 2018 17:06:03 +0200 Subject: [PATCH 1/3] integrationMode was always 1, very obscure code was most likely the leftover of old functionality --- src/crystallite.f90 | 70 ++++++++++++---------------------------- src/numerics.f90 | 11 +++---- src/plastic_nonlocal.f90 | 17 ++++------ 3 files changed, 31 insertions(+), 67 deletions(-) diff --git a/src/crystallite.f90 b/src/crystallite.f90 index aea4fb993..88782751d 100644 --- a/src/crystallite.f90 +++ b/src/crystallite.f90 @@ -114,6 +114,7 @@ module crystallite end enum integer(kind(undefined_ID)),dimension(:,:), allocatable, private :: & crystallite_outputID !< ID of each post result output + procedure(), pointer :: integrateState public :: & crystallite_init, & @@ -122,6 +123,7 @@ module crystallite crystallite_push33ToRef, & crystallite_postResults private :: & + integrateState, & crystallite_integrateStateFPI, & crystallite_integrateStateEuler, & crystallite_integrateStateAdaptiveEuler, & @@ -149,6 +151,7 @@ subroutine crystallite_init debug_crystallite, & debug_levelBasic use numerics, only: & + numerics_integrator, & worldrank, & usePingPong use math, only: & @@ -269,6 +272,20 @@ subroutine crystallite_init allocate(crystallite_sizePostResult(maxval(crystallite_Noutput), & size(config_crystallite)), source=0_pInt) + select case(numerics_integrator(1)) + case(1_pInt) + integrateState => crystallite_integrateStateFPI + case(2_pInt) + integrateState => crystallite_integrateStateEuler + case(3_pInt) + integrateState => crystallite_integrateStateAdaptiveEuler + case(4_pInt) + integrateState => crystallite_integrateStateRK4 + case(5_pInt) + integrateState => crystallite_integrateStateRKCK45 + end select + + do c = 1_pInt, size(config_crystallite) #if defined(__GFORTRAN__) @@ -494,9 +511,6 @@ subroutine crystallite_stressAndItsTangent(updateJaco) subStepMinCryst, & subStepSizeCryst, & stepIncreaseCryst, & - nCryst, & - numerics_integrator, & - numerics_integrationMode, & numerics_timeSyncing use debug, only: & debug_level, & @@ -648,7 +662,6 @@ subroutine crystallite_stressAndItsTangent(updateJaco) endif singleRun NiterationCrystallite = 0_pInt - numerics_integrationMode = 1_pInt cutbackLooping: do while (any(crystallite_todo(:,startIP:endIP,FEsolving_execELem(1):FEsolving_execElem(2)))) if (iand(debug_level(debug_crystallite),debug_levelExtensive) /= 0_pInt) & @@ -1026,25 +1039,7 @@ subroutine crystallite_stressAndItsTangent(updateJaco) ! --- integrate --- requires fully defined state array (basic + dependent state) - if (any(crystallite_todo)) then - if (iand(debug_level(debug_crystallite),debug_levelExtensive) /= 0_pInt) then - write(6,'(/,a,i3)') '<< CRYST >> using state integrator ',numerics_integrator(numerics_integrationMode) - flush(6) - endif - select case(numerics_integrator(numerics_integrationMode)) - case(1_pInt) - call crystallite_integrateStateFPI() - case(2_pInt) - call crystallite_integrateStateEuler() - case(3_pInt) - call crystallite_integrateStateAdaptiveEuler() - case(4_pInt) - call crystallite_integrateStateRK4() - case(5_pInt) - call crystallite_integrateStateRKCK45() - end select - endif - + if (any(crystallite_todo)) call integrateState() where(.not. crystallite_converged .and. crystallite_subStep > subStepMinCryst) & ! do not try non-converged & fully cutbacked any further crystallite_todo = .true. @@ -1215,8 +1210,6 @@ end subroutine crystallite_stressAndItsTangent subroutine crystallite_integrateStateRK4() use, intrinsic :: & IEEE_arithmetic - use numerics, only: & - numerics_integrationMode use debug, only: & #ifdef DEBUG debug_e, & @@ -1517,8 +1510,7 @@ subroutine crystallite_integrateStateRKCK45() debug_levelExtensive, & debug_levelSelective use numerics, only: & - rTol_crystalliteState, & - numerics_integrationMode + rTol_crystalliteState use FEsolving, only: & FEsolving_execElem, & FEsolving_execIP @@ -2013,8 +2005,7 @@ subroutine crystallite_integrateStateAdaptiveEuler() debug_levelExtensive, & debug_levelSelective use numerics, only: & - rTol_crystalliteState, & - numerics_integrationMode + rTol_crystalliteState use FEsolving, only: & FEsolving_execElem, & FEsolving_execIP @@ -2082,7 +2073,6 @@ subroutine crystallite_integrateStateAdaptiveEuler() sourceStateResiduum = 0.0_pReal relSourceStateResiduum = 0.0_pReal - integrationMode: if (numerics_integrationMode == 1_pInt) then !$OMP PARALLEL ! --- DOT STATE (EULER INTEGRATION) --- @@ -2182,7 +2172,6 @@ subroutine crystallite_integrateStateAdaptiveEuler() enddo; enddo; enddo !$OMP ENDDO !$OMP END PARALLEL - endif integrationMode ! --- STRESS INTEGRATION (EULER INTEGRATION) --- @@ -2202,9 +2191,6 @@ subroutine crystallite_integrateStateAdaptiveEuler() enddo; enddo; enddo !$OMP END PARALLEL DO - - if (numerics_integrationMode == 1_pInt) then - !$OMP PARALLEL ! --- DOT STATE (HEUN METHOD) --- @@ -2323,17 +2309,6 @@ subroutine crystallite_integrateStateAdaptiveEuler() !$OMP ENDDO !$OMP END PARALLEL - elseif (numerics_integrationMode > 1) then ! stiffness calculation - - !$OMP PARALLEL 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 - crystallite_converged(g,i,e) = crystallite_todo(g,i,e) .or. crystallite_converged(g,i,e) ! ... converged per definitionem - enddo; enddo; enddo - !$OMP END PARALLEL DO - - endif - - ! --- NONLOCAL CONVERGENCE CHECK --- @@ -2364,7 +2339,6 @@ subroutine crystallite_integrateStateEuler() debug_levelExtensive, & debug_levelSelective use numerics, only: & - numerics_integrationMode, & numerics_timeSyncing use FEsolving, only: & FEsolving_execElem, & @@ -2411,7 +2385,6 @@ eIter = FEsolving_execElem(1:2) singleRun = (eIter(1) == eIter(2) .and. iIter(1,eIter(1)) == iIter(2,eIter(2))) - if (numerics_integrationMode == 1_pInt) then !$OMP PARALLEL ! --- DOT STATE --- @@ -2516,7 +2489,6 @@ eIter = FEsolving_execElem(1:2) enddo; enddo; enddo !$OMP ENDDO !$OMP END PARALLEL - endif !$OMP PARALLEL @@ -2581,7 +2553,6 @@ subroutine crystallite_integrateStateFPI() debug_levelSelective use numerics, only: & nState, & - numerics_integrationMode, & rTol_crystalliteState use FEsolving, only: & FEsolving_execElem, & @@ -3156,7 +3127,6 @@ logical function crystallite_integrateStress(& aTol_crystalliteStress, & rTol_crystalliteStress, & iJacoLpresiduum, & - numerics_integrationMode, & subStepSizeLp, & subStepSizeLi use debug, only: debug_level, & diff --git a/src/numerics.f90 b/src/numerics.f90 index 8de664248..e73bc8321 100644 --- a/src/numerics.f90 +++ b/src/numerics.f90 @@ -27,9 +27,8 @@ module numerics worldsize = 0_pInt !< MPI worldsize (/=0 for MPI simulations only) integer(4), protected, public :: & DAMASK_NumThreadsInt = 0 !< value stored in environment variable DAMASK_NUM_THREADS, set to zero if no OpenMP directive - integer(pInt), public :: & - numerics_integrationMode = 0_pInt !< integrationMode 1 = central solution; integrationMode 2 = perturbation, Default 0: undefined, is not read from file - integer(pInt), dimension(2) , protected, public :: & +!> ToDo: numerics_integrator in an array for historical reasons, only element 1 is used! + integer(pInt), dimension(2), protected, public :: & numerics_integrator = 1_pInt !< method used for state integration (central & perturbed state), Default 1: fix-point iteration for both states real(pReal), protected, public :: & relevantStrain = 1.0e-7_pReal, & !< strain increment considered significant (used by crystallite to determine whether strain inc is considered significant) @@ -317,9 +316,7 @@ subroutine numerics_init case ('atol_crystallitestress') aTol_crystalliteStress = IO_floatValue(line,chunkPos,2_pInt) case ('integrator') - numerics_integrator(1) = IO_intValue(line,chunkPos,2_pInt) - case ('integratorstiffness') - numerics_integrator(2) = IO_intValue(line,chunkPos,2_pInt) + numerics_integrator = IO_intValue(line,chunkPos,2_pInt) case ('usepingpong') usepingpong = IO_intValue(line,chunkPos,2_pInt) > 0_pInt case ('timesyncing') @@ -531,7 +528,7 @@ subroutine numerics_init write(6,'(a24,1x,es8.1)') ' rTol_crystalliteState: ',rTol_crystalliteState write(6,'(a24,1x,es8.1)') ' rTol_crystalliteStress: ',rTol_crystalliteStress write(6,'(a24,1x,es8.1)') ' aTol_crystalliteStress: ',aTol_crystalliteStress - write(6,'(a24,2(1x,i8))') ' integrator: ',numerics_integrator + write(6,'(a24,1x,es8.1)') ' integtator: ',numerics_integrator write(6,'(a24,1x,L8)') ' timeSyncing: ',numerics_timeSyncing write(6,'(a24,1x,L8)') ' use ping pong scheme: ',usepingpong write(6,'(a24,1x,es8.1,/)')' unitlength: ',numerics_unitlength diff --git a/src/plastic_nonlocal.f90 b/src/plastic_nonlocal.f90 index 41666a34c..e1355da8f 100644 --- a/src/plastic_nonlocal.f90 +++ b/src/plastic_nonlocal.f90 @@ -2382,8 +2382,7 @@ use, intrinsic :: & use prec, only: dNeq0, & dNeq, & dEq0 -use numerics, only: numerics_integrationMode, & - numerics_timeSyncing +use numerics, only: numerics_timeSyncing use IO, only: IO_error use debug, only: debug_level, & debug_constitutive, & @@ -2942,14 +2941,12 @@ rhoDot = rhoDotFlux & + rhoDotAthermalAnnihilation & + rhoDotThermalAnnihilation -if (numerics_integrationMode == 1_pInt) then ! save rates for output if in central integration mode - rhoDotFluxOutput(1:ns,1:8,1_pInt,ip,el) = rhoDotFlux(1:ns,1:8) - rhoDotMultiplicationOutput(1:ns,1:2,1_pInt,ip,el) = rhoDotMultiplication(1:ns,[1,3]) - rhoDotSingle2DipoleGlideOutput(1:ns,1:2,1_pInt,ip,el) = rhoDotSingle2DipoleGlide(1:ns,9:10) - rhoDotAthermalAnnihilationOutput(1:ns,1:2,1_pInt,ip,el) = rhoDotAthermalAnnihilation(1:ns,9:10) - rhoDotThermalAnnihilationOutput(1:ns,1:2,1_pInt,ip,el) = rhoDotThermalAnnihilation(1:ns,9:10) - rhoDotEdgeJogsOutput(1:ns,1_pInt,ip,el) = 2.0_pReal * rhoDotThermalAnnihilation(1:ns,1) -endif +rhoDotFluxOutput(1:ns,1:8,1_pInt,ip,el) = rhoDotFlux(1:ns,1:8) +rhoDotMultiplicationOutput(1:ns,1:2,1_pInt,ip,el) = rhoDotMultiplication(1:ns,[1,3]) +rhoDotSingle2DipoleGlideOutput(1:ns,1:2,1_pInt,ip,el) = rhoDotSingle2DipoleGlide(1:ns,9:10) +rhoDotAthermalAnnihilationOutput(1:ns,1:2,1_pInt,ip,el) = rhoDotAthermalAnnihilation(1:ns,9:10) +rhoDotThermalAnnihilationOutput(1:ns,1:2,1_pInt,ip,el) = rhoDotThermalAnnihilation(1:ns,9:10) +rhoDotEdgeJogsOutput(1:ns,1_pInt,ip,el) = 2.0_pReal * rhoDotThermalAnnihilation(1:ns,1) #ifdef DEBUG From 67314fc204f010e955443eb63b524c9c2d680846 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 18 Aug 2018 16:40:13 +0200 Subject: [PATCH 2/3] leftover from reducing numerics_integrator to scalar postponed simplification to scalar because of heavy modification in plastic constitutive laws --- src/numerics.f90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/numerics.f90 b/src/numerics.f90 index 05abc0198..146961b86 100644 --- a/src/numerics.f90 +++ b/src/numerics.f90 @@ -26,7 +26,7 @@ module numerics worldsize = 0_pInt !< MPI worldsize (/=0 for MPI simulations only) integer(4), protected, public :: & DAMASK_NumThreadsInt = 0 !< value stored in environment variable DAMASK_NUM_THREADS, set to zero if no OpenMP directive -!> ToDo: numerics_integrator in an array for historical reasons, only element 1 is used! + !< ToDo: numerics_integrator is an array for historical reasons, only element 1 is used! integer(pInt), dimension(2), protected, public :: & numerics_integrator = 1_pInt !< method used for state integration (central & perturbed state), Default 1: fix-point iteration for both states real(pReal), protected, public :: & @@ -525,7 +525,7 @@ subroutine numerics_init write(6,'(a24,1x,es8.1)') ' rTol_crystalliteState: ',rTol_crystalliteState write(6,'(a24,1x,es8.1)') ' rTol_crystalliteStress: ',rTol_crystalliteStress write(6,'(a24,1x,es8.1)') ' aTol_crystalliteStress: ',aTol_crystalliteStress - write(6,'(a24,1x,es8.1)') ' integtator: ',numerics_integrator + write(6,'(a24,2(1x,i8))') ' integrator: ',numerics_integrator write(6,'(a24,1x,L8)') ' timeSyncing: ',numerics_timeSyncing write(6,'(a24,1x,L8)') ' use ping pong scheme: ',usepingpong write(6,'(a24,1x,es8.1,/)')' unitlength: ',numerics_unitlength From 3b3e0bc06851363062efce0d02619f2fa7d0cdf3 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 18 Aug 2018 17:11:50 +0200 Subject: [PATCH 3/3] forgotten use --- src/crystallite.f90 | 2 -- 1 file changed, 2 deletions(-) diff --git a/src/crystallite.f90 b/src/crystallite.f90 index cdb112a48..88782751d 100644 --- a/src/crystallite.f90 +++ b/src/crystallite.f90 @@ -511,8 +511,6 @@ subroutine crystallite_stressAndItsTangent(updateJaco) subStepMinCryst, & subStepSizeCryst, & stepIncreaseCryst, & - numerics_integrator, & - numerics_integrationMode, & numerics_timeSyncing use debug, only: & debug_level, &