From 5405686888b1e859dc26302495dbd492a89fe202 Mon Sep 17 00:00:00 2001
From: Antoine Hoffmann <antoine.hoffmann@epfl.ch>
Date: Wed, 22 Mar 2023 13:50:21 +0100
Subject: [PATCH] Optimizations and tests
---
matlab/profiler.m | 25 +++-
src/array_mod.F90 | 4 +-
src/closure_mod.F90 | 121 +++++++++---------
src/cosolver_interface_mod.F90 | 50 +++++---
src/ghosts_mod.F90 | 46 ++++---
src/inital.F90 | 2 -
src/memory.F90 | 2 +
src/moments_eq_rhs_mod.F90 | 22 ++--
src/stepon.F90 | 53 ++++----
src/utility_mod.F90 | 56 ++++----
testcases/cyclone_example/fort.90 | 92 +++++++++++++
testcases/cyclone_example/fort_00.90 | 81 ------------
.../gyacomo_debug | 0
testcases/linear_entropy_mode/fort.90 | 100 ---------------
testcases/linear_entropy_mode/fort_00.90 | 82 ------------
testcases/linear_entropy_mode/gyacomo_1_debug | 1 -
testcases/zpinch_example/fort.90 | 10 +-
testcases/zpinch_example/fort_00.90 | 83 ------------
wk/analysis_gyacomo.m | 2 +-
19 files changed, 307 insertions(+), 525 deletions(-)
delete mode 100644 testcases/cyclone_example/fort_00.90
rename testcases/{linear_entropy_mode => cyclone_example}/gyacomo_debug (100%)
delete mode 100644 testcases/linear_entropy_mode/fort.90
delete mode 100644 testcases/linear_entropy_mode/fort_00.90
delete mode 120000 testcases/linear_entropy_mode/gyacomo_1_debug
delete mode 100644 testcases/zpinch_example/fort_00.90
diff --git a/matlab/profiler.m b/matlab/profiler.m
index 4294f6da..9523e622 100644
--- a/matlab/profiler.m
+++ b/matlab/profiler.m
@@ -44,7 +44,23 @@ Sapj_Ta = mean(diff(Sapj_Tc));
checkfield_Ta = mean(diff(checkfield_Tc));
process_Ta = mean(diff(process_Tc));
diag_Ta = mean(diff(diag_Tc));
-
+miss_Ta = mean(diff(missing_Tc));
+total_Ta = mean(diff(total_Tc));
+names = {...
+ 'Mrhs';
+ 'Advf';
+ 'Ghst';
+ 'Clos';
+ 'Coll';
+ 'Pois';
+ 'Sapj';
+ 'Chck';
+ 'Diag';
+ 'Proc';
+ 'Miss';
+};
+Ts_A = [rhs_Tc(end) adv_field_Tc(end) ghost_Tc(end) clos_Tc(end) coll_Tc(end)...
+ poisson_Tc(end) Sapj_Tc(end) checkfield_Tc(end) diag_Tc(end) process_Tc(end) missing_Tc(end)];
NSTEP_PER_SAMP= mean(diff(Ts0D))/DT_SIM;
%% Plots
@@ -54,8 +70,11 @@ fig = figure;
% colors = rand(N_T,3);
colors = lines(N_T);
p1 = area(Ts0D(2:end),TIME_PER_FCT,'LineStyle','none');
-for i = 1:N_T; p1(i).FaceColor = colors(i,:); end;
-legend('Compute RHS','Adv. fields','ghosts comm', 'closure', 'collision','Poisson','Nonlin','Check+sym', 'Diagnos.', 'Process', 'Missing')
+for i = 1:N_T; p1(i).FaceColor = colors(i,:);
+ LEGEND{i} = sprintf('%s t=%1.1e[s] (%0.1f %s)',names{i},Ts_A(i),Ts_A(i)/total_Tc(end)*100,'\%');
+end;
+legend(LEGEND);
+% legend('Compute RHS','Adv. fields','ghosts comm', 'closure', 'collision','Poisson','Nonlin','Check+sym', 'Diagnos.', 'Process', 'Missing')
xlabel('Sim. Time [$\rho_s/c_s$]'); ylabel('Step Comp. Time [s]')
xlim([Ts0D(2),Ts0D(end)]);
title(sprintf('Proc. 1, total sim. time ~%.0f [h]',CPUTIME/3600))
diff --git a/src/array_mod.F90 b/src/array_mod.F90
index da4d7c9b..f4878080 100644
--- a/src/array_mod.F90
+++ b/src/array_mod.F90
@@ -17,10 +17,12 @@ MODULE array
! Non linear term array (ip,ij,iky,ikx,iz)
COMPLEX(dp), DIMENSION(:,:,:,:,:,:), ALLOCATABLE :: Sapj ! electron
- ! self collision matrix (ia,ip,ij,iky,ikx,iz)
+ ! a-a collision matrix (ia,ip,ij,iky,ikx,iz)
REAL(dp), DIMENSION(:,:,:,:,:,:), ALLOCATABLE :: Caa
! Test and field collision matrices (ia,ib,ip,ij,iky,ikx,iz)
REAL(dp), DIMENSION(:,:,:,:,:,:,:), ALLOCATABLE :: Cab_F, Cab_T
+ ! nu x self collision matrix nuCself = nuaa*Caa + sum_b_neq_a nu_ab Cab_T
+ REAL(dp), DIMENSION(:,:,:,:,:,:), ALLOCATABLE :: nuCself
! Collision term (ip,ij,iky,ikx,iz)
COMPLEX(dp), DIMENSION(:,:,:,:,:,:), ALLOCATABLE :: Capj
diff --git a/src/closure_mod.F90 b/src/closure_mod.F90
index c44a5436..4cff3e6e 100644
--- a/src/closure_mod.F90
+++ b/src/closure_mod.F90
@@ -28,21 +28,22 @@ SUBROUTINE apply_closure_model
! only Napj s.t. p+2j <= 3 are evolved
! -> (p,j) allowed are (0,0),(1,0),(0,1),(2,0),(1,1),(3,0)
! =>> Dmax = min(Pmax,2*Jmax+1)
- z: DO iz = 1,local_nz+ngz
- kx:DO ikx= 1,local_nkx
- ky:DO iky= 1,local_nky
+ ! z: DO iz = 1,local_nz+ngz
+ ! kx:DO ikx= 1,local_nkx
+ ! ky:DO iky= 1,local_nky
j: DO ij = 1,local_nj+ngj
p: DO ip = 1,local_np+ngp
IF ( parray(ip)+2*jarray(ij) .GT. dmax) THEN
- a:DO ia = 1,local_na
- moments(ia,ip,ij,iky,ikx,iz,updatetlevel) = 0._dp
- ENDDO a
+ ! a:DO ia = 1,local_na
+ ! moments(ia,ip,ij,iky,ikx,iz,updatetlevel) = 0._dp
+ moments(ia,ip,ij,:,:,:,updatetlevel) = 0._dp
+ ! ENDDO a
ENDIF
ENDDO p
ENDDO j
- ENDDO ky
- ENDDO kx
- ENDDO z
+ ! ENDDO ky
+ ! ENDDO kx
+ ! ENDDO z
ELSE
ERROR STOP '>> ERROR << Closure scheme not found '
ENDIF
@@ -66,35 +67,37 @@ SUBROUTINE ghosts_upper_truncation
! zero truncation, An+1=0 for n+1>nmax
! applies only for the processes that evolve the highest moment degree
IF(local_jmax .GE. Jmax) THEN
- DO iz = 1,local_nz+ngz
- DO ikx= 1,local_nkx
- DO iky= 1,local_nky
+ ! DO iz = 1,local_nz+ngz
+ ! DO ikx= 1,local_nkx
+ ! DO iky= 1,local_nky
DO ig = 1,ngj/2
- DO ip = 1,local_np+ngp
- DO ia = 1,local_na
- moments(ia,ip,local_nj+ngj/2+ig,iky,ikx,iz,updatetlevel) = 0._dp
- ENDDO
- ENDDO
- ENDDO
- ENDDO
- ENDDO
- ENDDO
+ ! DO ip = 1,local_np+ngp
+ ! DO ia = 1,local_na
+ ! moments(ia,ip,local_nj+ngj/2+ig,iky,ikx,iz,updatetlevel) = 0._dp
+ moments(:,:,local_nj+ngj/2+ig,:,:,:,updatetlevel) = 0._dp
+ ! ENDDO
+ ! ENDDO
+ ENDDO
+ ! ENDDO
+ ! ENDDO
+ ! ENDDO
ENDIF
! applies only for the process that has largest p index
IF(local_pmax .GE. Pmax) THEN
- DO iz = 1,local_nz+ngz
- DO ikx = 1,local_nkx
- DO iky = 1,local_nky
- DO ij = 1,local_nj+ngj
- DO ia = 1,local_na
- DO ig = 1,ngp/2
- moments(ia,local_np+ngp/2+ig,ij,iky,ikx,iz,updatetlevel) = 0._dp
- ENDDO
- ENDDO
- ENDDO
- ENDDO
- ENDDO
- ENDDO
+ ! DO iz = 1,local_nz+ngz
+ ! DO ikx = 1,local_nkx
+ ! DO iky = 1,local_nky
+ ! DO ij = 1,local_nj+ngj
+ DO ig = 1,ngp/2
+ ! DO ia = 1,local_na
+ ! moments(ia,local_np+ngp/2+ig,ij,iky,ikx,iz,updatetlevel) = 0._dp
+ moments(:,local_np+ngp/2+ig,:,:,:,:,updatetlevel) = 0._dp
+ ! ENDDO
+ ENDDO
+ ! ENDDO
+ ! ENDDO
+ ! ENDDO
+ ! ENDDO
ENDIF
END SUBROUTINE ghosts_upper_truncation
@@ -112,36 +115,38 @@ SUBROUTINE ghosts_lower_truncation
INTEGER :: iz,ikx,iky,ip,ia,ij,ig
! zero truncation, An=0 for n<0
IF(local_jmin .EQ. 0) THEN
- DO iz = 1,local_nz+ngz
- DO ikx = 1,local_nkx
- DO iky = 1,local_nky
+ ! DO iz = 1,local_nz+ngz
+ ! DO ikx = 1,local_nkx
+ ! DO iky = 1,local_nky
DO ig = 1,ngj/2
- DO ip = 1,local_np+ngp
- DO ia = 1,local_na
+ ! DO ip = 1,local_np+ngp
+ ! DO ia = 1,local_na
! set to zero the first ghosts cells
- moments(ia,ip,ig,iky,ikx,iz,updatetlevel) = 0._dp
- ENDDO
- ENDDO
- ENDDO
- ENDDO
- ENDDO
- ENDDO
+ ! moments(ia,ip,ig,iky,ikx,iz,updatetlevel) = 0._dp
+ moments(:,:,ig,:,:,:,updatetlevel) = 0._dp
+ ! ENDDO
+ ! ENDDO
+ ENDDO
+ ! ENDDO
+ ! ENDDO
+ ! ENDDO
ENDIF
! applies only for the process that has lowest p index
IF(local_pmin .EQ. 0) THEN
- DO iz = 1,local_nz+ngz
- DO ikx = 1,local_nkx
- DO iky = 1,local_nky
- DO ij = 1,local_nj+ngj
+ ! DO iz = 1,local_nz+ngz
+ ! DO ikx = 1,local_nkx
+ ! DO iky = 1,local_nky
+ ! DO ij = 1,local_nj+ngj
DO ig = 1,ngp/2
- DO ia = 1,local_na
- moments(ia,ig,ij,iky,ikx,iz,updatetlevel) = 0._dp
- ENDDO
- ENDDO
- ENDDO
- ENDDO
- ENDDO
- ENDDO
+ ! DO ia = 1,local_na
+ ! moments(ia,ig,ij,iky,ikx,iz,updatetlevel) = 0._dp
+ moments(:,ig,:,:,:,:,updatetlevel) = 0._dp
+ ENDDO
+ ! ENDDO
+ ! ENDDO
+ ! ENDDO
+ ! ENDDO
+ ! ENDDO
ENDIF
END SUBROUTINE ghosts_lower_truncation
diff --git a/src/cosolver_interface_mod.F90 b/src/cosolver_interface_mod.F90
index 888c722e..019de6f8 100644
--- a/src/cosolver_interface_mod.F90
+++ b/src/cosolver_interface_mod.F90
@@ -19,7 +19,7 @@ CONTAINS
USE MPI
IMPLICIT NONE
LOGICAL, INTENT(IN) :: GK_CO
- COMPLEX(dp), DIMENSION(total_np) :: local_sum, buffer
+ COMPLEX(dp), DIMENSION(total_np) :: local_coll, buffer
COMPLEX(dp), DIMENSION(local_np) :: TColl_distr
COMPLEX(dp) :: Tmp_
INTEGER :: iz,ikx,iky,ij,ip,ia,ikx_C,iky_C,iz_C
@@ -38,18 +38,18 @@ CONTAINS
ENDIF
!! Apply the cosolver collision matrix
CALL apply_cosolver_mat(ia,ip,ij,iky,ikx,iz,ikx_C,iky_C,iz_C,Tmp_)
- local_sum(ip) = Tmp_
+ local_coll(ip) = Tmp_
ENDDO p
IF (num_procs_p .GT. 1) THEN
! Reduce the local_sums to root = 0
- CALL MPI_REDUCE(local_sum, buffer, total_np, MPI_DOUBLE_COMPLEX, MPI_SUM, 0, comm_p, ierr)
+ CALL MPI_REDUCE(local_coll, buffer, total_np, MPI_DOUBLE_COMPLEX, MPI_SUM, 0, comm_p, ierr)
! buffer contains the entire collision term along p, we scatter it between
! the other processes (use of scatterv since Pmax/Np is not an integer)
CALL MPI_SCATTERV(buffer, rcv_p, dsp_p, MPI_DOUBLE_COMPLEX,&
TColl_distr, local_np, MPI_DOUBLE_COMPLEX, &
0, comm_p, ierr)
ELSE
- TColl_distr = local_sum
+ TColl_distr = local_coll
ENDIF
! Write in output variable
DO ip = 1,local_np
@@ -65,17 +65,17 @@ CONTAINS
!******************************************************************************!
!! compute the collision terms in a (Np x Nj x Nkx x Nky) matrix all at once
!******************************************************************************!
- SUBROUTINE apply_cosolver_mat(ia,ip,ij,iky,ikx,iz,ikx_C,iky_C,iz_C,local_sum)
+ SUBROUTINE apply_cosolver_mat(ia,ip,ij,iky,ikx,iz,ikx_C,iky_C,iz_C,local_coll)
USE grid, ONLY: &
- local_na, &
+ total_na, &
local_np, parray,parray_full, ngp,&
total_nj, jarray,jarray_full, dmax, ngj, bar, ngz
- USE array, ONLY: Caa, Cab_T, Cab_F, nadiab_moments
+ USE array, ONLY: nuCself, Cab_F, nadiab_moments
USE model, ONLY: CLOS
USE species, ONLY: nu_ab
IMPLICIT NONE
INTEGER, INTENT(IN) :: ia,ip,ij,iky,ikx,iz,ikx_C,iky_C,iz_C
- COMPLEX(dp), INTENT(OUT) :: local_sum
+ COMPLEX(dp), INTENT(OUT) :: local_coll
INTEGER :: ib,iq,il
INTEGER :: p_int,q_int,j_int,l_int
INTEGER :: izi, iqi, ili
@@ -83,24 +83,21 @@ CONTAINS
p_int = parray_full(ip)
j_int = jarray_full(ij)
!! Apply the cosolver collision matrix
- local_sum = 0._dp ! Initialization
+ local_coll = 0._dp ! Initialization
q:DO iq = 1,local_np
iqi = iq + ngp/2
q_int = parray(iqi)
l:DO il = 1,total_nj
ili = il + ngj/2
l_int = jarray(ili)
- ! self interaction
- local_sum = local_sum + nadiab_moments(ia,iqi,ili,iky,ikx,izi) &
- * nu_ab(ia,ia)*Caa(ia,bar(p_int,j_int), bar(q_int,l_int), iky_C, ikx_C, iz_C)
+ ! self interaction + test interaction
+ local_coll = local_coll + nadiab_moments(ia,iqi,ili,iky,ikx,izi) &
+ * nuCself(ia,bar(p_int,j_int), bar(q_int,l_int), iky_C, ikx_C, iz_C)
! sum the contribution over the other species
- b:DO ib = 1,local_na
- IF((CLOS .NE. 1) .OR. (q_int+2*l_int .LE. dmax)) THEN
- ! Test contribution
- local_sum = local_sum + nadiab_moments(ia,iqi,ili,iky,ikx,izi) &
- * nu_ab(ia,ib) * Cab_T(ia,ib,bar(p_int,j_int), bar(q_int,l_int), iky_C, ikx_C, iz_C)
+ b:DO ib = 1,total_na
+ IF(ib .NE. ia) THEN
! Field contribution
- local_sum = local_sum + nadiab_moments(ib,iqi,ili,iky,ikx,izi) &
+ local_coll = local_coll + nadiab_moments(ib,iqi,ili,iky,ikx,izi) &
* nu_ab(ia,ib) * Cab_F(ia,ib,bar(p_int,j_int), bar(q_int,l_int), iky_C, ikx_C, iz_C)
ENDIF
ENDDO b
@@ -115,14 +112,14 @@ CONTAINS
USE basic, ONLY: allocate_array, speak
USE parallel, ONLY: comm_p, my_id
USE grid, ONLY: &
- local_na,&
+ local_na, total_na, &
local_nkx,local_nky, kparray,&
local_nz, ngz, bar,&
pmax, jmax, ieven
- USE array, ONLY: Caa, Cab_T, Cab_F
+ USE array, ONLY: Caa, Cab_T, Cab_F, nuCself
USE MPI
USE model, ONLY: Na, LINEARITY
- USE species, ONLY: name
+ USE species, ONLY: name, nu_ab
USE futils
IMPLICIT NONE
! Input
@@ -289,6 +286,17 @@ CONTAINS
Cab_F = 0._dp;
Cab_T = 0._dp;
ENDIF
+ ! Build the self matrix
+ ! nuCself = nuaa*Caa + sum_b_neq_a nu_ab Cab_T
+ DO ia = 1,total_na
+ nuCself(ia,:,:,:,:,:) = nu_ab(ia,ia)*Caa(ia,:,:,:,:,:)
+ DO ib = 1,total_na
+ IF(ib .NE. ia) THEN
+ nuCself(ia,:,:,:,:,:) = nuCself(ia,:,:,:,:,:)&
+ +nu_ab(ia,ib)*Cab_T(ia,ib,:,:,:,:,:)
+ ENDIF
+ ENDDO
+ ENDDO
CALL speak('============DONE===========')
END SUBROUTINE load_COSOlver_mat
diff --git a/src/ghosts_mod.F90 b/src/ghosts_mod.F90
index c9861d1c..bd83c119 100644
--- a/src/ghosts_mod.F90
+++ b/src/ghosts_mod.F90
@@ -64,36 +64,44 @@ SUBROUTINE update_ghosts_p_mom
ngp,ngj,ngz
IMPLICIT NONE
COMPLEX(dp), DIMENSION(local_np+ngp) :: slice_p
- INTEGER :: ia,iz,ij,iy,ix
- DO iz = 1+ngz/2,local_nz+ngz/2
- DO ix = 1,local_nkx
- DO iy = 1,local_nky
- DO ij = 1+ngj/2,local_nj+ngj/2
- DO ia = 1,local_na
- slice_p = moments(ia,:,ij,iy,ix,iz,updatetlevel)
- CALL exchange_ghosts_1D(slice_p,nbr_L,nbr_R,local_np,ngp)
- moments(ia,:,ij,iy,ix,iz,updatetlevel) = slice_p
- ENDDO
- ENDDO
- ENDDO
- ENDDO
- ENDDO
- ! INTEGER :: ierr, first, last, ig, count
- ! first = 1 + ngp/2
- ! last = local_np + ngp/2
+ !! SUPER SLOW
+ ! INTEGER :: ia,iz,ij,iy,ix
+ ! DO iz = 1+ngz/2,local_nz+ngz/2
+ ! DO ix = 1,local_nkx
+ ! DO iy = 1,local_nky
+ ! DO ij = 1+ngj/2,local_nj+ngj/2
+ ! DO ia = 1,local_na
+ ! slice_p = moments(ia,:,ij,iy,ix,iz,updatetlevel)
+ ! CALL exchange_ghosts_1D(slice_p,nbr_L,nbr_R,local_np,ngp)
+ ! moments(ia,:,ij,iy,ix,iz,updatetlevel) = slice_p
+ ! ENDDO
+ ! ENDDO
+ ! ENDDO
+ ! ENDDO
+ ! ENDDO
+ INTEGER :: ierr, first, last, ig, count
+ first = 1 + ngp/2
+ last = local_np + ngp/2
! count = local_na*(local_nj+ngj)*local_nky*local_nkx*(local_nz+ngz) ! Number of elements to send
- ! !!!!!!!!!!! Send ghost to right neighbour !!!!!!!!!!!!!!!!!!!!!!
+ !!!!!!!!!!! Send ghost to right neighbour !!!!!!!!!!!!!!!!!!!!!!
! DO ig = 1,ngp/2
! CALL mpi_sendrecv(moments(:,last-(ig-1),:,:,:,:,updatetlevel), count, MPI_DOUBLE_COMPLEX, nbr_R, 14+ig, &
! moments(:,first-ig ,:,:,:,:,updatetlevel), count, MPI_DOUBLE_COMPLEX, nbr_L, 14+ig, &
! comm0, status, ierr)
! ENDDO
- ! !!!!!!!!!!! Send ghost to left neighbour !!!!!!!!!!!!!!!!!!!!!!
+ !!!!!!!!!!! Send ghost to left neighbour !!!!!!!!!!!!!!!!!!!!!!
! DO ig = 1,ngp/2
! CALL mpi_sendrecv(moments(:,first+(ig-1),:,:,:,:,updatetlevel), count, MPI_DOUBLE_COMPLEX, nbr_L, 16+ig, &
! moments(:,last + ig ,:,:,:,:,updatetlevel), count, MPI_DOUBLE_COMPLEX, nbr_R, 16+ig, &
! comm0, status, ierr)
! ENDDO
+ count = (ngp/2)*local_na*(local_nj+ngj)*local_nky*local_nkx*(local_nz+ngz) ! Number of elements to send
+ CALL mpi_sendrecv(moments(:,(last-(ngp/2-1)):(last),:,:,:,:,updatetlevel), count, MPI_DOUBLE_COMPLEX, nbr_R, 14, &
+ moments(:,(first-ngp/2):(first-1) ,:,:,:,:,updatetlevel), count, MPI_DOUBLE_COMPLEX, nbr_L, 14, &
+ comm0, status, ierr)
+ CALL mpi_sendrecv(moments(:,(first):(first+(ngp/2-1)),:,:,:,:,updatetlevel), count, MPI_DOUBLE_COMPLEX, nbr_L, 16, &
+ moments(:,(last+1):(last+ngp/2),:,:,:,:,updatetlevel), count, MPI_DOUBLE_COMPLEX, nbr_R, 16, &
+ comm0, status, ierr)
END SUBROUTINE update_ghosts_p_mom
!Communicate z+1, z+2 moments to left neighboor and z-1, z-2 moments to right one
diff --git a/src/inital.F90 b/src/inital.F90
index 9e9abed6..54262f73 100644
--- a/src/inital.F90
+++ b/src/inital.F90
@@ -101,7 +101,6 @@ SUBROUTINE init_moments
USE initial_par,ONLY: iseed, init_noiselvl, init_background
USE fields, ONLY: moments
USE prec_const, ONLY: dp
- USE utility, ONLY: checkfield
USE model, ONLY: LINEARITY
USE parallel, ONLY: my_id
IMPLICIT NONE
@@ -165,7 +164,6 @@ SUBROUTINE init_gyrodens
ngp, ngj, ngz, iky0, parray, jarray, contains_ky0, AA_x, AA_y
USE fields, ONLY: moments
USE prec_const, ONLY: dp
- USE utility, ONLY: checkfield
USE initial_par,ONLY: iseed, init_noiselvl, init_background
USE model, ONLY: LINEARITY
USE parallel, ONLY: my_id
diff --git a/src/memory.F90 b/src/memory.F90
index 2b27684f..0fccaf0f 100644
--- a/src/memory.F90
+++ b/src/memory.F90
@@ -70,9 +70,11 @@ SUBROUTINE memory
CALL allocate_array( Cab_F, 1,na, 1,na, 1,(pmax+1)*(jmax+1), 1,(pmax+1)*(jmax+1), 1,local_nky, 1,local_nkx, 1,local_nz)
CALL allocate_array( Cab_T, 1,na, 1,na, 1,(pmax+1)*(jmax+1), 1,(pmax+1)*(jmax+1), 1,local_nky, 1,local_nkx, 1,local_nz)
CALL allocate_array( Caa, 1,na, 1,(pmax+1)*(jmax+1), 1,(pmax+1)*(jmax+1), 1,local_nky, 1,local_nkx, 1,local_nz)
+ CALL allocate_array(nuCself, 1,na, 1,(pmax+1)*(jmax+1), 1,(pmax+1)*(jmax+1), 1,local_nky, 1,local_nkx, 1,local_nz)
ELSE !DK collision matrix (same for every k)
CALL allocate_array( Cab_F, 1,na, 1,na, 1,(pmax+1)*(jmax+1), 1,(pmax+1)*(jmax+1), 1,1, 1,1, 1,1)
CALL allocate_array( Cab_T, 1,na, 1,na, 1,(pmax+1)*(jmax+1), 1,(pmax+1)*(jmax+1), 1,1, 1,1, 1,1)
CALL allocate_array( Caa, 1,na, 1,(pmax+1)*(jmax+1), 1,(pmax+1)*(jmax+1), 1,1, 1,1, 1,1)
+ CALL allocate_array(nuCself, 1,na, 1,(pmax+1)*(jmax+1), 1,(pmax+1)*(jmax+1), 1,1, 1,1, 1,1)
ENDIF
END SUBROUTINE memory
\ No newline at end of file
diff --git a/src/moments_eq_rhs_mod.F90 b/src/moments_eq_rhs_mod.F90
index eecef868..65b68f1b 100644
--- a/src/moments_eq_rhs_mod.F90
+++ b/src/moments_eq_rhs_mod.F90
@@ -30,7 +30,7 @@ SUBROUTINE compute_moments_eq_rhs
COMPLEX(dp) :: Mperp, Mpara, Dphi, Dpsi
COMPLEX(dp) :: Unapm1j, Unapm1jp1, Unapm1jm1 ! Terms from mirror force with adiab moments_
COMPLEX(dp) :: i_kx,i_ky
- COMPLEX(dp) :: Napj, RHS
+ COMPLEX(dp) :: Napj, RHS, phikykxz, psikykxz
! Measuring execution time
CALL cpu_time(t0_rhs)
! Spatial loops
@@ -40,8 +40,10 @@ SUBROUTINE compute_moments_eq_rhs
kx = kxarray(ikx) ! radial wavevector
i_kx = imagu * kx ! radial derivative
y:DO iky = 1,local_nky
- ky = kyarray(iky) ! binormal wavevector
- i_ky = imagu * ky ! binormal derivative
+ ky = kyarray(iky) ! binormal wavevector
+ i_ky = imagu * ky ! binormal derivative
+ phikykxz = phi(iky,ikx,izi)
+ psikykxz = psi(iky,ikx,izi)
! Kinetic loops
j:DO ij = 1,local_nj ! This loop is from 1 to jmaxi+1
iji = ij+ngj/2
@@ -55,7 +57,7 @@ SUBROUTINE compute_moments_eq_rhs
a:DO ia = 1,local_na
Napj = moments(ia,ipi,iji,iky,ikx,izi,updatetlevel)
RHS = 0._dp
- IF((CLOS .NE. 1) .OR. (p_int+2*j_int .LE. dmax)) THEN ! for the closure scheme
+ IF((CLOS .NE. 1) .OR. (p_int +2*j_int .LE. dmax)) THEN ! for the closure scheme
!! Compute moments_ mixing terms
! Perpendicular dynamic
! term propto n^{p,j}
@@ -90,18 +92,18 @@ SUBROUTINE compute_moments_eq_rhs
! Parallel magnetic term (Landau damping and the mirror force)
Mpara = gradz_coeff(izi,eo)*(Ldamp + Fmir)
!! Electrical potential term
- IF ( p_int .LE. 2 ) THEN ! kronecker p0 p1 p2
+ IF ( p_int .LE. 2 ) THEN ! kronecker p0 p1 p2
Dphi =i_ky*( xphij (ia,ip,ij)*kernel(ia,iji ,iky,ikx,izi,eo) &
+xphijp1(ia,ip,ij)*kernel(ia,iji+1,iky,ikx,izi,eo) &
- +xphijm1(ia,ip,ij)*kernel(ia,iji-1,iky,ikx,izi,eo) )*phi(iky,ikx,izi)
+ +xphijm1(ia,ip,ij)*kernel(ia,iji-1,iky,ikx,izi,eo) )*phikykxz
ELSE
Dphi = 0._dp
ENDIF
!! Vector potential term
- IF ( (p_int .LE. 3) .AND. (p_int .GE. 1) ) THEN ! Kronecker p1 or p3
+ IF ( (p_int .LE. 3) .AND. (p_int .GE. 1) ) THEN ! Kronecker p1 or p3
Dpsi =-i_ky*( xpsij (ia,ip,ij)*kernel(ia,iji ,iky,ikx,izi,eo) &
+xpsijp1(ia,ip,ij)*kernel(ia,iji+1,iky,ikx,izi,eo) &
- +xpsijm1(ia,ip,ij)*kernel(ia,iji-1,iky,ikx,izi,eo))*psi(iky,ikx,izi)
+ +xpsijm1(ia,ip,ij)*kernel(ia,iji-1,iky,ikx,izi,eo))*psikykxz
ELSE
Dpsi = 0._dp
ENDIF
@@ -129,13 +131,13 @@ SUBROUTINE compute_moments_eq_rhs
!! 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) &
+ IF (p_int .GT. 2) &
RHS = RHS - mu_p*diff_p_coeff*p_int**6*Napj
IF (j_int .GT. 1) &
RHS = RHS - mu_j*diff_j_coeff*j_int**6*Napj
CASE('dvpar4')
! fourth order numerical diffusion in vpar
- IF(p_int .GE. 4) &
+ IF(p_int .GE. 4) &
! Numerical parallel velocity hyperdiffusion "+ dvpar4 g_a" see Pueschel 2010 (eq 33)
! (not used often so not parallelized)
RHS = RHS + mu_p*dv4_Hp_coeff(p_int)*moments(ia,ipi-4,iji,iky,ikx,izi,updatetlevel)
diff --git a/src/stepon.F90 b/src/stepon.F90
index 33f0c0da..1001468b 100644
--- a/src/stepon.F90
+++ b/src/stepon.F90
@@ -63,9 +63,9 @@ SUBROUTINE stepon
END SUBROUTINE assemble_RHS
SUBROUTINE checkfield_all ! Check all the fields for inf or nan
- USE utility,ONLY: checkelem
+ USE utility,ONLY: is_nan, is_inf
USE basic, ONLY: t0_checkfield, t1_checkfield, tc_checkfield
- USE fields, ONLY: phi, moments
+ USE fields, ONLY: phi
USE grid, ONLY: local_na,local_np,local_nj,local_nky,local_nkx,local_nz,&
ngp,ngj,ngz
USE MPI
@@ -74,6 +74,7 @@ SUBROUTINE stepon
IMPLICIT NONE
LOGICAL :: checkf_
INTEGER :: ia, ip, ij, iky, ikx, iz
+ REAL :: sum_
! Execution time start
CALL cpu_time(t0_checkfield)
@@ -82,23 +83,10 @@ SUBROUTINE stepon
mlend=.FALSE.
IF(.NOT.nlend) THEN
- z: DO iz = 1,local_nz+ngz
- kx:DO ikx= 1,local_nkx
- ky:DO iky=1,local_nky
- checkf_ = checkelem(phi(iky,ikx,iz),' phi')
- mlend= (mlend .or. checkf_)
- j: DO ij=1,local_nj+ngj
- p: DO ip=1,local_np+ngp
- a: DO ia=1,local_na
- checkf_ = checkelem(moments(ia,ip,ij,iky,ikx,iz,updatetlevel),' moments')
- mlend = (mlend .or. checkf_)
- ENDDO a
- ENDDO p
- ENDDO j
- ENDDO ky
- ENDDO kx
- ENDDO z
- CALL MPI_ALLREDUCE(mlend, nlend, 1, MPI_LOGICAL, MPI_LOR, MPI_COMM_WORLD, ierr)
+ sum_ = SUM(ABS(phi))
+ checkf_ = (is_nan(sum_,'phi') .OR. is_inf(sum_,'phi'))
+ mlend = (mlend .or. checkf_)
+ CALL MPI_ALLREDUCE(mlend, nlend, 1, MPI_LOGICAL, MPI_LOR, MPI_COMM_WORLD, ierr)
ENDIF
! Execution time end
CALL cpu_time(t1_checkfield)
@@ -107,6 +95,7 @@ SUBROUTINE stepon
SUBROUTINE anti_aliasing
USE fields, ONLY: moments
+ USE time_integration, ONLY: updatetlevel
USE grid, ONLY: local_na,local_np,local_nj,local_nky,local_nkx,local_nz,&
ngp,ngj,ngz, AA_x, AA_y
IMPLICIT NONE
@@ -117,7 +106,8 @@ SUBROUTINE stepon
j: DO ij=1,local_nj+ngj
p: DO ip=1,local_np+ngp
a: DO ia=1,local_na
- moments(ia,ip,ij,iky,ikx,iz,:) = AA_x(ikx)* AA_y(iky) * moments(ia,ip,ij,iky,ikx,iz,:)
+ moments(ia,ip,ij,iky,ikx,iz,updatetlevel) =&
+ AA_x(ikx)* AA_y(iky) * moments(ia,ip,ij,iky,ikx,iz,updatetlevel)
ENDDO a
ENDDO p
ENDDO j
@@ -128,25 +118,28 @@ SUBROUTINE stepon
SUBROUTINE enforce_symmetry ! Force X(k) = X(N-k)* complex conjugate symmetry
USE fields, ONLY: phi, psi, moments
+ USE time_integration, ONLY: updatetlevel
USE grid, ONLY: local_na,local_np,local_nj,total_nkx,local_nz,&
ngp,ngj,ngz, ikx0,iky0, contains_ky0
IMPLICIT NONE
INTEGER :: ia, ip, ij, ikx, iz
IF ( contains_ky0 ) THEN
! moments
- z: DO iz = 1,local_nz+ngz
- j: DO ij=1,local_nj+ngj
- p: DO ip=1,local_np+ngp
- a: DO ia=1,local_na
+ ! z: DO iz = 1,local_nz+ngz
+ ! j: DO ij=1,local_nj+ngj
+ ! p: DO ip=1,local_np+ngp
+ ! a: DO ia=1,local_na
DO ikx=2,total_nkx/2 !symmetry at ky = 0
- moments(ia,ip,ij,iky0,ikx,iz,:) = CONJG(moments(ia,ip,ij,iky0,total_nkx+2-ikx,iz,:))
+ moments(:,:,:,iky0,ikx,:,updatetlevel) = &
+ CONJG(moments(:,:,:,iky0,total_nkx+2-ikx,:,updatetlevel))
END DO
! must be real at origin and top right
- moments(ia,ip,ij, iky0,ikx0,iz,:) = REAL(moments(ia,ip,ij, iky0,ikx0,iz,:),dp)
- ENDDO a
- ENDDO p
- ENDDO j
- ENDDO z
+ moments(:,:,:, iky0,ikx0,:,updatetlevel) = &
+ REAL(moments(:,:,:, iky0,ikx0,:,updatetlevel),dp)
+ ! ENDDO a
+ ! ENDDO p
+ ! ENDDO j
+ ! ENDDO z
! Phi
DO ikx=2,total_nkx/2 !symmetry at ky = 0
phi(iky0,ikx,:) = phi(iky0,total_nkx+2-ikx,:)
diff --git a/src/utility_mod.F90 b/src/utility_mod.F90
index 6433db45..601b6b7f 100644
--- a/src/utility_mod.F90
+++ b/src/utility_mod.F90
@@ -1,6 +1,6 @@
MODULE utility
IMPLICIT NONE
- PUBLIC :: checkfield, checkelem
+ PUBLIC :: is_nan, is_inf!. checkfield, checkelem
CONTAINS
FUNCTION is_nan(x,str) RESULT(isn)
@@ -9,7 +9,7 @@ CONTAINS
USE prec_const, ONLY: dp, stdout
IMPLICIT NONE
- real(dp), INTENT(IN) :: x
+ real, INTENT(IN) :: x
CHARACTER(LEN=*), INTENT(IN) :: str
LOGICAL :: isn
@@ -28,7 +28,7 @@ CONTAINS
USE prec_const, ONLY: dp, stdout
IMPLICIT NONE
- real(dp), INTENT(IN) :: x
+ real, INTENT(IN) :: x
CHARACTER(LEN=*), INTENT(IN) :: str
LOGICAL :: isi
@@ -43,33 +43,33 @@ CONTAINS
END FUNCTION is_inf
- FUNCTION checkfield(n1,n2,n3,field,str) RESULT(mlend)
- use prec_const, ONLY: dp
- IMPLICIT NONE
- !! BUG found (or feature?)
- ! if one put the commented first line (commented) instead of the second one,
- ! no error will be risen by the compiler even if the rank of the array is not matching (should be 3D!)
- ! COMPLEX(dp), DIMENSION(ikys:ikye,ikxs:ikxe), INTENT(IN) :: field
- INTEGER, INTENT(in) :: n1,n2,n3
- COMPLEX(dp), DIMENSION(n1,n2,n3), INTENT(IN) :: field
- CHARACTER(LEN=*), INTENT(IN) :: str
- LOGICAL :: mlend
- COMPLEX(dp) :: sumfield
+ ! FUNCTION checkfield(n1,n2,n3,field,str) RESULT(mlend)
+ ! use prec_const, ONLY: dp
+ ! IMPLICIT NONE
+ ! !! BUG found (or feature?)
+ ! ! if one put the commented first line (commented) instead of the second one,
+ ! ! no error will be risen by the compiler even if the rank of the array is not matching (should be 3D!)
+ ! ! COMPLEX(dp), DIMENSION(ikys:ikye,ikxs:ikxe), INTENT(IN) :: field
+ ! INTEGER, INTENT(in) :: n1,n2,n3
+ ! COMPLEX(dp), DIMENSION(n1,n2,n3), INTENT(IN) :: field
+ ! CHARACTER(LEN=*), INTENT(IN) :: str
+ ! LOGICAL :: mlend
+ ! COMPLEX(dp) :: sumfield
- sumfield=SUM(field)
+ ! sumfield=SUM(field)
- mlend= is_nan( REAL(sumfield),str).OR.is_inf( REAL(sumfield),str) &
- .OR. is_nan(AIMAG(sumfield),str).OR.is_inf(AIMAG(sumfield),str)
- END FUNCTION checkfield
+ ! mlend= is_nan( REAL(sumfield),str).OR.is_inf( REAL(sumfield),str) &
+ ! .OR. is_nan(AIMAG(sumfield),str).OR.is_inf(AIMAG(sumfield),str)
+ ! END FUNCTION checkfield
- FUNCTION checkelem(elem,str) RESULT(mlend)
- use prec_const, ONLY: dp
- IMPLICIT NONE
- COMPLEX(dp), INTENT(IN) :: elem
- CHARACTER(LEN=*), INTENT(IN) :: str
- LOGICAL :: mlend
+ ! FUNCTION checkelem(elem,str) RESULT(mlend)
+ ! use prec_const, ONLY: dp
+ ! IMPLICIT NONE
+ ! COMPLEX(dp), INTENT(IN) :: elem
+ ! CHARACTER(LEN=*), INTENT(IN) :: str
+ ! LOGICAL :: mlend
- mlend= is_nan( REAL(elem),str).OR.is_inf( REAL(elem),str) &
- .OR. is_nan(AIMAG(elem),str).OR.is_inf(AIMAG(elem),str)
- END FUNCTION checkelem
+ ! mlend= is_nan( REAL(elem),str).OR.is_inf( REAL(elem),str) &
+ ! .OR. is_nan(AIMAG(elem),str).OR.is_inf(AIMAG(elem),str)
+ ! END FUNCTION checkelem
END MODULE utility
diff --git a/testcases/cyclone_example/fort.90 b/testcases/cyclone_example/fort.90
index e69de29b..badc2ce7 100644
--- a/testcases/cyclone_example/fort.90
+++ b/testcases/cyclone_example/fort.90
@@ -0,0 +1,92 @@
+&BASIC
+ nrun = 99999999
+ dt = 0.01
+ tmax = 25
+ maxruntime = 356400
+ job2load = -1
+/
+&GRID
+ pmax = 4
+ jmax = 2
+ Nx = 128
+ Lx = 120
+ Ny = 48
+ Ly = 120
+ Nz = 16
+ SG = .f.
+ Nexc = 0
+/
+&GEOMETRY
+ geom = 's-alpha'
+ q0 = 1.4
+ shear = 0.8
+ eps = 0.18
+ kappa = 1.0
+ s_kappa= 0.0
+ delta = 0.0
+ s_delta= 0.0
+ zeta = 0.0
+ s_zeta = 0.0
+ parallel_bc = 'dirichlet'
+ shift_y= 0.0
+/
+&OUTPUT_PAR
+ dtsave_0d = 0.1
+ dtsave_1d = -1
+ dtsave_2d = -1
+ dtsave_3d = 1
+ dtsave_5d = 10
+ write_doubleprecision = .f.
+ write_gamma = .t.
+ write_hf = .t.
+ write_phi = .t.
+ write_Na00 = .t.
+ write_Napj = .t.
+ write_dens = .t.
+ write_fvel = .t.
+ write_temp = .t.
+/
+&MODEL_PAR
+ ! Collisionality
+ CLOS = 0
+ NL_CLOS = 0
+ LINEARITY = 'nonlinear'
+ Na = 1 ! number of species
+ mu_x = 1.0
+ mu_y = 1.0
+ N_HD = 4
+ mu_z = 2.0
+ HYP_V = 'hypcoll'
+ mu_p = 0.0
+ mu_j = 0.0
+ nu = 0.001
+ beta = 0.0
+ ADIAB_E = .t.
+ tau_e = 1.0
+/
+&SPECIES
+ ! ions
+ name_ = 'ions'
+ tau_ = 1.0
+ sigma_= 1.0
+ q_ = 1.0
+ k_N_ = 2.22
+ k_T_ = 6.96
+/
+
+&COLLISION_PAR
+ collision_model = 'DG' !DG/SG/PA/LD (dougherty, sugama, pitch angle, landau)
+ GK_CO = .f.
+ INTERSPECIES = .true.
+ mat_file = 'gk_sugama_P_20_J_10_N_150_kpm_8.0.h5'
+/
+&INITIAL_CON
+ INIT_OPT = 'blob'
+ ACT_ON_MODES = 'donothing'
+ init_background = 0.0
+ init_noiselvl = 0.005
+ iseed = 42
+/
+&TIME_INTEGRATION_PAR
+ numerical_scheme = 'RK4'
+/
diff --git a/testcases/cyclone_example/fort_00.90 b/testcases/cyclone_example/fort_00.90
deleted file mode 100644
index 6ca55404..00000000
--- a/testcases/cyclone_example/fort_00.90
+++ /dev/null
@@ -1,81 +0,0 @@
-&BASIC
- nrun = 100000000
- dt = 0.01
- tmax = 50
- maxruntime = 356400
-/
-&GRID
- pmaxe = 4
- jmaxe = 2
- pmaxi = 4
- jmaxi = 2
- Nx = 128
- Lx = 120
- Ny = 64
- Ly = 160
- Nz = 16
- Nexc = 0
- SG = .true.
-/
-&GEOMETRY
- geom = 's-alpha'
- q0 = 1.4
- shear = 0.8
- eps = 0.18
- parallel_BC = 'dirichlet'
-/
-&OUTPUT_PAR
- nsave_0d = 50
- nsave_1d = -1
- nsave_2d = -1
- nsave_3d = 100
- nsave_5d = 500
- write_doubleprecision = .false.
- write_gamma = .true.
- write_hf = .true.
- write_phi = .true.
- write_Na00 = .false.
- write_Napj = .true.
- write_Sapj = .false.
- write_dens = .true.
- write_temp = .true.
- job2load = -1
-/
-&MODEL_PAR
- ! Collisionality
- CLOS = 0
- NL_CLOS = 0
- LINEARITY = 'nonlinear'
- KIN_E = .false.
- mu_x = 0.1
- mu_y = 0.1
- N_HD = 2
- mu_z = 2.0
- mu_p = 0
- mu_j = 0
- nu = 0.05
- tau_e = 1
- tau_i = 1
- sigma_e = 0.023338
- sigma_i = 1
- q_e = -1
- q_i = 1
- K_Ni = 2.22
- K_Ti = 6.92
-/
-&COLLISION_PAR
- collision_model = 'DG'
- GK_CO = .false.
- INTERSPECIES = .true.
- mat_file = 'null'
-/
-&INITIAL_CON
- INIT_OPT = 'blob'
- ACT_ON_MODES = 'donothing'
- init_background = 0
- init_noiselvl = 1e-3
- iseed = 42
-/
-&TIME_INTEGRATION_PAR
- numerical_scheme = 'RK4'
-/
diff --git a/testcases/linear_entropy_mode/gyacomo_debug b/testcases/cyclone_example/gyacomo_debug
similarity index 100%
rename from testcases/linear_entropy_mode/gyacomo_debug
rename to testcases/cyclone_example/gyacomo_debug
diff --git a/testcases/linear_entropy_mode/fort.90 b/testcases/linear_entropy_mode/fort.90
deleted file mode 100644
index 173b1626..00000000
--- a/testcases/linear_entropy_mode/fort.90
+++ /dev/null
@@ -1,100 +0,0 @@
-&BASIC
- nrun = 1
- dt = 0.01
- tmax = 5
- maxruntime = 356400
- job2load = -1
-/
-&GRID
- pmax = 4
- jmax = 2
- Nx = 2
- Lx = 200
- Ny = 2
- Ly = 60
- Nz = 1
- SG = .f.
- Nexc = 1
-/
-&GEOMETRY
- geom = 'zpinch'
- q0 = 1.4
- shear = 0.0
- eps = 0.18
- kappa = 1.0
- s_kappa= 0.0
- delta = 0.0
- s_delta= 0.0
- zeta = 0.0
- s_zeta = 0.0
- parallel_bc = 'dirichlet'
- shift_y= 0.0
-/
-&OUTPUT_PAR
- dtsave_0d = 0.01
- dtsave_1d = -1
- dtsave_2d = -1
- dtsave_3d = 0.01
- dtsave_5d = 0.01
- write_doubleprecision = .f.
- write_gamma = .t.
- write_hf = .t.
- write_phi = .t.
- write_Na00 = .t.
- write_Napj = .t.
- write_dens = .t.
- write_fvel = .t.
- write_temp = .t.
-/
-&MODEL_PAR
- ! Collisionality
- CLOS = 0
- NL_CLOS = -1
- LINEARITY = 'linear'
- Na = 2 ! number of species
- mu_x = 0
- mu_y = 0
- N_HD = 4
- mu_z = 0
- mu_p = 0
- mu_j = 0
- nu = 0
- beta = 0
- ADIAB_E = .f.
- tau_e = 1.0
-/
-&SPECIES
- ! ions
- name_ = 'ions'
- tau_ = 1.0
- sigma_= 1.0
- q_ = 1.0
- k_N_ = 1.0!2.22
- k_T_ = 0.0!6.96
-/
-&SPECIES
- ! electrons
- name_ = 'electrons'
- tau_ = 1.0
- sigma_= 0.023338
- q_ = -1.0
- k_N_ = 0.0!2.22
- k_T_ = 0.0!6.96
-/
-
-&COLLISION_PAR
- collision_model = 'DG' !DG/SG/PA/LD (dougherty, sugama, pitch angle, landau)
- GK_CO = .f.
- INTERSPECIES = .true.
- !mat_file = 'gk_sugama_P_20_J_10_N_150_kpm_8.0.h5'
-/
-&INITIAL_CON
- INIT_OPT = 'mom00'
- ACT_ON_MODES = 'donothing'
- init_background = 1.0
- init_noiselvl = 0.0
- iseed = 42
-/
-&TIME_INTEGRATION_PAR
- numerical_scheme = 'RK4'
-/
diff --git a/testcases/linear_entropy_mode/fort_00.90 b/testcases/linear_entropy_mode/fort_00.90
deleted file mode 100644
index 158f88d3..00000000
--- a/testcases/linear_entropy_mode/fort_00.90
+++ /dev/null
@@ -1,82 +0,0 @@
-&BASIC
- nrun = 1
- dt = 0.01
- tmax = 5
- maxruntime = 356400
-/
-&GRID
- pmaxe = 4
- jmaxe = 2
- pmaxi = 4
- jmaxi = 2
- Nx = 2
- Lx = 200
- Ny = 2
- Ly = 60
- Nz = 1
- SG = .f.
-/
-&GEOMETRY
- geom = 'zpinch'
- q0 = 1.4
- shear = 0.0
- eps = 0.18
- parallel_bc = 'dirichlet'
-/
-&OUTPUT_PAR
- nsave_0d = 1
- nsave_1d = -1
- nsave_2d = -1
- nsave_3d = 1
- nsave_5d = 1
- write_doubleprecision = .f.
- write_gamma = .t.
- write_hf = .t.
- write_phi = .t.
- write_Na00 = .t.
- write_Napj = .t.
- write_dens = .t.
- write_temp = .t.
- job2load = -1
-/
-&MODEL_PAR
- ! Collisionality
- CLOS = 0
- NL_CLOS = -1
- LINEARITY = 'linear'
- KIN_E = .t.
- mu_x = 0.0
- mu_y = 0.0
- N_HD = 4
- mu_z = 0!0.1
- mu_p = 0
- mu_j = 0
- nu = 0!1.0
- tau_e = 1
- tau_i = 1
- sigma_e = 0.023338
- sigma_i = 1
- q_e = -1
- q_i = 1
- K_Ne = 0!6.96
- K_Te = 0!2.22
- K_Ni = 1!6.96
- K_Ti = 0!2.22
- beta = 0
-/
-&COLLISION_PAR
- collision_model = 'DG' !DG/SG/PA/LD (dougherty, sugama, pitch angle, landau)
- gyrokin_CO = .f.
- interspecies = .true.
- !mat_file = 'gk_sugama_P_20_J_10_N_150_kpm_8.0.h5'
-/
-&INITIAL_CON
- INIT_OPT = 'mom00'
- ACT_ON_MODES = 'donothing'
- init_background = 1.0
- init_noiselvl = 0.0
- iseed = 42
-/
-&TIME_INTEGRATION_PAR
- numerical_scheme = 'RK4'
-/
diff --git a/testcases/linear_entropy_mode/gyacomo_1_debug b/testcases/linear_entropy_mode/gyacomo_1_debug
deleted file mode 120000
index f0c8bbfd..00000000
--- a/testcases/linear_entropy_mode/gyacomo_1_debug
+++ /dev/null
@@ -1 +0,0 @@
-../../../gyacomo_1/bin/gyacomo_debug
\ No newline at end of file
diff --git a/testcases/zpinch_example/fort.90 b/testcases/zpinch_example/fort.90
index 5c5bd290..ee9320ef 100644
--- a/testcases/zpinch_example/fort.90
+++ b/testcases/zpinch_example/fort.90
@@ -1,7 +1,7 @@
&BASIC
nrun = 99999999
dt = 0.01
- tmax = 50
+ tmax = 10
maxruntime = 356400
job2load = -1
/
@@ -60,7 +60,7 @@
mu_p = 0.0
mu_j = 0.0
nu = 0.1
- beta = 0.0
+ beta = 0.01
ADIAB_E = .f.
tau_e = 1.0
/
@@ -84,7 +84,7 @@
/
&COLLISION_PAR
- collision_model = 'DG' !DG/SG/PA/LD (dougherty, sugama, pitch angle, landau)
+ collision_model = 'SG' !DG/SG/PA/LD (dougherty, sugama, pitch angle, landau)
GK_CO = .t.
INTERSPECIES = .true.
mat_file = 'gk_sugama_P_20_J_10_N_150_kpm_8.0.h5'
@@ -92,8 +92,8 @@
&INITIAL_CON
INIT_OPT = 'phi'
ACT_ON_MODES = 'donothing'
- init_background = 0.0
- init_noiselvl = 0.005
+ init_background = 0.005
+ init_noiselvl = 0.00
iseed = 42
/
&TIME_INTEGRATION_PAR
diff --git a/testcases/zpinch_example/fort_00.90 b/testcases/zpinch_example/fort_00.90
deleted file mode 100644
index 5864634f..00000000
--- a/testcases/zpinch_example/fort_00.90
+++ /dev/null
@@ -1,83 +0,0 @@
-&BASIC
- nrun = 100000000
- dt = 0.01
- tmax = 50
- maxruntime = 356400
-/
-&GRID
- pmaxe = 4
- jmaxe = 2
- pmaxi = 4
- jmaxi = 2
- Nx = 128
- Lx = 200
- Ny = 48
- Ly = 60
- Nz = 1
- SG = .f.
-/
-&GEOMETRY
- geom = 'Z-pinch'
- q0 = 0
- shear = 0
- eps = 0
- parallel_bc = 'shearless'
-/
-&OUTPUT_PAR
- nsave_0d = 10
- nsave_1d = -1
- nsave_2d = -1
- nsave_3d = 100
- nsave_5d = 1000
- write_doubleprecision = .f.
- write_gamma = .t.
- write_hf = .t.
- write_phi = .t.
- write_Na00 = .f.
- write_Napj = .t.
- write_Sapj = .f.
- write_dens = .t.
- write_temp = .t.
- job2load = -1
-/
-&MODEL_PAR
- ! Collisionality
- CLOS = 0
- NL_CLOS = -1
- LINEARITY = 'nonlinear'
- KIN_E = .t.
- mu_x = 1.0
- mu_y = 1.0
- N_HD = 4
- mu_z = 0.0
- mu_p = 0
- mu_j = 0
- nu = 0.1
- tau_e = 1
- tau_i = 1
- sigma_e = 0.023338
- sigma_i = 1
- q_e = -1
- q_i = 1
- K_Ne = 2.0
- K_Te = 0.4
- K_Ni = 2.0
- K_Ti = 0.4
- beta = 0
-/
-&COLLISION_PAR
- collision_model = 'DG' !DG/SG/PA/LD (dougherty, sugama, pitch angle, landau)
- GK_CO = .true.
- INTERSPECIES = .true.
- !mat_file = 'gk_sugama_P_20_J_10_N_150_kpm_8.0.h5'
-/
-&INITIAL_CON
- INIT_OPT = 'phi'
- ACT_ON_MODES = 'donothing'
- init_background = 0
- init_noiselvl = 0.005
- iseed = 42
-/
-&TIME_INTEGRATION_PAR
- numerical_scheme = 'RK4'
-/
diff --git a/wk/analysis_gyacomo.m b/wk/analysis_gyacomo.m
index 1ff98e2e..f0e73003 100644
--- a/wk/analysis_gyacomo.m
+++ b/wk/analysis_gyacomo.m
@@ -93,7 +93,7 @@ options.RESOLUTION = 256;
create_film(data,options,'.gif')
end
-if 1
+if 0
%% fields snapshots
% Options
options.INTERP = 0;
--
GitLab