Permalink
Browse files

This update introduces the routine factorize in the initialization

procedure of the code. This routine is used to check the input of the
horizontal grid resolution and checks whether the input grid resolution
can be factorized in small prime factors. If not the simulation is
aborted. For now the threshold is set at 7. In the main a line is added
to show what the shown output values represent (removed double printing
of velocity divergence).
  • Loading branch information...
stevensrjam committed Jan 7, 2016
1 parent 8320ebc commit 73672a673d31615aaeff71807ce8e0fcd1895199
Showing with 58 additions and 3 deletions.
  1. +2 −2 Makefile.am
  2. +2 −0 QuitRoutine.F90
  3. +34 −0 factorize.F90
  4. +20 −1 main.F90
View
4 Makefile.am 100644 → 100755
@@ -10,7 +10,7 @@ CheckDivergence.F90 ExplicitTermsVY.F90 ImplicitAndUpdateVX.F90
CorrectPressure.F90 ExplicitTermsVZ.F90 ImplicitAndUpdateVY.F90 StatReadReduceWrite.F90 \
CorrectVelocity.F90 ImplicitAndUpdateVZ.F90 param.F90 StatRoutines.F90 \
CreateGrid.F90 QuitRoutine.F90 TimeMarcher.F90 \
CreateInitialConditions.F90 GlobalQuantities.F90 InitTimeMarchScheme.F90 ReadInputFile.F90 ReadFlowField.F90 HdfReadContinua.F90
CreateInitialConditions.F90 GlobalQuantities.F90 InitTimeMarchScheme.F90 ReadInputFile.F90 factorize.F90 ReadFlowField.F90 HdfReadContinua.F90
afid_LDADD = $(FFTW3_LIBS) $(LAPACK_LIBS) $(BLAS_LIBS) $(HDF5_LIBS)
@@ -77,4 +77,4 @@ transpose_y_to_x.o: transpose_y_to_x.F90
ImplicitAndUpdateTemp.o: ImplicitAndUpdateTemp.F90 param.o
decomp_2d.o: decomp_2d.F90 transpose_z_to_x.F90 transpose_x_to_z.F90 transpose_x_to_y.F90 transpose_y_to_x.F90 transpose_y_to_z.F90 transpose_z_to_y.F90
interp.o: interp.F90 param.o
factorize.o: factorize.F90
View
@@ -67,6 +67,8 @@ subroutine NotifyError(errorcode)
else if(errorcode.eq.334) then
write(*,*) "walltime greater than walltimemax"
write(*,*) "statistics and continuation updated"
else if(errorcode.eq.444) then
write(*,*) "FFT size in ny or nz is not efficient"
else
write(*,*) "Maximum number of timesteps reached"
end if
View
@@ -0,0 +1,34 @@
! ---------------------------------------------------------------
! This routine factorizes the grid resolution in the horizontal
! directions to make sure that efficient FFT's are used. The input
! is (In2) is "ny or ny" and the output (Out2) is the largest factor
! A loop in the main routine checks whether the largest factor is
! smaller than 7. If not the simulation is aborted.
! ---------------------------------------------------------------
subroutine Factorize(In2,Out2)
implicit none
integer,intent(in) :: In2 ! Horizontal grid resolution
integer,intent(out) :: Out2 ! Largest factor
integer :: Input,Divisor
Input=In2
! Remove all factors of 2
do
if (mod(Input,2) /= 0 .or. Input == 1) exit
Input = Input / 2 ! remove this factor
enddo
! consider odd factors
Divisor = 3
do ! try 3, 5, 7, etc
if (Divisor > Input) exit ! if a factor is too large, exit and done
do ! try this factor repeatedly
if (mod(Input,Divisor) /= 0 .or. Input == 1) exit
Input = Input / Divisor ! remove this factor
enddo
Divisor = Divisor + 2 ! Next odd number
enddo
Out2=Divisor ! Largest factor. Should be small for effcient FFT
return
end
View
@@ -13,6 +13,7 @@ program AFiD
real :: ti(2), tin(3), minwtdt
real :: ts
integer :: prow=0,pcol=0
integer :: lfactor,lfactor2
character(100) :: arg
!*******************************************************
@@ -151,6 +152,18 @@ program AFiD
! ********* starts the time dependent calculation ***
errorcode = 0 !EP set errocode to 0 (OK)
minwtdt = huge(0.0d0) !EP initialize minimum time step walltime
! Check input for efficient FFT
! factorize input FFT directions. The largest factor should
! be relatively small to have efficient FFT's
lfactor=2 ! initialize
call Factorize(nym,lfactor2) ! check nym
lfactor=max(lfactor,lfactor2)
call Factorize(nzm,lfactor2)
lfactor=max(lfactor,lfactor2)
! if largest factor larger than 7 quit the simulation
if (lfactor>7) errorcode=444
do ntime=0,ntst
ti(1) = MPI_WTIME()
@@ -215,7 +228,9 @@ program AFiD
if(mod(time,tout).lt.dt) then
if(ismaster) then
write(6,*) 'Maximum divergence = ', dmax
write(6,*)ntime,time,vmax(1),vmax(2),vmax(3),dmax,tempm,tempmax,tempmin
write(6,*)'ntime - time - vmax(1) - vmax(2) - vmax(3) -&
tempm - tempmax - tempmin'
write(6,*)ntime,time,vmax(1),vmax(2),vmax(3),tempm,tempmax,tempmin
write(6,'(a,f8.3,a)') 'Minimum Iteration Time = ', minwtdt, &
' sec.'
endif
@@ -247,6 +262,10 @@ program AFiD
!EP walltime exceeded walltimemax, no error; normal quit
if(errorcode.eq.334) call QuitRoutine(tin,.true.,errorcode)
!RS FFT input not correct
if(errorcode.eq.444) call QuitRoutine(tin,.true.,errorcode)
errorcode = 100 !EP already finalized
exit

0 comments on commit 73672a6

Please sign in to comment.