Skip to content
Snippets Groups Projects
Commit 358d3684 authored by Antoine Cyril David Hoffmann's avatar Antoine Cyril David Hoffmann
Browse files

Kernel are now z-dependant

parent fba92bd9
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
src/closure_mod.F90
+ 8
− 8
View file @ 358d3684
......@@ -32,8 +32,8 @@ SUBROUTINE apply_closure_model
moments_e(ipeg_e-1,ij,ikx,iky,iz,updatetlevel) = 0._dp
moments_e(ipeg_e ,ij,ikx,iky,iz,updatetlevel) = 0._dp
ENDDO
kernel_e(ijsg_e,ikx,iky) = 0._dp
kernel_e(ijeg_e,ikx,iky) = 0._dp
kernel_e(ijsg_e,ikx,iky,iz) = 0._dp
kernel_e(ijeg_e,ikx,iky,iz) = 0._dp
DO ip = ipsg_i,ipeg_i
moments_i(ip,ijsg_i,ikx,iky,iz,updatetlevel) = 0._dp
......@@ -45,8 +45,8 @@ SUBROUTINE apply_closure_model
moments_i(ipeg_i-1,ij,ikx,iky,iz,updatetlevel) = 0._dp
moments_i(ipeg_i ,ij,ikx,iky,iz,updatetlevel) = 0._dp
ENDDO
kernel_i(ijsg_i,ikx,iky) = 0._dp
kernel_i(ijeg_i,ikx,iky) = 0._dp
kernel_i(ijsg_i,ikx,iky,iz) = 0._dp
kernel_i(ijeg_i,ikx,iky,iz) = 0._dp
ENDDO
ENDDO
ENDDO
......@@ -66,8 +66,8 @@ SUBROUTINE apply_closure_model
moments_e(ipeg_e-1,ij,ikx,iky,iz,updatetlevel) = 0._dp
moments_e(ipeg_e ,ij,ikx,iky,iz,updatetlevel) = 0._dp
ENDDO
kernel_e(ijsg_e,ikx,iky) = 0._dp
kernel_e(ijeg_e,ikx,iky) = 0._dp
kernel_e(ijsg_e,ikx,iky,iz) = 0._dp
kernel_e(ijeg_e,ikx,iky,iz) = 0._dp
DO ip = ipsg_i,ipeg_i
moments_i(ip,ijsg_i,ikx,iky,iz,updatetlevel) = 0._dp
......@@ -79,8 +79,8 @@ SUBROUTINE apply_closure_model
moments_i(ipeg_i-1,ij,ikx,iky,iz,updatetlevel) = 0._dp
moments_i(ipeg_i ,ij,ikx,iky,iz,updatetlevel) = 0._dp
ENDDO
kernel_i(ijsg_i,ikx,iky) = 0._dp
kernel_i(ijeg_i,ikx,iky) = 0._dp
kernel_i(ijsg_i,ikx,iky,iz) = 0._dp
kernel_i(ijeg_i,ikx,iky,iz) = 0._dp
ENDDO
ENDDO
ENDDO
......
src/collision_mod.F90
+ 28
− 28
View file @ 358d3684
......@@ -43,7 +43,7 @@ CONTAINS
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! Non zero term for p = 0 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
IF( p_dp .EQ. 0 ) THEN ! Kronecker p0
! Get adiabatic moment
TColl_ = TColl_ - (p_dp + 2._dp*j_dp + 2._dp*be_2) * qe_taue * Kernel_e(ij_,ikx_,iky_)*phi(ikx_,iky_,iz_)
TColl_ = TColl_ - (p_dp + 2._dp*j_dp + 2._dp*be_2) * qe_taue * Kernel_e(ij_,ikx_,iky_,iz_)*phi(ikx_,iky_,iz_)
!** build required fluid moments **
n_ = 0._dp
upar_ = 0._dp; uperp_ = 0._dp
......@@ -51,9 +51,9 @@ CONTAINS
DO in_ = 1,jmaxe+1
n_dp = REAL(in_-1,dp)
! Store the kernels for sparing readings
Knp0 = Kernel_e(in_,ikx_,iky_)
Knp1 = Kernel_e(in_+1,ikx_,iky_)
Knm1 = Kernel_e(in_-1,ikx_,iky_)
Knp0 = Kernel_e(in_ ,ikx_,iky_,iz_)
Knp1 = Kernel_e(in_+1,ikx_,iky_,iz_)
Knm1 = Kernel_e(in_-1,ikx_,iky_,iz_)
! Nonadiabatic moments (only different from moments when p=0)
nadiab_moment_0j = moments_e(1,in_ ,ikx_,iky_,iz_,updatetlevel) + qe_taue*Knp0*phi(ikx_,iky_,iz_)
! Density
......@@ -67,10 +67,10 @@ CONTAINS
ENDDO
T_ = (Tpar_ + 2._dp*Tperp_)/3._dp - n_
! Add energy restoring term
TColl_ = TColl_ + T_* 4._dp * j_dp * Kernel_e(ij_ ,ikx_,iky_)
TColl_ = TColl_ - T_* 2._dp * (j_dp + 1._dp) * Kernel_e(ij_+1,ikx_,iky_)
TColl_ = TColl_ - T_* 2._dp * j_dp * Kernel_e(ij_-1,ikx_,iky_)
TColl_ = TColl_ + uperp_*be_* (Kernel_e(ij_,ikx_,iky_) - Kernel_e(ij_-1,ikx_,iky_))
TColl_ = TColl_ + T_* 4._dp * j_dp * Kernel_e(ij_ ,ikx_,iky_,iz_)
TColl_ = TColl_ - T_* 2._dp * (j_dp + 1._dp) * Kernel_e(ij_+1,ikx_,iky_,iz_)
TColl_ = TColl_ - T_* 2._dp * j_dp * Kernel_e(ij_-1,ikx_,iky_,iz_)
TColl_ = TColl_ + uperp_*be_* (Kernel_e(ij_,ikx_,iky_,iz_) - Kernel_e(ij_-1,ikx_,iky_,iz_))
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! Non zero term for p = 1 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
......@@ -79,9 +79,9 @@ CONTAINS
upar_ = 0._dp
DO in_ = 1,jmaxe+1
! Parallel velocity
upar_ = upar_ + Kernel_e(in_,ikx_,iky_) * moments_e(2,in_,ikx_,iky_,iz_,updatetlevel)
upar_ = upar_ + Kernel_e(in_,ikx_,iky_,iz_) * moments_e(2,in_,ikx_,iky_,iz_,updatetlevel)
ENDDO
TColl_ = TColl_ + upar_*Kernel_e(ij_,ikx_,iky_)
TColl_ = TColl_ + upar_*Kernel_e(ij_,ikx_,iky_,iz_)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! Non zero term for p = 2 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
......@@ -93,9 +93,9 @@ CONTAINS
DO in_ = 1,jmaxe+1
n_dp = REAL(in_-1,dp)
! Store the kernels for sparing readings
Knp0 = Kernel_e(in_ ,ikx_,iky_)
Knp1 = Kernel_e(in_+1,ikx_,iky_)
Knm1 = Kernel_e(in_-1,ikx_,iky_)
Knp0 = Kernel_e(in_ ,ikx_,iky_,iz_)
Knp1 = Kernel_e(in_+1,ikx_,iky_,iz_)
Knm1 = Kernel_e(in_-1,ikx_,iky_,iz_)
! Nonadiabatic moments (only different from moments when p=0)
nadiab_moment_0j = moments_e(1,in_,ikx_,iky_,iz_,updatetlevel) + qe_taue*Knp0*phi(ikx_,iky_,iz_)
! Density
......@@ -106,7 +106,7 @@ CONTAINS
Tperp_ = Tperp_ + ((2._dp*n_dp+1._dp)*Knp0 - (n_dp+1._dp)*Knp1 - n_dp*Knm1)*nadiab_moment_0j
ENDDO
T_ = (Tpar_ + 2._dp*Tperp_)/3._dp - n_
TColl_ = TColl_ + T_*SQRT2*Kernel_e(ij_,ikx_,iky_)
TColl_ = TColl_ + T_*SQRT2*Kernel_e(ij_,ikx_,iky_,iz_)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
ENDIF
......@@ -149,7 +149,7 @@ CONTAINS
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! Non zero term for p = 0 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
IF( p_dp .EQ. 0 ) THEN ! Kronecker p0
! Get adiabatic moment
TColl_ = TColl_ - (p_dp + 2._dp*j_dp + 2._dp*bi_2) * qi_taui * Kernel_i(ij_,ikx_,iky_)*phi(ikx_,iky_,iz_)
TColl_ = TColl_ - (p_dp + 2._dp*j_dp + 2._dp*bi_2) * qi_taui * Kernel_i(ij_,ikx_,iky_,iz_)*phi(ikx_,iky_,iz_)
!** build required fluid moments **
n_ = 0._dp
upar_ = 0._dp; uperp_ = 0._dp
......@@ -157,9 +157,9 @@ CONTAINS
DO in_ = 1,jmaxi+1
n_dp = REAL(in_-1,dp)
! Store the kernels for sparing readings
Knp0 = Kernel_i(in_,ikx_,iky_)
Knp1 = Kernel_i(in_+1,ikx_,iky_)
Knm1 = Kernel_i(in_-1,ikx_,iky_)
Knp0 = Kernel_i(in_ ,ikx_,iky_,iz_)
Knp1 = Kernel_i(in_+1,ikx_,iky_,iz_)
Knm1 = Kernel_i(in_-1,ikx_,iky_,iz_)
! Nonadiabatic moments (only different from moments when p=0)
nadiab_moment_0j = moments_i(1,in_ ,ikx_,iky_,iz_,updatetlevel) + qi_taui*Knp0*phi(ikx_,iky_,iz_)
! Density
......@@ -173,10 +173,10 @@ CONTAINS
ENDDO
T_ = (Tpar_ + 2._dp*Tperp_)/3._dp - n_
! Add energy restoring term
TColl_ = TColl_ + T_* 4._dp * j_dp * Kernel_i(ij_ ,ikx_,iky_)
TColl_ = TColl_ - T_* 2._dp * (j_dp + 1._dp) * Kernel_i(ij_+1,ikx_,iky_)
TColl_ = TColl_ - T_* 2._dp * j_dp * Kernel_i(ij_-1,ikx_,iky_)
TColl_ = TColl_ + uperp_*bi_* (Kernel_i(ij_,ikx_,iky_) - Kernel_i(ij_-1,ikx_,iky_))
TColl_ = TColl_ + T_* 4._dp * j_dp * Kernel_i(ij_ ,ikx_,iky_,iz_)
TColl_ = TColl_ - T_* 2._dp * (j_dp + 1._dp) * Kernel_i(ij_+1,ikx_,iky_,iz_)
TColl_ = TColl_ - T_* 2._dp * j_dp * Kernel_i(ij_-1,ikx_,iky_,iz_)
TColl_ = TColl_ + uperp_*bi_* (Kernel_i(ij_,ikx_,iky_,iz_) - Kernel_i(ij_-1,ikx_,iky_,iz_))
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! Non zero term for p = 1 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
......@@ -185,9 +185,9 @@ CONTAINS
upar_ = 0._dp
DO in_ = 1,jmaxi+1
! Parallel velocity
upar_ = upar_ + Kernel_i(in_,ikx_,iky_) * moments_i(2,in_,ikx_,iky_,iz_,updatetlevel)
upar_ = upar_ + Kernel_i(in_,ikx_,iky_,iz_) * moments_i(2,in_,ikx_,iky_,iz_,updatetlevel)
ENDDO
TColl_ = TColl_ + upar_*Kernel_i(ij_,ikx_,iky_)
TColl_ = TColl_ + upar_*Kernel_i(ij_,ikx_,iky_,iz_)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! Non zero term for p = 2 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
......@@ -199,9 +199,9 @@ CONTAINS
DO in_ = 1,jmaxi+1
n_dp = REAL(in_-1,dp)
! Store the kernels for sparing readings
Knp0 = Kernel_i(in_ ,ikx_,iky_)
Knp1 = Kernel_i(in_+1,ikx_,iky_)
Knm1 = Kernel_i(in_-1,ikx_,iky_)
Knp0 = Kernel_i(in_ ,ikx_,iky_,iz_)
Knp1 = Kernel_i(in_+1,ikx_,iky_,iz_)
Knm1 = Kernel_i(in_-1,ikx_,iky_,iz_)
! Nonadiabatic moments (only different from moments when p=0)
nadiab_moment_0j = moments_i(1,in_,ikx_,iky_,iz_,updatetlevel) + qi_taui*Knp0*phi(ikx_,iky_,iz_)
! Density
......@@ -212,7 +212,7 @@ CONTAINS
Tperp_ = Tperp_ + ((2._dp*n_dp+1._dp)*Knp0 - (n_dp+1._dp)*Knp1 - n_dp*Knm1)*nadiab_moment_0j
ENDDO
T_ = (Tpar_ + 2._dp*Tperp_)/3._dp - n_
TColl_ = TColl_ + T_*SQRT2*Kernel_i(ij_,ikx_,iky_)
TColl_ = TColl_ + T_*SQRT2*Kernel_i(ij_,ikx_,iky_,iz_)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
ENDIF
......
src/compute_Sapj.F90
+ 2
− 2
View file @ 358d3684
......@@ -54,7 +54,7 @@ zloop: DO iz = izs,ize
kyloope: DO iky = ikys,ikye ! Loop over ky
kx = kxarray(ikx)
ky = kyarray(iky)
kerneln = kernel_e(in, ikx, iky)
kerneln = kernel_e(in, ikx, iky, iz)
! First convolution terms
Fx_cmpx(ikx,iky) = imagu*kx* phi(ikx,iky,iz) * kerneln
......@@ -145,7 +145,7 @@ zloop: DO iz = izs,ize
kyloopi: DO iky = ikys,ikye ! Loop over ky
kx = kxarray(ikx)
ky = kyarray(iky)
kerneln = kernel_i(in, ikx, iky)
kerneln = kernel_i(in, ikx, iky, iz)
! First convolution terms
Fx_cmpx(ikx,iky) = imagu*kx* phi(ikx,iky,iz) * kerneln
......
src/poisson.F90
+ 2
− 2
View file @ 358d3684
......@@ -37,7 +37,7 @@ SUBROUTINE poisson
sum_kernel2_e = 0._dp
! loop over n only if the max polynomial degree
DO ine=1,jmaxe+1 ! ine = n+1
Kne = kernel_e(ine, ikx, iky)
Kne = kernel_e(ine,ikx,iky,iz)
sum_kernel_mom_e = sum_kernel_mom_e + Kne * moments_e(1, ine, ikx, iky, iz, updatetlevel)
sum_kernel2_e = sum_kernel2_e + Kne**2 ! ... sum recursively ...
END DO
......@@ -47,7 +47,7 @@ SUBROUTINE poisson
sum_kernel2_i = 0._dp
! loop over n only if the max polynomial degree
DO ini=1,jmaxi+1
Kni = kernel_i(ini, ikx, iky)
Kni = kernel_i(ini,ikx,iky,iz)
sum_kernel_mom_i = sum_kernel_mom_i + Kni * moments_i(1, ini, ikx, iky, iz, updatetlevel)
sum_kernel2_i = sum_kernel2_i + Kni**2 ! ... sum recursively ...
END DO
......
src/processing_mod.F90
+ 11
− 11
View file @ 358d3684
......@@ -34,7 +34,7 @@ SUBROUTINE compute_radial_ion_transport
imagu * ky_ * moments_i(1,1,ikx,iky,iz,updatetlevel) * CONJG(phi(ikx,iky,iz))
DO ij = ijs_i, ije_i
pflux_local = pflux_local - &
imagu * ky_ * kernel_i(ij,ikx,iky) * moments_i(1,ij,ikx,iky,iz,updatetlevel) * CONJG(phi(ikx,iky,iz))
imagu * ky_ * kernel_i(ij,ikx,iky,iz) * moments_i(1,ij,ikx,iky,iz,updatetlevel) * CONJG(phi(ikx,iky,iz))
ENDDO
ENDDO
ENDDO
......@@ -81,14 +81,14 @@ SUBROUTINE compute_density
! electron density
dens_e(ikx,iky,iz) = 0._dp
DO ij = ijs_e, ije_e
dens_e(ikx,iky,iz) = dens_e(ikx,iky,iz) + kernel_e(ij,ikx,iky) * &
(moments_e(1,ij,ikx,iky,iz,updatetlevel)+q_e/tau_e*kernel_e(ij,ikx,iky)*phi(ikx,iky,iz))
dens_e(ikx,iky,iz) = dens_e(ikx,iky,iz) + kernel_e(ij,ikx,iky,iz) * &
(moments_e(1,ij,ikx,iky,iz,updatetlevel)+q_e/tau_e*kernel_e(ij,ikx,iky,iz)*phi(ikx,iky,iz))
ENDDO
! ion density
dens_i(ikx,iky,iz) = 0._dp
DO ij = ijs_i, ije_i
dens_i(ikx,iky,iz) = dens_i(ikx,iky,iz) + kernel_i(ij,ikx,iky) * &
(moments_i(1,ij,ikx,iky,iz,updatetlevel)+q_i/tau_i*kernel_i(ij,ikx,iky)*phi(ikx,iky,iz))
dens_i(ikx,iky,iz) = dens_i(ikx,iky,iz) + kernel_i(ij,ikx,iky,iz) * &
(moments_i(1,ij,ikx,iky,iz,updatetlevel)+q_i/tau_i*kernel_i(ij,ikx,iky,iz)*phi(ikx,iky,iz))
ENDDO
ENDDO
ENDDO
......@@ -118,9 +118,9 @@ SUBROUTINE compute_temperature
DO ij = ijs_e, ije_e
j_dp = REAL(ij-1,dp)
temp_e(ikx,iky,iz) = temp_e(ikx,iky,iz) + &
2._dp/3._dp * (2._dp*j_dp*kernel_e(ij,ikx,iky) - (j_dp+1)*kernel_e(ij+1,ikx,iky) - j_dp*kernel_e(ij-1,ikx,iky))&
* (moments_e(1,ij,ikx,iky,iz,updatetlevel)+q_e/tau_e*kernel_e(ij,ikx,iky)*phi(ikx,iky,iz)) + &
SQRT2/3._dp * kernel_e(ij,ikx,iky) * moments_e(3,ij,ikx,iky,iz,updatetlevel)
2._dp/3._dp * (2._dp*j_dp*kernel_e(ij,ikx,iky,iz) - (j_dp+1)*kernel_e(ij+1,ikx,iky,iz) - j_dp*kernel_e(ij-1,ikx,iky,iz))&
* (moments_e(1,ij,ikx,iky,iz,updatetlevel)+q_e/tau_e*kernel_e(ij,ikx,iky,iz)*phi(ikx,iky,iz)) &
+ SQRT2/3._dp * kernel_e(ij,ikx,iky,iz) * moments_e(3,ij,ikx,iky,iz,updatetlevel)
ENDDO
! ion temperature
......@@ -128,9 +128,9 @@ SUBROUTINE compute_temperature
DO ij = ijs_i, ije_i
j_dp = REAL(ij-1,dp)
temp_i(ikx,iky,iz) = temp_i(ikx,iky,iz) + &
2._dp/3._dp * (2._dp*j_dp*kernel_i(ij,ikx,iky) - (j_dp+1)*kernel_i(ij+1,ikx,iky) - j_dp*kernel_i(ij-1,ikx,iky))&
* (moments_i(1,ij,ikx,iky,iz,updatetlevel)+q_i/tau_i*kernel_i(ij,ikx,iky)*phi(ikx,iky,iz)) + &
SQRT2/3._dp * kernel_i(ij,ikx,iky) * moments_i(3,ij,ikx,iky,iz,updatetlevel)
2._dp/3._dp * (2._dp*j_dp*kernel_i(ij,ikx,iky,iz) - (j_dp+1)*kernel_i(ij+1,ikx,iky,iz) - j_dp*kernel_i(ij-1,ikx,iky,iz))&
* (moments_i(1,ij,ikx,iky,iz,updatetlevel)+q_i/tau_i*kernel_i(ij,ikx,iky,iz)*phi(ikx,iky,iz)) &
+ SQRT2/3._dp * kernel_i(ij,ikx,iky,iz) * moments_i(3,ij,ikx,iky,iz,updatetlevel)
ENDDO
ENDDO
ENDDO
......
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