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

add dv4 framework

parent f6452117
Branches
Tags
No related merge requests found
Expand all Hide whitespace changes
Inline Side-by-side
src/auxval.F90
+ 2
− 0
View file @ c082a7a6
......@@ -50,6 +50,8 @@ subroutine auxval
CALL build_dnjs_table ! precompute the Laguerre nonlin product coeffs
ENDIF
CALL build_dv4Hp_table ! precompute the hermite fourth derivative table
!! Display parallel settings
DO i_ = 0,num_procs-1
CALL mpi_barrier(MPI_COMM_WORLD, ierr)
......
src/collision_mod.F90
+ 137
− 236
View file @ c082a7a6
This diff is collapsed.
src/model_mod.F90
+ 10
− 4
View file @ c082a7a6
MODULE model
! Module for diagnostic parameters
USE prec_const
IMPLICIT NONE
PRIVATE
......@@ -11,12 +10,15 @@ MODULE model
CHARACTER(len=16), &
PUBLIC, PROTECTED ::LINEARITY= 'linear' ! To turn on non linear bracket term
LOGICAL, PUBLIC, PROTECTED :: KIN_E = .true. ! Simulate kinetic electron (adiabatic otherwise)
INTEGER, PUBLIC, PROTECTED :: N_HD = 4 ! order of numerical spatial diffusion
REAL(dp), PUBLIC, PROTECTED :: mu_x = 0._dp ! spatial x-Hyperdiffusivity coefficient (for num. stability)
REAL(dp), PUBLIC, PROTECTED :: mu_y = 0._dp ! spatial y-Hyperdiffusivity coefficient (for num. stability)
LOGICAL, PUBLIC, PROTECTED :: HDz_h = .false. ! to apply z-hyperdiffusion on non adiab part
REAL(dp), PUBLIC, PROTECTED :: mu_z = 0._dp ! spatial z-Hyperdiffusivity coefficient (for num. stability)
CHARACTER(len=16), &
PUBLIC, PROTECTED :: HYP_V = 'hypcoll' ! hyperdiffusion model for velocity space ('none','hypcoll','dvpar4')
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)
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 !
......@@ -59,17 +61,21 @@ CONTAINS
SUBROUTINE model_readinputs
! Read the input parameters
USE basic, ONLY : lu_in, my_id
USE basic, ONLY : lu_in, my_id, num_procs_p
USE prec_const
IMPLICIT NONE
NAMELIST /MODEL_PAR/ CLOS, NL_CLOS, KERN, LINEARITY, KIN_E, &
mu_x, mu_y, N_HD, mu_z, mu_p, mu_j, nu,&
mu_x, mu_y, N_HD, HDz_h, mu_z, mu_p, mu_j, HYP_V, nu,&
tau_e, tau_i, sigma_e, sigma_i, q_e, q_i,&
k_Ne, k_Ni, k_Te, k_Ti, k_gB, k_cB, lambdaD, beta
READ(lu_in,model_par)
IF((HYP_V .EQ. 'dvpar4') .AND. (num_procs_p .GT. 1)) THEN
ERROR STOP '>> ERROR << dvpar4 velocity dissipation is not compatible with current p parallelization'
ENDIF
IF(.NOT. KIN_E) THEN
IF(my_id.EQ.0) print*, 'Adiabatic electron model -> beta = 0'
beta = 0._dp
......
src/moments_eq_rhs_mod.F90
+ 25
− 17
View file @ c082a7a6
......@@ -24,7 +24,7 @@ SUBROUTINE compute_moments_eq_rhs
xphij_i, xphijp1_i, xphijm1_i, xpsij_i, xpsijp1_i, xpsijm1_i,&
kernel_i, nadiab_moments_i, ddz_nipj, interp_nipj, Sipj,&
moments_i(ipgs_i:ipge_i,ijgs_i:ijge_i,ikys:ikye,ikxs:ikxe,izgs:izge,updatetlevel),&
TColl_i, ddzND_nipj, diff_pi_coeff, diff_ji_coeff,&
TColl_i, ddzND_Nipj, diff_pi_coeff, diff_ji_coeff,&
moments_rhs_i(ips_i:ipe_i,ijs_i:ije_i,ikys:ikye,ikxs:ikxe,izs:ize,updatetlevel))
!compute ion moments_eq_rhs
......@@ -36,7 +36,7 @@ SUBROUTINE compute_moments_eq_rhs
xphij_e, xphijp1_e, xphijm1_e, xpsij_e, xpsijp1_e, xpsijm1_e,&
kernel_e, nadiab_moments_e, ddz_nepj, interp_nepj, Sepj,&
moments_e(ipgs_e:ipge_e,ijgs_e:ijge_e,ikys:ikye,ikxs:ikxe,izgs:izge,updatetlevel),&
TColl_e, ddzND_nepj, diff_pe_coeff, diff_je_coeff,&
TColl_e, ddzND_Nepj, diff_pe_coeff, diff_je_coeff,&
moments_rhs_e(ips_e:ipe_e,ijs_e:ije_e,ikys:ikye,ikxs:ikxe,izs:ize,updatetlevel))
CONTAINS
......@@ -193,21 +193,30 @@ SUBROUTINE compute_moments_eq_rhs
-mu_y*diff_ky_coeff*ky**N_HD*moments_(ip,ij,iky,ikx,iz) &
! Numerical parallel hyperdiffusion "mu_z*ddz**4" see Pueschel 2010 (eq 25)
-mu_z*diff_dz_coeff*ddzND_napj_(ip,ij,iky,ikx,iz)
! GX like Hermite hypercollisions see Mandell et al. 2023 (eq 3.23), unadvised to use it
IF (p_int .GT. 2) &
moments_rhs_(ip,ij,iky,ikx,iz) = &
moments_rhs_(ip,ij,iky,ikx,iz) - mu_p*diff_pe_coeff*p_int**6*moments_(ip,ij,iky,ikx,iz)
IF (j_int .GT. 1) &
moments_rhs_(ip,ij,iky,ikx,iz) = &
moments_rhs_(ip,ij,iky,ikx,iz) - mu_j*diff_je_coeff*j_int**6*moments_(ip,ij,iky,ikx,iz)
! fourth order numerical diffusion in vpar
! IF( (ip-4 .GT. 0) .AND. (num_procs_p .EQ. 1) ) &
! ! Numerical parallel velocity hyperdiffusion "+ dvpar4 g_a" see Pueschel 2010 (eq 33)
! ! (not used often so not parallelized)
! moments_rhs_(ip,ij,iky,ikx,iz) = &
! moments_rhs_(ip,ij,iky,ikx,iz) &
! + mu_p * moments_(ip-4,ij,iky,ikx,iz)
!! Velocity space dissipation (should be implemented somewhere else)
SELECT CASE(HYP_V)
CASE('hypcoll') ! GX like Hermite hypercollisions see Mandell et al. 2023 (eq 3.23), unadvised to use it
IF (p_int .GT. 2) &
moments_rhs_(ip,ij,iky,ikx,iz) = &
moments_rhs_(ip,ij,iky,ikx,iz) - mu_p*diff_pe_coeff*p_int**6*moments_(ip,ij,iky,ikx,iz)
IF (j_int .GT. 1) &
moments_rhs_(ip,ij,iky,ikx,iz) = &
moments_rhs_(ip,ij,iky,ikx,iz) - mu_j*diff_je_coeff*j_int**6*moments_(ip,ij,iky,ikx,iz)
CASE('dvpar4')
! fourth order numerical diffusion in vpar
IF(ip-4 .GT. 0) &
! Numerical parallel velocity hyperdiffusion "+ dvpar4 g_a" see Pueschel 2010 (eq 33)
! (not used often so not parallelized)
moments_rhs_(ip,ij,iky,ikx,iz) = &
moments_rhs_(ip,ij,iky,ikx,iz) &
+ mu_p*dv4_Hp_coeff(p_int)*moments_(ip-4,ij,iky,ikx,iz)
! + dummy Laguerre diff
IF (j_int .GT. 1) &
moments_rhs_(ip,ij,iky,ikx,iz) = &
moments_rhs_(ip,ij,iky,ikx,iz) - mu_j*diff_je_coeff*j_int**6*moments_(ip,ij,iky,ikx,iz)
CASE DEFAULT
END SELECT
ELSE
moments_rhs_(ip,ij,iky,ikx,iz) = 0._dp
ENDIF
......@@ -216,7 +225,6 @@ SUBROUTINE compute_moments_eq_rhs
END DO kyloop
END DO kxloop
END DO zloop
! Execution time end
CALL cpu_time(t1_rhs)
tc_rhs = tc_rhs + (t1_rhs-t0_rhs)
......
src/processing_mod.F90
+ 14
− 8
View file @ c082a7a6
......@@ -405,11 +405,11 @@ SUBROUTINE compute_nadiab_moments_z_gradients_and_interp
!
USE fields, ONLY : moments_i, moments_e, phi, psi
USE array, ONLY : kernel_e, kernel_i, nadiab_moments_e, nadiab_moments_i, &
ddz_nepj, ddzND_nepj, interp_nepj,&
ddz_nipj, ddzND_nipj, interp_nipj!, ddz_phi
ddz_nepj, ddzND_Nepj, interp_nepj,&
ddz_nipj, ddzND_Nipj, interp_nipj!, ddz_phi
USE time_integration, ONLY : updatetlevel
USE model, ONLY : qe_taue, qi_taui,q_o_sqrt_tau_sigma_e, q_o_sqrt_tau_sigma_i, &
KIN_E, CLOS, beta
KIN_E, CLOS, beta, HDz_h
USE calculus, ONLY : grad_z, grad_z2, grad_z4, interp_z
IMPLICIT NONE
INTEGER :: eo, p_int, j_int
......@@ -489,9 +489,12 @@ SUBROUTINE compute_nadiab_moments_z_gradients_and_interp
! Compute z derivatives
CALL grad_z(eo,nadiab_moments_e(ip,ij,iky,ikx,izgs:izge), ddz_nepj(ip,ij,iky,ikx,izs:ize))
! Parallel hyperdiffusion
CALL grad_z4(moments_e(ip,ij,iky,ikx,izgs:izge,updatetlevel),ddzND_nepj(ip,ij,iky,ikx,izs:ize))
! CALL grad_z2(nadiab_moments_e(ip,ij,iky,ikx,izgs:izge),ddzND_nepj(ip,ij,iky,ikx,izs:ize))
! Compute even odd grids interpolation
IF (HDz_h) THEN
CALL grad_z4(nadiab_moments_e(ip,ij,iky,ikx,izgs:izge),ddzND_Nepj(ip,ij,iky,ikx,izs:ize))
ELSE
CALL grad_z4(moments_e(ip,ij,iky,ikx,izgs:izge,updatetlevel),ddzND_Nepj(ip,ij,iky,ikx,izs:ize))
ENDIF
! Compute even odd grids interpolation
CALL interp_z(eo,nadiab_moments_e(ip,ij,iky,ikx,izgs:izge), interp_nepj(ip,ij,iky,ikx,izs:ize))
ENDDO
ENDDO
......@@ -508,8 +511,11 @@ SUBROUTINE compute_nadiab_moments_z_gradients_and_interp
! Compute z first derivative
CALL grad_z(eo,nadiab_moments_i(ip,ij,iky,ikx,izgs:izge), ddz_nipj(ip,ij,iky,ikx,izs:ize))
! Parallel numerical diffusion
CALL grad_z4(moments_i(ip,ij,iky,ikx,izgs:izge,updatetlevel),ddzND_nipj(ip,ij,iky,ikx,izs:ize))
! CALL grad_z2(nadiab_moments_i(ip,ij,iky,ikx,izgs:izge),ddzND_nipj(ip,ij,iky,ikx,izs:ize))
IF (HDz_h) THEN
CALL grad_z4(nadiab_moments_i(ip,ij,iky,ikx,izgs:izge),ddzND_Nipj(ip,ij,iky,ikx,izs:ize))
ELSE
CALL grad_z4(moments_i(ip,ij,iky,ikx,izgs:izge,updatetlevel),ddzND_Nipj(ip,ij,iky,ikx,izs:ize))
ENDIF
! Compute even odd grids interpolation
CALL interp_z(eo,nadiab_moments_i(ip,ij,iky,ikx,izgs:izge), interp_nipj(ip,ij,iky,ikx,izs:ize))
ENDDO
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment