diff --git a/src/auxval.F90 b/src/auxval.F90
index 136e4bc0980540b783a4d6b78691374e51cb394c..1a043a312da7893b2dc9b42d44719747599d03ce 100644
--- a/src/auxval.F90
+++ b/src/auxval.F90
@@ -16,36 +16,34 @@ subroutine auxval
USE DLRA, ONLY: init_DLRA
#endif
IMPLICIT NONE
-
INTEGER :: i_, ierr
- IF (my_id .EQ. 0) WRITE(*,*) '=== Set auxiliary values ==='
-
+ CALL speak('=== Set auxiliary values ===')
! Init the grids
- CALL set_grids(shear,Npol,LINEARITY,N_HD,EM,Na) ! radial modes (MPI distributed by FFTW)
-
- CALL memory ! Allocate memory for global arrays
-
+ CALL set_grids(shear,Npol,LINEARITY,N_HD,EM,Na)
+ ! Allocate memory for global arrays
+ CALL memory
+ ! Initialize displacement and receive arrays
CALL init_parallel_var(local_np,total_np,local_nky,total_nky,local_nz)
-
+ ! Initialize heatflux variables
CALL init_process
-
- CALL eval_magnetic_geometry ! precompute coeff for lin equation
-
- CALL compute_lin_coeff ! precompute coeff for lin equation and geometry
-
- CALL evaluate_kernels ! precompute the kernels
-
- CALL evaluate_EM_op ! compute inverse of poisson and ampere operators
-
- IF ( LINEARITY .NE. 'linear' ) THEN;
- CALL build_dnjs_table ! precompute the Laguerre nonlin product coeffs
- ENDIF
-
- CALL build_dv4Hp_table ! precompute the hermite fourth derivative table
-
- CALL set_closure_model ! set the closure scheme in use
+ ! precompute coeff for lin equation
+ CALL eval_magnetic_geometry
+ ! precompute coeff for lin equation and geometry
+ CALL compute_lin_coeff
+ ! precompute the kernels
+ CALL evaluate_kernels
+ ! compute inverse of poisson and ampere operators
+ CALL evaluate_EM_op
+ ! precompute the Laguerre nonlin product coeffs
+ IF ( LINEARITY .NE. 'linear' ) &
+ CALL build_dnjs_table
+ ! precompute the hermite fourth derivative table
+ CALL build_dv4Hp_table
+ ! set the closure scheme in use
+ CALL set_closure_model
#ifdef TEST_SVD
+ ! If we want to test SVD decomposition etc.
CALL init_DLRA(local_nky,local_np*local_nj)
#endif
diff --git a/src/basic_mod.F90 b/src/basic_mod.F90
index 2aa6dfd8f5b26cec7d3e1382068ff4b522bbc4d6..7119944d027eeffb748ba82090298dfd895d5336 100644
--- a/src/basic_mod.F90
+++ b/src/basic_mod.F90
@@ -4,7 +4,6 @@ MODULE basic
use prec_const, ONLY : xp
IMPLICIT none
PRIVATE
- ! INCLUDE 'fftw3-mpi.f03'
! INPUT PARAMETERS
INTEGER, PUBLIC, PROTECTED :: nrun = 1 ! Number of time steps to run
real(xp), PUBLIC, PROTECTED :: tmax = 100000.0 ! Maximum simulation time
@@ -19,37 +18,35 @@ MODULE basic
INTEGER, PUBLIC, PROTECTED :: cstep = 0 ! Current step number (Init from restart file)
LOGICAL, PUBLIC :: nlend = .FALSE. ! Signal end of run
LOGICAL, PUBLIC :: crashed = .FALSE. ! Signal end of crashed run
-
INTEGER, PUBLIC :: iframe0d ! counting the number of times 0d datasets are outputed (for diagnose)
INTEGER, PUBLIC :: iframe1d ! counting the number of times 1d datasets are outputed (for diagnose)
INTEGER, PUBLIC :: iframe2d ! counting the number of times 2d datasets are outputed (for diagnose)
INTEGER, PUBLIC :: iframe3d ! counting the number of times 3d datasets are outputed (for diagnose)
INTEGER, PUBLIC :: iframe5d ! counting the number of times 5d datasets are outputed (for diagnose)
-
! List of logical file units
INTEGER, PUBLIC, PROTECTED :: lu_in = 90 ! File duplicated from STDIN
INTEGER, PUBLIC, PROTECTED :: lu_stop = 91 ! stop file, see subroutine TESEND
-
! To measure computation time
type :: chrono
real(xp) :: tstart !start of the chrono
real(xp) :: tstop !stop
real(xp) :: ttot !cumulative time
end type chrono
-
+ ! Define the chronos for each relevant routines
type(chrono), PUBLIC, PROTECTED :: chrono_runt, chrono_mrhs, chrono_advf, chrono_pois, chrono_sapj,&
chrono_diag, chrono_chck, chrono_step, chrono_clos, chrono_ghst, chrono_coll, chrono_napj, chrono_grad
-
#ifdef TEST_SVD
+ ! A chrono for SVD tests
type(chrono), PUBLIC, PROTECTED :: chrono_DLRA
#endif
-
+ ! This sets if the outputs is done through a large gather or using parallelization from futils
+ ! it is recommended to set it to .true.
LOGICAL, PUBLIC, PROTECTED :: GATHERV_OUTPUT = .true.
-
+ ! Routines interfaces
PUBLIC :: allocate_array, basic_outputinputs,basic_data,&
speak, str, increase_step, increase_cstep, increase_time, display_h_min_s,&
set_basic_cp, daytim, start_chrono, stop_chrono
-
+ ! Interface for allocating arrays, these routines allocate and initialize directly to zero
INTERFACE allocate_array
MODULE PROCEDURE allocate_array_xp1,allocate_array_xp2,allocate_array_xp3, &
allocate_array_xp4, allocate_array_xp5, allocate_array_xp6, allocate_array_xp7
diff --git a/src/closure_mod.F90 b/src/closure_mod.F90
index 7baace2eecd84016bf6d25c2d7de65e2549c4eb3..32e2c1374aa0670bc0a6ae7fd9b7310aba12e32f 100644
--- a/src/closure_mod.F90
+++ b/src/closure_mod.F90
@@ -22,6 +22,7 @@ SUBROUTINE closure_readinputs
END SUBROUTINE
SUBROUTINE set_closure_model
+ USE basic, ONLY: speak
USE grid, ONLY: local_np, ngp, local_nj, ngj, parray, jarray,&
pmax, jmax
IMPLICIT NONE
@@ -57,7 +58,10 @@ SUBROUTINE set_closure_model
SELECT CASE(nonlinear_closure)
CASE('truncation')
IF(nmax .LT. 0) THEN
- ERROR STOP "cannot truncate the sum with a number smaller than 0"
+ CALL speak("Set nonlinear truncation to anti Laguerre aliasing")
+ DO ij = 1,local_nj
+ nmaxarray(ij) = jmax - jarray(ij+ngj/2)
+ ENDDO
ELSE
nmaxarray(:) = nmax
ENDIF
@@ -96,83 +100,6 @@ SUBROUTINE apply_closure_model
END SELECT
END SUBROUTINE apply_closure_model
-! ! Positive Oob indices are approximated with a model
-! SUBROUTINE apply_closure_model
-! USE prec_const, ONLY: xp
-! USE grid, ONLY: local_nj,ngj, jarray,&
-! local_np,ngp, parray
-! USE fields, ONLY: moments
-! USE time_integration, ONLY: updatetlevel
-! IMPLICIT NONE
-! INTEGER ::ij,ip,ia
-! SELECT CASE (hierarchy_closure)
-! CASE('truncation')
-! ! zero truncation, An+1=0 for n+1>nmax only
-! CALL ghosts_upper_truncation
-! CASE('max_degree')
-! ! truncation at highest fully represented kinetic moment
-! ! e.g. Dmax = 3 means
-! ! only Napj s.t. p+2j <= 3 are evolved
-! ! -> (p,j) allowed are (0,0),(1,0),(0,1),(2,0),(1,1),(3,0)
-! ! =>> Dmax = min(Pmax,2*Jmax+1)
-! j: DO ij = 1,local_nj+ngj
-! p: DO ip = 1,local_np+ngp
-! IF ( parray(ip)+2*jarray(ij) .GT. dmax) THEN
-! moments(ia,ip,ij,:,:,:,updatetlevel) = 0._xp
-! ENDIF
-! ENDDO p
-! ENDDO j
-! END SELECT
-! CALL ghosts_lower_truncation
-! END SUBROUTINE apply_closure_model
-
-! ! Positive Oob indices are approximated with a model
-! SUBROUTINE ghosts_upper_truncation
-! USE prec_const, ONLY: xp
-! USE grid, ONLY: local_np,ngp,local_pmax, pmax,&
-! local_nj,ngj,local_jmax, jmax
-! USE fields, ONLY: moments
-! USE time_integration, ONLY: updatetlevel
-! IMPLICIT NONE
-! INTEGER ::ig
-! ! zero truncation, An+1=0 for n+1>nmax
-! ! applies only for the processes that evolve the highest moment degree
-! IF(local_jmax .GE. Jmax) THEN
-! DO ig = 1,ngj/2
-! moments(:,:,local_nj+ngj/2+ig,:,:,:,updatetlevel) = 0._xp
-! ENDDO
-! ENDIF
-! ! applies only for the process that has largest p index
-! IF(local_pmax .GE. Pmax) THEN
-! DO ig = 1,ngp/2
-! moments(:,local_np+ngp/2+ig,:,:,:,:,updatetlevel) = 0._xp
-! ENDDO
-! ENDIF
-! END SUBROUTINE ghosts_upper_truncation
-
-! ! Negative OoB indices are 0
-! SUBROUTINE ghosts_lower_truncation
-! USE prec_const, ONLY: xp
-! USE grid, ONLY: ngp,ngj,local_pmin,local_jmin
-! USE fields, ONLY: moments
-! USE time_integration, ONLY: updatetlevel
-! IMPLICIT NONE
-! INTEGER :: ig
-! ! zero truncation, An=0 for n<0
-! IF(local_jmin .EQ. 0) THEN
-! DO ig = 1,ngj/2
-! moments(:,:,ig,:,:,:,updatetlevel) = 0._xp
-! ENDDO
-! ENDIF
-! ! applies only for the process that has lowest p index
-! IF(local_pmin .EQ. 0) THEN
-! DO ig = 1,ngp/2
-! moments(:,ig,:,:,:,:,updatetlevel) = 0._xp
-! ENDDO
-! ENDIF
-
-! END SUBROUTINE ghosts_lower_truncation
-
SUBROUTINE closure_outputinputs(fid)
! Write the input parameters to the results_xx.h5 file
diff --git a/src/grid_mod.F90 b/src/grid_mod.F90
index b759c4e64ca32a31576fafd4435baeeea8dd8cec..4817c7b84040025ce0aa01558defc117478a6329 100644
--- a/src/grid_mod.F90
+++ b/src/grid_mod.F90
@@ -10,7 +10,7 @@ MODULE grid
! GRID Input
INTEGER, PUBLIC, PROTECTED :: pmax = 1 ! The maximal Hermite-moment computed
INTEGER, PUBLIC, PROTECTED :: jmax = 1 ! The maximal Laguerre-moment computed
- INTEGER, PUBLIC, PROTECTED :: maxj = 1 ! The maximal Laguerre-moment
+ INTEGER, PUBLIC, PROTECTED :: maxj = 1 ! The maximal Laguerre-moment
INTEGER, PUBLIC, PROTECTED :: dmax = 1 ! The maximal full GF set of i-moments v^dmax
INTEGER, PUBLIC, PROTECTED :: Nx = 4 ! Number of total internal grid points in x
REAL(xp), PUBLIC, PROTECTED :: Lx = 120_xp ! horizontal length of the spatial box
@@ -563,12 +563,6 @@ CONTAINS
! Find local extrema
local_zmax(eo) = zarray(local_nz+ngz/2,eo)
local_zmin(eo) = zarray(1+ngz/2,eo)
- ! Fill the ghosts
- ! Continue angles
- ! DO ig = 1,ngz/2
- ! zarray(ig,eo) = local_zmin(eo)-REAL(ngz/2-(ig-1),xp)*deltaz
- ! zarray(local_nz+ngz/2+ig,eo) = local_zmax(eo)+REAL(ig,xp)*deltaz
- ! ENDDO
! Periodic z \in (-pi pi-dz)
DO ig = 1,ngz/2 ! first ghost cells
iglob = ig+local_nz_offset-ngz/2
diff --git a/src/parallel_mod.F90 b/src/parallel_mod.F90
index b4e398b6e7e13b67ae36524cd8b29279bd02166d..0739c3a2a9ab5b8db638e57f9fbc92c5d384aa5e 100644
--- a/src/parallel_mod.F90
+++ b/src/parallel_mod.F90
@@ -107,73 +107,11 @@ CONTAINS
CALL MPI_CART_SHIFT(comm0, 2, 1, nbr_D, nbr_U, ierr) !down up neighbours
END SUBROUTINE ppinit
+ ! Here we intialize arrays that describes the data parallelization for each processes. These are used then in gather calls
SUBROUTINE init_parallel_var(np_loc,np_tot,nky_loc,nky_tot,nz_loc)
IMPLICIT NONE
INTEGER, INTENT(IN) :: np_loc,np_tot,nky_loc,nky_tot,nz_loc
INTEGER :: i_
- ! !! P reduction at constant x,y,z,j
- ! ALLOCATE(rcv_p(0:num_procs_p-1),dsp_p(0:num_procs_p-1)) !Displacement sizes for balance diagnostic
- ! ! all processes share their local number of points
- ! CALL MPI_ALLGATHER(np_loc,1,MPI_INTEGER,rcv_p,1,MPI_INTEGER,comm_p,ierr)
- ! ! the displacement array can be build from each np_loc as
- ! dsp_p(0)=0
- ! DO i_=1,num_procs_p-1
- ! dsp_p(i_) =dsp_p(i_-1) + rcv_p(i_-1)
- ! END DO
- ! !!!!!! XYZ gather variables
- ! !! Y reduction at constant x and z
- ! ! number of points recieved and displacement for the y reduction
- ! ALLOCATE(rcv_y(0:num_procs_ky-1),dsp_y(0:num_procs_ky-1)) !Displacement sizes for balance diagnostic
- ! ! all processes share their local number of points
- ! CALL MPI_ALLGATHER(nky_loc,1,MPI_INTEGER,rcv_y,1,MPI_INTEGER,comm_ky,ierr)
- ! ! the displacement array can be build from each np_loc as
- ! dsp_y(0)=0
- ! DO i_=1,num_procs_ky-1
- ! dsp_y(i_) =dsp_y(i_-1) + rcv_y(i_-1)
- ! END DO
- ! !! Z reduction for full slices of y data but constant x
- ! ! number of points recieved and displacement for the z reduction
- ! ALLOCATE(rcv_zy(0:num_procs_z-1),dsp_zy(0:num_procs_z-1)) !Displacement sizes for balance diagnostic
- ! ! all processes share their local number of points
- ! CALL MPI_ALLGATHER(nz_loc*nky_tot,1,MPI_INTEGER,rcv_zy,1,MPI_INTEGER,comm_z,ierr)
- ! ! the displacement array can be build from each np_loc as
- ! dsp_zy(0)=0
- ! DO i_=1,num_procs_z-1
- ! dsp_zy(i_) =dsp_zy(i_-1) + rcv_zy(i_-1)
- ! END DO
- ! !!!!! PJZ gather variables
- ! !! P reduction at constant j and z is already done in module GRID
- ! !! Z reduction for full slices of p data but constant j
- ! ! number of points recieved and displacement for the z reduction
- ! ALLOCATE(rcv_zp(0:num_procs_z-1),dsp_zp(0:num_procs_z-1)) !Displacement sizes for balance diagnostic
- ! ! all processes share their local number of points
- ! CALL MPI_ALLGATHER(nz_loc*np_tot,1,MPI_INTEGER,rcv_zp,1,MPI_INTEGER,comm_z,ierr)
- ! ! the displacement array can be build from each np_loc as
- ! dsp_zp(0)=0
- ! DO i_=1,num_procs_z-1
- ! dsp_zp(i_) =dsp_zp(i_-1) + rcv_zp(i_-1)
- ! END DO
- ! !!!!! PJXYZ gather variables
- ! !! Y reduction for full slices of p data but constant j
- ! ! number of points recieved and displacement for the y reduction
- ! ALLOCATE(rcv_yp(0:num_procs_ky-1),dsp_yp(0:num_procs_ky-1)) !Displacement sizes for balance diagnostic
- ! ! all processes share their local number of points
- ! CALL MPI_ALLGATHER(nky_loc*np_tot,1,MPI_INTEGER,rcv_yp,1,MPI_INTEGER,comm_ky,ierr)
- ! ! the displacement array can be build from each np_loc as
- ! dsp_yp(0)=0
- ! DO i_=1,num_procs_ky-1
- ! dsp_yp(i_) =dsp_yp(i_-1) + rcv_yp(i_-1)
- ! END DO
- ! !! Z reduction for full slices of py data but constant j
- ! ! number of points recieved and displacement for the z reduction
- ! ALLOCATE(rcv_zyp(0:num_procs_z-1),dsp_zyp(0:num_procs_z-1)) !Displacement sizes for balance diagnostic
- ! ! all processes share their local number of points
- ! CALL MPI_ALLGATHER(nz_loc*np_tot*nky_tot,1,MPI_INTEGER,rcv_zyp,1,MPI_INTEGER,comm_z,ierr)
- ! ! the displacement array can be build from each np_loc as
- ! dsp_zyp(0)=0
- ! DO i_=1,num_procs_z-1
- ! dsp_zyp(i_) =dsp_zyp(i_-1) + rcv_zyp(i_-1)
- ! END DO
! P reduction at constant x,y,z,j
ALLOCATE(rcv_p(num_procs_p),dsp_p(num_procs_p)) !Displacement sizes for balance diagnostic
CALL MPI_ALLGATHER(np_loc,1,MPI_INTEGER,rcv_p,1,MPI_INTEGER,comm_p,ierr)
@@ -316,7 +254,6 @@ CONTAINS
END SUBROUTINE gather_pjz
!!!!! Gather a field in kinetic + spatial coordinates on rank 0 !!!!!
- !!!!! Gather a field in spatial coordinates on rank 0 !!!!!
SUBROUTINE gather_pjxyz(field_loc,field_tot,na_tot,np_loc,np_tot,nj_tot,nky_loc,nky_tot,nkx_tot,nz_loc,nz_tot)
IMPLICIT NONE
INTEGER, INTENT(IN) :: na_tot,np_loc,np_tot,nj_tot,nky_loc,nky_tot,nkx_tot,nz_loc,nz_tot