From 46d880bc2a65b752d80162da7715ae6bb6a08a43 Mon Sep 17 00:00:00 2001
From: Antoine Cyril David Hoffmann <ahoffman@spcpc606.epfl.ch>
Date: Mon, 5 Oct 2020 09:59:55 +0200
Subject: [PATCH] non-linear term written more explicitely
---
src/compute_Sapj.F90 | 94 +++++++++++++++++++++++++++++++++-----------
1 file changed, 70 insertions(+), 24 deletions(-)
diff --git a/src/compute_Sapj.F90 b/src/compute_Sapj.F90
index 8a382789..a7a40a8c 100644
--- a/src/compute_Sapj.F90
+++ b/src/compute_Sapj.F90
@@ -14,11 +14,11 @@ SUBROUTINE compute_Sapj
! COMPLEX(dp), DIMENSION(ips_e:ipe_e, ijs_e:ije_e, Nkr, Nkz) :: Sepj
! COMPLEX(dp), DIMENSION(ips_i:ipe_i, ijs_i:ije_i, Nkr, Nkz) :: Sipj
INTEGER :: in, is
- REAL(dp):: kr, kz, kernel, be_2, bi_2, factn
+ REAL(dp):: kr, kz, kerneln, be_2, bi_2, factn
REAL(dp):: sigmae2_taue_o2, sigmai2_taui_o2
!!!!!!!!!!!!!!!!!!!! ELECTRON non linear term computation (Sepj)!!!!!!!!!!
- sigmae2_taue_o2 = sigma_e**2 * tau_e/2._dp ! factor of the Kernel argument
+ sigmae2_taue_o2 = sigma_e**2 * tau_e/2._dp ! factor of the kerneln argument
ploope: DO ip = ips_e,ipe_e ! Loop over Hermite moments
jloope: DO ij = ijs_e, ije_e ! Loop over Laguerre moments
Sepj(ip,ij,:,:) = 0._dp
@@ -26,30 +26,53 @@ SUBROUTINE compute_Sapj
nloope: DO in = 1,jmaxe+1 ! Loop over laguerre for the sum
- krloope: DO ikr = 1,Nkr ! Loop over kr
- kzloope: DO ikz = 1,Nkz ! Loop over kz
+ krloope1: DO ikr = 1,Nkr ! Loop over kr
+ kzloope1: DO ikz = 1,Nkz ! Loop over kz
kr = krarray(ikr)
kz = kzarray(ikz)
be_2 = sigmae2_taue_o2 * (kr**2 + kz**2)
- kernel = be_2**(in-1)/factn * EXP(-be_2)
+ kerneln = be_2**(in-1)/factn * EXP(-be_2)
! First convolution term
- F_(ikr,ikz) = (kz - kr) * phi(ikr,ikz) !* kernel
+ F_(ikr,ikz) = imagu*kz* phi(ikr,ikz) * kerneln
! Second convolution term
G_(ikr,ikz) = 0._dp ! initialization of the sum
DO is = 1, MIN( in+ij-1, jmaxe+1 ) ! sum truncation on number of moments
G_(ikr,ikz) = G_(ikr,ikz) + &
dnjs(in,ij,is) * moments_e(ip,is,ikr,ikz,updatetlevel)
ENDDO
- G_(ikr,ikz) = (kz - kr) * G_(ikr,ikz)
- ENDDO kzloope
- ENDDO krloope
+ G_(ikr,ikz) = imagu*kr*G_(ikr,ikz)
- CALL convolve_2D_F2F( F_, G_, CONV ) ! Convolve and go back to Fourier space
- !CALL convolve_2D_F2R( F_, G_, CONV ) ! .. or convolve and keep the results in real space**
+ ENDDO kzloope1
+ ENDDO krloope1
+ CALL convolve_2D_F2F( F_, G_, CONV )
Sepj(ip,ij,:,:) = Sepj(ip,ij,:,:) + CONV ! Add it to Sepj (Real space here)
+ krloope2: DO ikr = 1,Nkr ! Loop over kr
+ kzloope2: DO ikz = 1,Nkz ! Loop over kz
+ kr = krarray(ikr)
+ kz = kzarray(ikz)
+ be_2 = sigmae2_taue_o2 * (kr**2 + kz**2)
+ kerneln = be_2**(in-1)/factn * EXP(-be_2)
+
+ ! First convolution term
+ F_(ikr,ikz) = imagu*kr* phi(ikr,ikz) * kerneln
+ ! Second convolution term
+ G_(ikr,ikz) = 0._dp ! initialization of the sum
+ DO is = 1, MIN( in+ij-1, jmaxe+1 ) ! sum truncation on number of moments
+ G_(ikr,ikz) = G_(ikr,ikz) + &
+ dnjs(in,ij,is) * moments_e(ip,is,ikr,ikz,updatetlevel)
+ ENDDO
+ G_(ikr,ikz) = imagu*kz*G_(ikr,ikz)
+
+ ENDDO kzloope2
+ ENDDO krloope2
+
+
+ CALL convolve_2D_F2F( F_, G_, CONV )
+ Sepj(ip,ij,:,:) = Sepj(ip,ij,:,:) - CONV ! Add it to Sepj (Real space here)
+
IF ( in + 1 .LE. jmaxe+1 ) THEN
factn = real(in,dp) * factn ! factorial(n+1)
ENDIF
@@ -71,29 +94,52 @@ SUBROUTINE compute_Sapj
nloopi: DO in = 1,jmaxi+1 ! Loop over laguerre for the sum
- krloopi: DO ikr = 1,Nkr ! Loop over kr
- kzloopi: DO ikz = 1,Nkz ! Loop over kz
- kr = krarray(ikr)
- kz = kzarray(ikz)
- bi_2 = sigmai2_taui_o2 * (kr**2 + kz**2)
- kernel = bi_2**(in-1)/factn * EXP(-bi_2)
-
- F_(ikr,ikz) = (kz - kr) * phi(ikr,ikz) !* kernel
+ krloopi1: DO ikr = 1,Nkr ! Loop over kr
+ kzloopi1: DO ikz = 1,Nkz ! Loop over kz
+ kr = krarray(ikr)
+ kz = kzarray(ikz)
+ bi_2 = sigmai2_taui_o2 * (kr**2 + kz**2)
+ kerneln = bi_2**(in-1)/factn * EXP(-bi_2)
+ ! First convolution term
+ F_(ikr,ikz) = imagu*kz*phi(ikr,ikz) * kerneln
+ ! Second convolution term
G_(ikr,ikz) = 0._dp ! initialization of the sum
DO is = 1, MIN( in+ij-1, jmaxi+1 )
G_(ikr,ikz) = G_(ikr,ikz) + &
dnjs(in,ij,is) * moments_i(ip,is,ikr,ikz,updatetlevel)
ENDDO
- G_(ikr,ikz) = (kz - kr) * G_(ikr,ikz)
- ENDDO kzloopi
- ENDDO krloopi
+ G_(ikr,ikz) = imagu*kr*G_(ikr,ikz)
- CALL convolve_2D_F2F( F_, G_, CONV ) ! Convolve and back to Fourier
- !CALL convolve_2D_F2R( F_, G_, CONV ) ! or convolve and keep the results in real space**
+ ENDDO kzloopi1
+ ENDDO krloopi1
+ CALL convolve_2D_F2F( F_, G_, CONV ) ! Convolve and back to Fourier
Sipj(ip,ij,:,:) = Sipj(ip,ij,:,:) + CONV ! Add it to Sipj (Real space here)
+ krloopi2: DO ikr = 1,Nkr ! Loop over kr
+ kzloopi2: DO ikz = 1,Nkz ! Loop over kz
+ kr = krarray(ikr)
+ kz = kzarray(ikz)
+ bi_2 = sigmai2_taui_o2 * (kr**2 + kz**2)
+ kerneln = bi_2**(in-1)/factn * EXP(-bi_2)
+
+ ! First convolution term
+ F_(ikr,ikz) = imagu*kr*phi(ikr,ikz) * kerneln
+ ! Second convolution term
+ G_(ikr,ikz) = 0._dp ! initialization of the sum
+ DO is = 1, MIN( in+ij-1, jmaxi+1 )
+ G_(ikr,ikz) = G_(ikr,ikz) + &
+ dnjs(in,ij,is) * moments_i(ip,is,ikr,ikz,updatetlevel)
+ ENDDO
+ G_(ikr,ikz) = imagu*kz*G_(ikr,ikz)
+
+ ENDDO kzloopi2
+ ENDDO krloopi2
+
+ CALL convolve_2D_F2F( F_, G_, CONV ) ! Convolve and back to Fourier
+ Sipj(ip,ij,:,:) = Sipj(ip,ij,:,:) - CONV ! Add it to Sipj (Real space here)
+
IF ( in + 1 .LE. jmaxi+1 ) THEN
factn = real(in,dp) * factn ! factorial(n+1)
ENDIF
--
GitLab