diff --git a/PRIVATE b/PRIVATE index 84c4973a3..4532b2772 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit 84c4973a378814b91f6c3525db76d8afe6bc84b7 +Subproject commit 4532b27728f7316c79d76d44368da19c27a66769 diff --git a/src/numerics.f90 b/src/numerics.f90 index ac6ce9a16..abeaff480 100644 --- a/src/numerics.f90 +++ b/src/numerics.f90 @@ -98,8 +98,7 @@ module numerics &-thermal_snes_type ngmres ', & petsc_options = '' integer(pInt), protected, public :: & - fftw_planner_flag = 32_pInt, & !< conversion of fftw_plan_mode to integer, basically what is usually done in the include file of fftw - divergence_correction = 2_pInt !< correct divergence calculation in fourier space 0: no correction, 1: size scaled to 1, 2: size scaled to Npoints + fftw_planner_flag = 32_pInt !< conversion of fftw_plan_mode to integer, basically what is usually done in the include file of fftw logical, protected, public :: & continueCalculation = .false., & !< false:exit if BVP solver does not converge, true: continue calculation despite BVP solver not converging memory_efficient = .true., & !< for fast execution (pre calculation of gamma_hat), Default .true.: do not precalculate @@ -332,8 +331,6 @@ subroutine numerics_init fftw_plan_mode = IO_lc(IO_stringValue(line,chunkPos,2_pInt)) case ('spectralderivative') spectral_derivative = IO_lc(IO_stringValue(line,chunkPos,2_pInt)) - case ('divergence_correction') - divergence_correction = IO_intValue(line,chunkPos,2_pInt) case ('update_gamma') update_gamma = IO_intValue(line,chunkPos,2_pInt) > 0_pInt case ('petsc_options') @@ -460,7 +457,6 @@ subroutine numerics_init ! spectral parameters #ifdef Grid write(6,'(a24,1x,L8)') ' continueCalculation: ',continueCalculation - write(6,'(a24,1x,i8)') ' divergence_correction: ',divergence_correction write(6,'(a24,1x,a)') ' spectral_derivative: ',trim(spectral_derivative) write(6,'(a24,1x,a)') ' fftw_plan_mode: ',trim(fftw_plan_mode) write(6,'(a24,1x,i8)') ' fftw_planner_flag: ',fftw_planner_flag @@ -538,8 +534,6 @@ subroutine numerics_init if (err_damage_tolabs <= 0.0_pReal) call IO_error(301_pInt,ext_msg='err_damage_tolabs') if (err_damage_tolrel <= 0.0_pReal) call IO_error(301_pInt,ext_msg='err_damage_tolrel') #ifdef Grid - if (divergence_correction < 0_pInt .or. & - divergence_correction > 2_pInt) call IO_error(301_pInt,ext_msg='divergence_correction') if (update_gamma .and. & .not. memory_efficient) call IO_error(error_ID = 847_pInt) if (err_stress_tolrel <= 0.0_pReal) call IO_error(301_pInt,ext_msg='err_stress_tolRel') diff --git a/src/spectral_utilities.f90 b/src/spectral_utilities.f90 index 10b57f0ae..41987f932 100644 --- a/src/spectral_utilities.f90 +++ b/src/spectral_utilities.f90 @@ -105,13 +105,14 @@ module spectral_utilities type, private :: tNumerics logical :: & memory_efficient + integer :: & + divergence_correction !< correct divergence calculation in fourier space 0: no correction, 1: size scaled to 1, 2: size scaled to Npoints real(pReal) :: & spectral_derivative, & fftw_planner_flag, & FFTW_timelimit, & !< timelimit for FFTW plan creation for FFTW, see www.fftw.org petsc_defaultOptions, & - petsc_options, & - divergence_correction + petsc_options end type tNumerics type(tNumerics) :: num ! numerics parameters. Better name? @@ -171,8 +172,7 @@ subroutine utilities_init spectral_derivative, & fftw_planner_flag, & petsc_defaultOptions, & - petsc_options, & - divergence_correction + petsc_options use debug, only: & debug_level, & debug_SPECTRAL, & @@ -247,8 +247,12 @@ subroutine utilities_init write(6,'(/,a,3(i12 ))') ' grid a b c: ', grid write(6,'(a,3(es12.5))') ' size x y z: ', geomSize - num%memory_efficient = config_numerics%getInt ('memory_efficient',defaultVal=1) > 0 - num%FFTW_timelimit = config_numerics%getFloat('fftw_timelimit', defaultVal=-1.0) + num%memory_efficient = config_numerics%getInt ('memory_efficient', defaultVal=1) > 0 + num%FFTW_timelimit = config_numerics%getFloat('fftw_timelimit', defaultVal=-1.0) + num%divergence_correction = config_numerics%getInt ('divergence_correction', defaultVal=2) + + if (num%divergence_correction < 0 .or. num%divergence_correction > 2) & + call IO_error(301_pInt,ext_msg='divergence_correction') select case (spectral_derivative) case ('continuous') @@ -264,12 +268,12 @@ subroutine utilities_init !-------------------------------------------------------------------------------------------------- ! scale dimension to calculate either uncorrected, dimension-independent, or dimension- and ! resolution-independent divergence - if (divergence_correction == 1_pInt) then + if (num%divergence_correction == 1) then do j = 1_pInt, 3_pInt if (j /= minloc(geomSize,1) .and. j /= maxloc(geomSize,1)) & scaledGeomSize = geomSize/geomSize(j) enddo - elseif (divergence_correction == 2_pInt) then + elseif (num%divergence_correction == 2) then do j = 1_pInt, 3_pInt if ( j /= int(minloc(geomSize/real(grid,pReal),1),pInt) & .and. j /= int(maxloc(geomSize/real(grid,pReal),1),pInt)) & @@ -362,7 +366,7 @@ subroutine utilities_init if(j > grid(2)/2_pInt + 1_pInt) k_s(2) = k_s(2) - grid(2) ! running from 0,1,...,N/2,N/2+1,-N/2,-N/2+1,...,-1 do i = 1_pInt, grid1Red k_s(1) = i - 1_pInt ! symmetry, junst running from 0,1,...,N/2,N/2+1 - xi2nd(1:3,i,j,k-grid3Offset) = utilities_getFreqDerivative(k_s) ! if divergence_correction is set, frequencies are calculated on unit length + xi2nd(1:3,i,j,k-grid3Offset) = utilities_getFreqDerivative(k_s) where(mod(grid,2)==0 .and. [i,j,k] == grid/2+1 .and. & spectral_derivative_ID == DERIVATIVE_CONTINUOUS_ID) ! for even grids, set the Nyquist Freq component to 0.0 xi1st(1:3,i,j,k-grid3Offset) = cmplx(0.0_pReal,0.0_pReal,pReal)