Skip to content
Snippets Groups Projects
Commit fa51343e authored by Antoine Cyril David Hoffmann's avatar Antoine Cyril David Hoffmann :seedling:
Browse files

Added new parameter to control the order of numerical dissipation

parent 9dcecc98
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
src/grid_mod.F90
+ 4
− 4
View file @ fa51343e
......@@ -280,7 +280,7 @@ CONTAINS
SUBROUTINE set_kygrid
USE prec_const
USE model, ONLY: LINEARITY
USE model, ONLY: LINEARITY, N_HD
IMPLICIT NONE
INTEGER :: i_, in, istart, iend
Nky = Ny/2+1 ! Defined only on positive kx since fields are real
......@@ -293,7 +293,7 @@ CONTAINS
deltaky = 2._dp*PI/Ly
ky_max = Nky*deltaky
ky_min = deltaky
diff_ky_coeff= (1._dp/ky_max)**4
diff_ky_coeff= (1._dp/ky_max)**N_HD
ENDIF
! Build the full grids on process 0 to diagnose it without comm
ALLOCATE(kyarray_full(1:Nky))
......@@ -345,7 +345,7 @@ CONTAINS
SUBROUTINE set_kxgrid(shear)
USE prec_const
USE model, ONLY: LINEARITY
USE model, ONLY: LINEARITY, N_HD
IMPLICIT NONE
REAL(dp), INTENT(IN) :: shear
INTEGER :: i_, counter
......@@ -375,7 +375,7 @@ CONTAINS
deltakx = 2._dp*PI/Lx
kx_max = (Nkx/2)*deltakx
kx_min = deltakx
diff_kx_coeff= (1._dp/kx_max)**4
diff_kx_coeff= (1._dp/kx_max)**N_HD
! Creating a grid ordered as dk*(0 1 2 3 -2 -1)
local_kxmax = 0._dp
DO ikx = ikxs,ikxe
......
src/model_mod.F90
+ 9
− 6
View file @ fa51343e
......@@ -16,10 +16,11 @@ MODULE model
REAL(dp), PUBLIC, PROTECTED :: mu_z = 0._dp ! spatial z-Hyperdiffusivity coefficient (for num. stability)
REAL(dp), PUBLIC, PROTECTED :: mu_p = 0._dp ! kinetic para hyperdiffusivity coefficient (for num. stability)
REAL(dp), PUBLIC, PROTECTED :: mu_j = 0._dp ! kinetic perp hyperdiffusivity coefficient (for num. stability)
REAL(dp), PUBLIC, PROTECTED :: nu = 1._dp ! Collision frequency
INTEGER, PUBLIC, PROTECTED :: N_HD = 4 ! order of numerical spatial diffusion
REAL(dp), PUBLIC, PROTECTED :: nu = 0._dp ! Collision frequency
REAL(dp), PUBLIC, PROTECTED :: tau_e = 1._dp ! Temperature
REAL(dp), PUBLIC, PROTECTED :: tau_i = 1._dp !
REAL(dp), PUBLIC, PROTECTED :: sigma_e = 1._dp ! Mass
REAL(dp), PUBLIC, PROTECTED :: sigma_e = 0.023338_dp! sqrt of electron ion mass ratio
REAL(dp), PUBLIC, PROTECTED :: sigma_i = 1._dp !
REAL(dp), PUBLIC, PROTECTED :: q_e = -1._dp ! Charge
REAL(dp), PUBLIC, PROTECTED :: q_i = 1._dp !
......@@ -57,7 +58,7 @@ CONTAINS
IMPLICIT NONE
NAMELIST /MODEL_PAR/ CLOS, NL_CLOS, KERN, LINEARITY, KIN_E, &
mu_x, mu_y, mu_z, mu_p, mu_j, nu,&
mu_x, mu_y, N_HD, mu_z, mu_p, mu_j, nu,&
tau_e, tau_i, sigma_e, sigma_i, q_e, q_i,&
K_n, K_T, K_E, GradB, CurvB, lambdaD
......@@ -109,9 +110,11 @@ CONTAINS
CALL attach(fidres, TRIM(str), "KIN_E", KIN_E)
CALL attach(fidres, TRIM(str), "nu", nu)
CALL attach(fidres, TRIM(str), "mu", 0)
CALL attach(fidres, TRIM(str), "mu_x", mu_x)
CALL attach(fidres, TRIM(str), "mu_y", mu_y)
CALL attach(fidres, TRIM(str), "mu_z", mu_z)
CALL attach(fidres, TRIM(str), "mu_x", mu_x)
CALL attach(fidres, TRIM(str), "mu_y", mu_y)
CALL attach(fidres, TRIM(str), "mu_z", mu_z)
CALL attach(fidres, TRIM(str), "mu_p", mu_p)
CALL attach(fidres, TRIM(str), "mu_j", mu_j)
CALL attach(fidres, TRIM(str), "tau_e", tau_e)
CALL attach(fidres, TRIM(str), "tau_i", tau_i)
CALL attach(fidres, TRIM(str), "sigma_e", sigma_e)
......
src/moments_eq_rhs_mod.F90
+ 9
− 9
View file @ fa51343e
......@@ -37,7 +37,7 @@ SUBROUTINE moments_eq_rhs_e
COMPLEX(dp) :: Tnepp1j, Tnepm1j, Tnepp1jm1, Tnepm1jm1 ! Terms from mirror force with non adiab moments
COMPLEX(dp) :: Tperp, Tpar, Tmir, Tphi
COMPLEX(dp) :: Unepm1j, Unepm1jp1, Unepm1jm1 ! Terms from mirror force with adiab moments
COMPLEX(dp) :: i_ky,phikykxz
COMPLEX(dp) :: i_kx,i_ky,phikykxz
! Measuring execution time
CALL cpu_time(t0_rhs)
......@@ -46,6 +46,7 @@ SUBROUTINE moments_eq_rhs_e
zloope : DO iz = izs,ize
kxloope : DO ikx = ikxs,ikxe
kx = kxarray(ikx) ! radial wavevector
i_kx = imagu * kx ! toroidal derivative
kyloope : DO iky = ikys,ikye
ky = kyarray(iky) ! toroidal wavevector
......@@ -110,11 +111,10 @@ SUBROUTINE moments_eq_rhs_e
! Drives (density + temperature gradients)
- i_ky * Tphi &
! Numerical perpendicular hyperdiffusion (totally artificial, for stability purpose)
- mu_x*diff_kx_coeff*kx**4*moments_e(ip,ij,iky,ikx,iz,updatetlevel) &
- mu_y*diff_ky_coeff*ky**4*moments_e(ip,ij,iky,ikx,iz,updatetlevel) &
! - (mu_x*kx**2 + mu_y*ky**2)*moments_e(ip,ij,iky,ikx,iz,updatetlevel) &
- mu_x*diff_kx_coeff*kx**N_HD*moments_e(ip,ij,iky,ikx,iz,updatetlevel) &
- mu_y*diff_ky_coeff*ky**N_HD*moments_e(ip,ij,iky,ikx,iz,updatetlevel) &
! Numerical parallel hyperdiffusion "+ (mu_z*kz**4)" see Pueschel 2010 (eq 25)
+ mu_z * diff_dz_coeff * ddz4_Nepj(ip,ij,iky,ikx,iz) &
- mu_z * diff_dz_coeff * ddz4_Nepj(ip,ij,iky,ikx,iz) &
! Collision term
+ TColl_e(ip,ij,iky,ikx,iz) &
! Nonlinear term
......@@ -166,7 +166,7 @@ SUBROUTINE moments_eq_rhs_i
COMPLEX(dp) :: Tnipp1j, Tnipm1j, Tnipp1jm1, Tnipm1jm1 ! Terms from mirror force with non adiab moments
COMPLEX(dp) :: Unipm1j, Unipm1jp1, Unipm1jm1 ! Terms from mirror force with adiab moments
COMPLEX(dp) :: Tperp, Tpar, Tmir, Tphi
COMPLEX(dp) :: i_ky, phikykxz
COMPLEX(dp) :: i_kx, i_ky, phikykxz
! Measuring execution time
CALL cpu_time(t0_rhs)
......@@ -176,7 +176,7 @@ SUBROUTINE moments_eq_rhs_i
zloopi : DO iz = izs,ize
kxloopi : DO ikx = ikxs,ikxe
kx = kxarray(ikx) ! radial wavevector
i_kx = imagu * kx ! toroidal derivative
kyloopi : DO iky = ikys,ikye
ky = kyarray(iky) ! toroidal wavevector
i_ky = imagu * ky ! toroidal derivative
......@@ -242,8 +242,8 @@ SUBROUTINE moments_eq_rhs_i
! Drives (density + temperature gradients)
- i_ky * Tphi &
! Numerical hyperdiffusion (totally artificial, for stability purpose)
- (mu_x*kx**4 + mu_y*ky**4)*moments_i(ip,ij,iky,ikx,iz,updatetlevel) &
! - (mu_x*kx**2 + mu_y*ky**2)*moments_i(ip,ij,iky,ikx,iz,updatetlevel) &
+ mu_x*diff_kx_coeff*i_kx**N_HD*moments_i(ip,ij,iky,ikx,iz,updatetlevel) &
+ mu_y*diff_ky_coeff*i_ky**N_HD*moments_i(ip,ij,iky,ikx,iz,updatetlevel) &
! Numerical parallel hyperdiffusion "+ (mu_z*kz**4)"
+ mu_z * diff_dz_coeff * ddz4_Nipj(ip,ij,iky,ikx,iz) &
! Collision term
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment