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

Corrections of Dougherty operator

parent 690e5640
Branches
Tags
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Showing
with 155 additions and 112 deletions
src/collision_mod.F90
+ 155
− 112
View file @ 5ed2781d
...@@ -22,82 +22,103 @@ CONTAINS ...@@ -22,82 +22,103 @@ CONTAINS
USE grid, ONLY: jmaxe, parray_e, jarray_e, krarray, kzarray USE grid, ONLY: jmaxe, parray_e, jarray_e, krarray, kzarray
USE prec_const USE prec_const
USE time_integration, ONLY: updatetlevel USE time_integration, ONLY: updatetlevel
USE model, ONLY: sigmae2_taue_o2, tau_e, nu_e USE model, ONLY: sigmae2_taue_o2, tau_e, nu_e, q_e
IMPLICIT NONE IMPLICIT NONE
INTEGER, INTENT(IN) :: ip,ij,ikr,ikz INTEGER, INTENT(IN) :: ip,ij,ikr,ikz
COMPLEX(dp), INTENT(INOUT) :: TColl COMPLEX(dp), INTENT(INOUT) :: TColl
COMPLEX(dp) :: n_,upar_,uperp_,Tpar_, Tperp_ COMPLEX(dp) :: n_,upar_,uperp_,Tpar_, Tperp_
COMPLEX(dp) :: Dpj, Ppj, T_, ibe_ COMPLEX(dp) :: Dpj, Ppj, T_, ibe_
COMPLEX(dp) :: nadiab_moment_0j, nadiab_moment_0jp1, nadiab_moment_0jm1 COMPLEX(dp) :: nadiab_moment_0j, nadiab_moment_0jp1, nadiab_moment_0jm1
INTEGER :: in_ INTEGER :: in_
REAL :: n_dp, j_dp, p_dp, be_2 REAL(dp) :: n_dp, j_dp, p_dp, be_2, q_e_tau_e
!** Auxiliary variables ** !** Auxiliary variables **
p_dp = REAL(parray_e(ij),dp) p_dp = REAL(parray_e(ip),dp)
j_dp = REAL(jarray_e(ij),dp) j_dp = REAL(jarray_e(ij),dp)
be_2 = (krarray(ikr)**2 + kzarray(ikz)**2) * sigmae2_taue_o2 be_2 = (krarray(ikr)**2 + kzarray(ikz)**2) * sigmae2_taue_o2
ibe_ = imagu*2._dp*SQRT(be_2) ibe_ = imagu*2._dp*SQRT(be_2)
q_e_tau_e = q_e/tau_e
!** build fluid moments used for collison operator **
n_ = 0._dp
upar_ = 0._dp; uperp_ = 0._dp
Tpar_ = 0._dp; Tperp_ = 0._dp
DO in_ = 1,jmaxe+1
n_dp = REAL(in_-1,dp)
! Nonadiabatic moments (only different from moments when p=0)
nadiab_moment_0j = moments_e(1,in_,ikr,ikz,updatetlevel) + kernel_e(in_,ikr,ikz)*phi(ikr,ikz)/tau_e
nadiab_moment_0jp1 = 0._dp; nadiab_moment_0jm1 = 0._dp
IF (n_dp+1 .LE. jmaxe)& ! Truncation
nadiab_moment_0jp1 = moments_e(1,in_+1,ikr,ikz,updatetlevel) + kernel_e(in_+1,ikr,ikz)*phi(ikr,ikz)/tau_e
IF (n_dp-1 .GE. 0)& ! Truncation
nadiab_moment_0jm1 = moments_e(1,in_-1,ikr,ikz,updatetlevel) + kernel_e(in_-1,ikr,ikz)*phi(ikr,ikz)/tau_e
! Density
n_ = n_ + Kernel_e(in_,ikr,ikz) * nadiab_moment_0j
! Parallel velocity
upar_ = upar_ + Kernel_e(in_,ikr,ikz) * moments_e(2,in_,ikr,ikz,updatetlevel)
! Perpendicular velocity
uperp_ = uperp_ + ibe_*0.5_dp*Kernel_e(in_,ikr,ikz) * (nadiab_moment_0j -nadiab_moment_0jp1)
! Parallel temperature
Tpar_ = Tpar_ + Kernel_e(in_,ikr,ikz) * (SQRT2*moments_e(3,in_,ikr,ikz,updatetlevel) + nadiab_moment_0j)
! Perpendicular temperature
Tperp_ = Tperp_ + Kernel_e(in_,ikr,ikz) * ((2._dp*n_dp+1._dp)*nadiab_moment_0j - n_dp*nadiab_moment_0jm1 - (n_dp+1)*nadiab_moment_0jp1)
ENDDO
T_ = (Tpar_ + 2._dp*Tperp_)/3._dp - n_
!** Assembling collison operator ** !** Assembling collison operator **
! Velocity-space diffusion ! Velocity-space diffusion (similar to Lenhard Bernstein)
TColl = -(2._dp*p_dp + j_dp + be_2)*moments_e(ip,ij,ikr,ikz,updatetlevel) TColl = -(p_dp + 2._dp*j_dp + be_2)*moments_e(ip,ij,ikr,ikz,updatetlevel)
IF( p_dp .EQ. 0 ) & ! Get adiabatic moment
TColl = TColl - (2._dp*p_dp + j_dp + be_2) * Kernel_e(ij,ikr,ikz)*phi(ikr,ikz)/tau_e !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! Non zero term for p = 0 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
IF( p_dp .EQ. 0 ) THEN ! Kronecker p0
! Add energy restoring term ! Get adiabatic moment
IF( p_dp .eq. 0 ) THEN ! kronecker p0 TColl = TColl - (p_dp + 2._dp*j_dp + be_2) * q_e_tau_e * Kernel_e(ij,ikr,ikz)*phi(ikr,ikz)
TColl = TColl + T_* 4._dp*j_dp * Kernel_e(ij,ikr,ikz) !** build required fluid moments **
IF( j_dp+1 .LE. jmaxe ) & ! Truncation n_ = 0._dp
TColl = TColl - T_*2._dp*(j_dp + 1._dp)*Kernel_e(ij+1,ikr,ikz) upar_ = 0._dp; uperp_ = 0._dp
IF( j_dp-1 .GE. 0 ) & ! Truncation Tpar_ = 0._dp; Tperp_ = 0._dp
TColl = TColl - T_*2._dp*j_dp*Kernel_e(ij-1,ikr,ikz) DO in_ = 1,jmaxe+1
ENDIF n_dp = REAL(in_-1,dp)
IF( p_dp .eq. 2 ) & ! kronecker p2 ! Nonadiabatic moments (only different from moments when p=0)
TColl = TColl + T_*SQRT2*Kernel_e(ij,ikr,ikz) nadiab_moment_0j = moments_e(1,in_ ,ikr,ikz,updatetlevel) + q_e_tau_e * kernel_e(in_ ,ikr,ikz)*phi(ikr,ikz)
nadiab_moment_0jp1 = moments_e(1,in_+1,ikr,ikz,updatetlevel) + q_e_tau_e * kernel_e(in_+1,ikr,ikz)*phi(ikr,ikz)
!Add momentum restoring term nadiab_moment_0jm1 = moments_e(1,in_-1,ikr,ikz,updatetlevel) + q_e_tau_e * kernel_e(in_-1,ikr,ikz)*phi(ikr,ikz)
IF( p_dp .eq. 1 ) & ! kronecker p1 ! Density
n_ = n_ + Kernel_e(in_,ikr,ikz) * nadiab_moment_0j
! Perpendicular velocity
uperp_ = uperp_ + ibe_*0.5_dp*Kernel_e(in_,ikr,ikz) * (nadiab_moment_0j - nadiab_moment_0jp1)
! Parallel temperature
Tpar_ = Tpar_ + Kernel_e(in_,ikr,ikz) * (SQRT2*moments_e(3,in_,ikr,ikz,updatetlevel) + nadiab_moment_0j)
! Perpendicular temperature
Tperp_ = Tperp_ + Kernel_e(in_,ikr,ikz) * ((2._dp*n_dp+1._dp)* nadiab_moment_0j &
- n_dp * nadiab_moment_0jm1 &
- (n_dp+1)* nadiab_moment_0jp1)
ENDDO
T_ = (Tpar_ + 2._dp*Tperp_)/3._dp - n_
! Add energy restoring term
TColl = TColl + T_* 4._dp * j_dp * Kernel_e(ij ,ikr,ikz)
TColl = TColl - T_* 2._dp * (j_dp + 1._dp) * Kernel_e(ij+1,ikr,ikz)
TColl = TColl - T_* 2._dp * j_dp * Kernel_e(ij-1,ikr,ikz)
TColl = TColl + uperp_*ibe_*((j_dp + 1._dp)* Kernel_e(ij ,ikr,ikz) &
-j_dp * Kernel_e(ij-1,ikr,ikz))
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! Non zero term for p = 1 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
ELSEIF( p_dp .eq. 1 ) THEN ! kronecker p1
!** build required fluid moments **
upar_ = 0._dp
DO in_ = 1,jmaxe+1
! Parallel velocity
upar_ = upar_ + Kernel_e(in_,ikr,ikz) * moments_e(2,in_,ikr,ikz,updatetlevel)
ENDDO
TColl = TColl + upar_*Kernel_e(ij,ikr,ikz) TColl = TColl + upar_*Kernel_e(ij,ikr,ikz)
IF( p_dp .eq. 0 ) & ! kronecker p0 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
TColl = TColl + uperp_*ibe_*( (j_dp + 1._dp)*Kernel_e(ij,ikr,ikz) - j_dp*Kernel_e(ij-1,ikr,ikz))
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! Non zero term for p = 2 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
ELSEIF( p_dp .eq. 2 ) THEN ! kronecker p2
!** build required fluid moments **
n_ = 0._dp
Tpar_ = 0._dp; Tperp_ = 0._dp
DO in_ = 1,jmaxe+1
n_dp = REAL(in_-1,dp)
! Nonadiabatic moments (only different from moments when p=0)
nadiab_moment_0j = moments_e(1,in_ ,ikr,ikz,updatetlevel) + q_e_tau_e * kernel_e(in_ ,ikr,ikz)*phi(ikr,ikz)/tau_e
nadiab_moment_0jp1 = moments_e(1,in_+1,ikr,ikz,updatetlevel) + q_e_tau_e * kernel_e(in_+1,ikr,ikz)*phi(ikr,ikz)/tau_e
nadiab_moment_0jm1 = moments_e(1,in_-1,ikr,ikz,updatetlevel) + q_e_tau_e * kernel_e(in_-1,ikr,ikz)*phi(ikr,ikz)/tau_e
! Density
n_ = n_ + Kernel_e(in_,ikr,ikz) * nadiab_moment_0j
! Parallel temperature
Tpar_ = Tpar_ + Kernel_e(in_,ikr,ikz) * (SQRT2*moments_e(3,in_,ikr,ikz,updatetlevel) + nadiab_moment_0j)
! Perpendicular temperature
Tperp_ = Tperp_ + Kernel_e(in_,ikr,ikz) * ((2._dp*n_dp+1._dp)*nadiab_moment_0j - n_dp*nadiab_moment_0jm1 - (n_dp+1)*nadiab_moment_0jp1)
ENDDO
T_ = (Tpar_ + 2._dp*Tperp_)/3._dp - n_
TColl = TColl + T_*SQRT2*Kernel_e(ij,ikr,ikz)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
ENDIF
! Multiply by electron-ion collision coefficient ! Multiply by electron-ion collision coefficient
TColl = nu_e * TColl TColl = nu_e * TColl
END SUBROUTINE DoughertyGK_e END SUBROUTINE DoughertyGK_e
!******************************************************************************! !******************************************************************************!
!! Doughtery gyrokinetic collision operator for ions !! Doughtery gyrokinetic collision operator for electrons
SUBROUTINE DoughertyGK_i(ip,ij,ikr,ikz,TColl) SUBROUTINE DoughertyGK_i(ip,ij,ikr,ikz,TColl)
USE fields, ONLY: moments_i, phi USE fields, ONLY: moments_i, phi
USE array, ONLY: Kernel_i USE array, ONLY: Kernel_i
...@@ -105,74 +126,96 @@ CONTAINS ...@@ -105,74 +126,96 @@ CONTAINS
USE grid, ONLY: jmaxi, parray_i, jarray_i, krarray, kzarray USE grid, ONLY: jmaxi, parray_i, jarray_i, krarray, kzarray
USE prec_const USE prec_const
USE time_integration, ONLY: updatetlevel USE time_integration, ONLY: updatetlevel
USE model, ONLY: sigmai2_taui_o2, tau_i, nu_i USE model, ONLY: sigmai2_taui_o2, tau_i, nu_i, q_i
IMPLICIT NONE IMPLICIT NONE
INTEGER, INTENT(IN) :: ip,ij,ikr,ikz INTEGER, INTENT(IN) :: ip,ij,ikr,ikz
COMPLEX(dp), INTENT(INOUT) :: TColl COMPLEX(dp), INTENT(INOUT) :: TColl
COMPLEX(dp) :: n_,upar_,uperp_,Tpar_, Tperp_ COMPLEX(dp) :: n_,upar_,uperp_,Tpar_, Tperp_
COMPLEX(dp) :: Dpj, Ppj, T_, ibi_ COMPLEX(dp) :: Dpj, Ppj, T_, ibi_
COMPLEX(dp) :: nadiab_moment_0j, nadiab_moment_0jp1, nadiab_moment_0jm1 COMPLEX(dp) :: nadiab_moment_0j, nadiab_moment_0jp1, nadiab_moment_0jm1
INTEGER :: in_ INTEGER :: in_
REAL :: n_dp, j_dp, p_dp, bi_2 REAL(dp) :: n_dp, j_dp, p_dp, bi_2, q_i_tau_i
!** Auxiliary variables ** !** Auxiliary variables **
p_dp = REAL(parray_i(ij),dp) p_dp = REAL(parray_i(ip),dp)
j_dp = REAL(jarray_i(ij),dp) j_dp = REAL(jarray_i(ij),dp)
bi_2 = (krarray(ikr)**2 + kzarray(ikz)**2) * sigmai2_taui_o2 bi_2 = (krarray(ikr)**2 + kzarray(ikz)**2) * sigmai2_taui_o2
ibi_ = imagu*2._dp*SQRT(bi_2) ibi_ = imagu*2._dp*SQRT(bi_2)
q_i_tau_i = q_i/tau_i
!** build fluid moments used for collison operator **
n_ = 0._dp
upar_ = 0._dp; uperp_ = 0._dp
Tpar_ = 0._dp; Tperp_ = 0._dp
DO in_ = 1,jmaxi+1
n_dp = REAL(in_-1,dp)
! Adiabatic moments (only different from moments when j=0)
nadiab_moment_0j = moments_i(1,in_,ikr,ikz,updatetlevel) + Kernel_i(in_ ,ikr,ikz)*phi(ikr,ikz)/tau_i
nadiab_moment_0jp1 = 0._dp; nadiab_moment_0jm1 = 0._dp
IF (n_dp+1 .LE. jmaxi)&
nadiab_moment_0jp1 = moments_i(1,in_+1,ikr,ikz,updatetlevel) + Kernel_i(in_+1,ikr,ikz)*phi(ikr,ikz)/tau_i
IF (n_dp-1 .GE. 0)&
nadiab_moment_0jm1 = moments_i(1,in_-1,ikr,ikz,updatetlevel) + Kernel_i(in_-1,ikr,ikz)*phi(ikr,ikz)/tau_i
! Density
n_ = n_ + Kernel_i(in_,ikr,ikz) * nadiab_moment_0j
! Parallel velocity
upar_ = upar_ + Kernel_i(in_,ikr,ikz) * moments_i(2,in_,ikr,ikz,updatetlevel)
! Perpendicular velocity
uperp_ = uperp_ + ibi_*0.5_dp*Kernel_i(in_,ikr,ikz) * (nadiab_moment_0j -nadiab_moment_0jp1)
! Parallel temperature
Tpar_ = Tpar_ + Kernel_i(in_,ikr,ikz) * (SQRT2*moments_i(3,in_,ikr,ikz,updatetlevel) + nadiab_moment_0j)
! Perpendicular temperature
Tperp_ = Tperp_ + Kernel_i(in_,ikr,ikz) * ((2._dp*n_dp+1._dp)*nadiab_moment_0j - n_dp*nadiab_moment_0jm1 - (n_dp+1)*nadiab_moment_0jp1)
ENDDO
T_ = (Tpar_ + 2._dp*Tperp_)/3._dp - n_
!** Assembling collison operator ** !** Assembling collison operator **
! Velocity-space diffusion ! Velocity-space diffusion (similar to Lenhard Bernstein)
TColl = -(2._dp*p_dp + j_dp + bi_2)*moments_i(ip,ij,ikr,ikz,updatetlevel) TColl = -(p_dp + 2._dp*j_dp + bi_2)*moments_i(ip,ij,ikr,ikz,updatetlevel)
IF( p_dp .EQ. 0 ) & ! Get adiabatic moment
TColl = TColl - (2._dp*p_dp + j_dp + bi_2) * Kernel_i(ij,ikr,ikz)*phi(ikr,ikz)/tau_i !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! Non zero term for p = 0 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
IF( p_dp .EQ. 0 ) THEN ! Kronecker p0
! Add energy restoring term ! Get adiabatic moment
IF( p_dp .EQ. 0 ) THEN TColl = TColl - (p_dp + 2._dp*j_dp + bi_2) * q_i_tau_i * Kernel_i(ij,ikr,ikz)*phi(ikr,ikz)
TColl = TColl + T_* 4._dp*j_dp * Kernel_i(ij,ikr,ikz) !** build required fluid moments **
IF( j_dp+1 .LE. jmaxi ) & n_ = 0._dp
TColl = TColl - T_*2._dp*(j_dp + 1._dp)*Kernel_i(ij+1,ikr,ikz) upar_ = 0._dp; uperp_ = 0._dp
IF( j_dp-1 .GE. 0 ) & Tpar_ = 0._dp; Tperp_ = 0._dp
TColl = TColl - T_*2._dp*j_dp*Kernel_i(ij-1,ikr,ikz) DO in_ = 1,jmaxi+1
ENDIF n_dp = REAL(in_-1,dp)
IF( p_dp .eq. 2 ) & ! Nonadiabatic moments (only different from moments when p=0)
nadiab_moment_0j = moments_i(1,in_ ,ikr,ikz,updatetlevel) + q_i_tau_i * kernel_i(in_ ,ikr,ikz)*phi(ikr,ikz)
nadiab_moment_0jp1 = moments_i(1,in_+1,ikr,ikz,updatetlevel) + q_i_tau_i * kernel_i(in_+1,ikr,ikz)*phi(ikr,ikz)
nadiab_moment_0jm1 = moments_i(1,in_-1,ikr,ikz,updatetlevel) + q_i_tau_i * kernel_i(in_-1,ikr,ikz)*phi(ikr,ikz)
! Density
n_ = n_ + Kernel_i(in_,ikr,ikz) * nadiab_moment_0j
! Perpendicular velocity
uperp_ = uperp_ + ibi_*0.5_dp*Kernel_i(in_,ikr,ikz) * (nadiab_moment_0j -nadiab_moment_0jp1)
! Parallel temperature
Tpar_ = Tpar_ + Kernel_i(in_,ikr,ikz) * (SQRT2*moments_i(3,in_,ikr,ikz,updatetlevel) + nadiab_moment_0j)
! Perpendicular temperature
Tperp_ = Tperp_ + Kernel_i(in_,ikr,ikz) * ((2._dp*n_dp+1._dp)* nadiab_moment_0j &
- n_dp * nadiab_moment_0jm1 &
- (n_dp+1)* nadiab_moment_0jp1)
ENDDO
T_ = (Tpar_ + 2._dp*Tperp_)/3._dp - n_
! Add energy restoring term
TColl = TColl + T_* 4._dp * j_dp * Kernel_i(ij ,ikr,ikz)
TColl = TColl - T_* 2._dp * (j_dp + 1._dp) * Kernel_i(ij+1,ikr,ikz)
TColl = TColl - T_* 2._dp * j_dp * Kernel_i(ij-1,ikr,ikz)
TColl = TColl + uperp_*ibi_*((j_dp + 1._dp)* Kernel_i(ij ,ikr,ikz) &
-j_dp * Kernel_i(ij-1,ikr,ikz))
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! Non zero term for p = 1 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
ELSEIF( p_dp .eq. 1 ) THEN ! kronecker p1
!** build required fluid moments **
upar_ = 0._dp
DO in_ = 1,jmaxi+1
! Parallel velocity
upar_ = upar_ + Kernel_i(in_,ikr,ikz) * moments_i(2,in_,ikr,ikz,updatetlevel)
ENDDO
TColl = TColl + upar_*Kernel_i(ij,ikr,ikz)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! Non zero term for p = 2 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
ELSEIF( p_dp .eq. 2 ) THEN ! kronecker p2
!** build required fluid moments **
n_ = 0._dp
Tpar_ = 0._dp; Tperp_ = 0._dp
DO in_ = 1,jmaxi+1
n_dp = REAL(in_-1,dp)
! Nonadiabatic moments (only different from moments when p=0)
nadiab_moment_0j = moments_i(1,in_ ,ikr,ikz,updatetlevel) + q_i_tau_i * kernel_i(in_ ,ikr,ikz)*phi(ikr,ikz)
nadiab_moment_0jp1 = moments_i(1,in_+1,ikr,ikz,updatetlevel) + q_i_tau_i * kernel_i(in_+1,ikr,ikz)*phi(ikr,ikz)
nadiab_moment_0jm1 = moments_i(1,in_-1,ikr,ikz,updatetlevel) + q_i_tau_i * kernel_i(in_-1,ikr,ikz)*phi(ikr,ikz)
! Density
n_ = n_ + Kernel_i(in_,ikr,ikz) * nadiab_moment_0j
! Parallel temperature
Tpar_ = Tpar_ + Kernel_i(in_,ikr,ikz) * (SQRT2*moments_i(3,in_,ikr,ikz,updatetlevel) + nadiab_moment_0j)
! Perpendicular temperature
Tperp_ = Tperp_ + Kernel_i(in_,ikr,ikz) * ((2._dp*n_dp+1._dp)*nadiab_moment_0j - n_dp*nadiab_moment_0jm1 - (n_dp+1)*nadiab_moment_0jp1)
ENDDO
T_ = (Tpar_ + 2._dp*Tperp_)/3._dp - n_
TColl = TColl + T_*SQRT2*Kernel_i(ij,ikr,ikz) TColl = TColl + T_*SQRT2*Kernel_i(ij,ikr,ikz)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! Add momentum restoring term ENDIF
IF( p_dp .eq. 1 ) &
TColl = TColl + upar_*Kernel_i(ij,ikr,ikz)
IF( p_dp .eq. 0 ) &
TColl = TColl + uperp_*ibi_*( (j_dp + 1._dp)*Kernel_i(ij,ikr,ikz) - j_dp*Kernel_i(ij-1,ikr,ikz))
! Multiply by ion-ion collision coefficient ! Multiply by ion-ion collision coefficient
TColl = nu_i * TColl TColl = nu_i * TColl
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment