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

COSOlver matrices are asked up to 2/3 kperpmax (since AA on larger k)

parent e7b01d2a
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
src/collision_mod.F90
+ 173
− 162
View file @ 607caa60
...@@ -425,175 +425,186 @@ CONTAINS ...@@ -425,175 +425,186 @@ CONTAINS
CALL allocate_array( CiepjF_kp, 1,(pmaxe+1)*(jmaxe+1), 1,(pmaxe+1)*(jmaxe+1), ikps_C,ikpe_C) CALL allocate_array( CiepjF_kp, 1,(pmaxe+1)*(jmaxe+1), 1,(pmaxe+1)*(jmaxe+1), ikps_C,ikpe_C)
DO ikp = ikps_C,ikpe_C ! Loop over everz kperp values DO ikp = ikps_C,ikpe_C ! Loop over everz kperp values
write(kperp_string,'(f6.4)') kparray(ikp) ! we put zeros if kp>2/3kpmax because thoses frequenvies are filtered through AA
NFLR = MIN(25,MAX(5,CEILING(kparray(ikp)**2))) IF(kparray(ikp) .GT. two_third_kpmax .AND. NON_LIN) THEN
write(NFLR_string,'(i2.1)') NFLR CiepjT_kp(:,:,ikp) = 0._dp
!!!!!!!!!!!! Electron matrices !!!!!!!!!!!! CiepjF_kp(:,:,ikp) = 0._dp
! get the self electron colision matrix CeipjT_kp(:,:,ikp) = 0._dp
IF (CO .GT. 0) THEN CeipjF_kp(:,:,ikp) = 0._dp
WRITE(mat_filename,'(a,a,a,a,a,a)') TRIM(selfmat_file),& Ceepj__kp(:,:,ikp) = 0._dp
'NFLR_',TRIM(ADJUSTL(NFLR_string)),'_kperp_',TRIM(ADJUSTL(kperp_string)),'.h5' Ciipj__kp(:,:,ikp) = 0._dp
ELSE ELSE
mat_filename = selfmat_file write(kperp_string,'(f6.4)') kparray(ikp)
ENDIF NFLR = MIN(25,MAX(5,CEILING(kparray(ikp)**2)))
write(NFLR_string,'(i2.1)') NFLR
CALL openf(mat_filename,fid1, 'r', 'D', mpicomm=comm_p); !!!!!!!!!!!! Electron matrices !!!!!!!!!!!!
CALL getatt(fid1,'/Caapj/Ceepj/','Pmaxe',pdime) ! get the self electron colision matrix
CALL getatt(fid1,'/Caapj/Ceepj/','Jmaxe',jdime) IF (CO .GT. 0) THEN
IF ( ((pdime .LT. pmaxe) .OR. (jdime .LT. jmaxe)) .AND. (my_id .EQ. 0)) WRITE(*,*) '!! P,J Matrix too small !!' WRITE(mat_filename,'(a,a,a,a,a,a)') TRIM(selfmat_file),&
'NFLR_',TRIM(ADJUSTL(NFLR_string)),'_kperp_',TRIM(ADJUSTL(kperp_string)),'.h5'
CALL allocate_array( Ceepj_full, 1,(pdime+1)*(jdime+1), 1,(pdime+1)*(jdime+1)) ELSE
CALL getarr(fid1, '/Caapj/Ceepj', Ceepj_full) ! get array (moli format) mat_filename = selfmat_file
! Fill sub array with only usefull polynmials degree ENDIF
DO ip_e = 0,pmaxe ! Loop over rows ! write(*,*) 'loading ', mat_filename
DO ij_e = 0,jmaxe
irow_sub = (jmaxe +1)*ip_e + ij_e +1 CALL openf(mat_filename,fid1, 'r', 'D', mpicomm=comm_p);
irow_full = (jdime +1)*ip_e + ij_e +1 CALL getatt(fid1,'/Caapj/Ceepj/','Pmaxe',pdime)
DO il_e = 0,pmaxe ! Loop over columns CALL getatt(fid1,'/Caapj/Ceepj/','Jmaxe',jdime)
DO ik_e = 0,jmaxe IF ( ((pdime .LT. pmaxe) .OR. (jdime .LT. jmaxe)) .AND. (my_id .EQ. 0)) WRITE(*,*) '!! P,J Matrix too small !!'
icol_sub = (jmaxe +1)*il_e + ik_e +1
icol_full = (jdime +1)*il_e + ik_e +1 CALL allocate_array( Ceepj_full, 1,(pdime+1)*(jdime+1), 1,(pdime+1)*(jdime+1))
Ceepj__kp (irow_sub,icol_sub,ikp) = Ceepj_full (irow_full,icol_full) CALL getarr(fid1, '/Caapj/Ceepj', Ceepj_full) ! get array (moli format)
ENDDO ! Fill sub array with only usefull polynmials degree
ENDDO DO ip_e = 0,pmaxe ! Loop over rows
ENDDO DO ij_e = 0,jmaxe
ENDDO irow_sub = (jmaxe +1)*ip_e + ij_e +1
irow_full = (jdime +1)*ip_e + ij_e +1
CALL closef(fid1) DO il_e = 0,pmaxe ! Loop over columns
DEALLOCATE(Ceepj_full) DO ik_e = 0,jmaxe
icol_sub = (jmaxe +1)*il_e + ik_e +1
! get the Test and Back field electron ion collision matrix icol_full = (jdime +1)*il_e + ik_e +1
IF (CO .GT. 0) THEN Ceepj__kp (irow_sub,icol_sub,ikp) = Ceepj_full (irow_full,icol_full)
WRITE(mat_filename,'(a,a,a,a,a,a,a,a)') TRIM(eimat_file),& ENDDO
'NFLRe_',TRIM(ADJUSTL(NFLR_string)),'_NFLRi_',TRIM(ADJUSTL(NFLR_string)),& ENDDO
'_kperp_',TRIM(ADJUSTL(kperp_string)),'.h5' ENDDO
ELSE ENDDO
mat_filename = eimat_file
ENDIF
! WRITE(*,*) 'loading : ', mat_filename
CALL openf(mat_filename,fid2, 'r', 'D', mpicomm=comm_p);
CALL getatt(fid2,'/Ceipj/CeipjT/','Pmaxi',pdimi)
CALL getatt(fid2,'/Ceipj/CeipjT/','Jmaxi',jdimi)
IF ( (pdimi .LT. pmaxi) .OR. (jdimi .LT. jmaxi) ) WRITE(*,*) '!! Sugama Matrix too small !!'
CALL allocate_array( CeipjT_full, 1,(pdime+1)*(jdime+1), 1,(pdime+1)*(jdime+1))
CALL allocate_array( CeipjF_full, 1,(pdime+1)*(jdime+1), 1,(pdimi+1)*(jdimi+1))
CALL getarr(fid2, '/Ceipj/CeipjT', CeipjT_full)
CALL getarr(fid2, '/Ceipj/CeipjF', CeipjF_full)
! Fill sub array with only usefull polynmials degree
DO ip_e = 0,pmaxe ! Loop over rows
DO ij_e = 0,jmaxe
irow_sub = (jmaxe +1)*ip_e + ij_e +1
irow_full = (jdime +1)*ip_e + ij_e +1
DO il_e = 0,pmaxe ! Loop over columns
DO ik_e = 0,jmaxe
icol_sub = (jmaxe +1)*il_e + ik_e +1
icol_full = (jdime +1)*il_e + ik_e +1
CeipjT_kp(irow_sub,icol_sub,ikp) = CeipjT_full(irow_full,icol_full)
ENDDO
ENDDO
DO il_i = 0,pmaxi ! Loop over columns
DO ik_i = 0,jmaxi
icol_sub = (Jmaxi +1)*il_i + ik_i +1
icol_full = (jdimi +1)*il_i + ik_i +1
CeipjF_kp(irow_sub,icol_sub,ikp) = CeipjF_full(irow_full,icol_full)
ENDDO
ENDDO
ENDDO
ENDDO
CALL closef(fid2)
DEALLOCATE(CeipjF_full)
DEALLOCATE(CeipjT_full)
!!!!!!!!!!!!!!! Ion matrices !!!!!!!!!!!!!!
! get the self electron colision matrix
IF (CO .GT. 0) THEN
WRITE(mat_filename,'(a,a,a,a,a,a)') TRIM(selfmat_file),&
'NFLR_',TRIM(ADJUSTL(NFLR_string)),'_kperp_',TRIM(ADJUSTL(kperp_string)),'.h5'
ELSE
mat_filename = selfmat_file
ENDIF
! WRITE(*,*) 'loading : ', mat_filename
CALL openf(mat_filename, fid3, 'r', 'D', mpicomm=comm_p);
CALL allocate_array( Ciipj_full, 1,(pdimi+1)*(jdimi+1), 1,(pdimi+1)*(jdimi+1)) CALL closef(fid1)
DEALLOCATE(Ceepj_full)
IF ( (pmaxe .EQ. pmaxi) .AND. (jmaxe .EQ. jmaxi) ) THEN ! if same degrees ion and electron matrices are alike so load Ceepj ! get the Test and Back field electron ion collision matrix
CALL getarr(fid3, '/Caapj/Ceepj', Ciipj_full) ! get array (moli format) IF (CO .GT. 0) THEN
ELSE WRITE(mat_filename,'(a,a,a,a,a,a,a,a)') TRIM(eimat_file),&
CALL getarr(fid3, '/Caapj/Ciipj', Ciipj_full) ! get array (moli format) 'NFLRe_',TRIM(ADJUSTL(NFLR_string)),'_NFLRi_',TRIM(ADJUSTL(NFLR_string)),&
ENDIF '_kperp_',TRIM(ADJUSTL(kperp_string)),'.h5'
ELSE
mat_filename = eimat_file
ENDIF
! WRITE(*,*) 'loading : ', mat_filename
CALL openf(mat_filename,fid2, 'r', 'D', mpicomm=comm_p);
CALL getatt(fid2,'/Ceipj/CeipjT/','Pmaxi',pdimi)
CALL getatt(fid2,'/Ceipj/CeipjT/','Jmaxi',jdimi)
IF ( (pdimi .LT. pmaxi) .OR. (jdimi .LT. jmaxi) ) WRITE(*,*) '!! Sugama Matrix too small !!'
CALL allocate_array( CeipjT_full, 1,(pdime+1)*(jdime+1), 1,(pdime+1)*(jdime+1))
CALL allocate_array( CeipjF_full, 1,(pdime+1)*(jdime+1), 1,(pdimi+1)*(jdimi+1))
CALL getarr(fid2, '/Ceipj/CeipjT', CeipjT_full)
CALL getarr(fid2, '/Ceipj/CeipjF', CeipjF_full)
! Fill sub array with only usefull polynmials degree
DO ip_e = 0,pmaxe ! Loop over rows
DO ij_e = 0,jmaxe
irow_sub = (jmaxe +1)*ip_e + ij_e +1
irow_full = (jdime +1)*ip_e + ij_e +1
DO il_e = 0,pmaxe ! Loop over columns
DO ik_e = 0,jmaxe
icol_sub = (jmaxe +1)*il_e + ik_e +1
icol_full = (jdime +1)*il_e + ik_e +1
CeipjT_kp(irow_sub,icol_sub,ikp) = CeipjT_full(irow_full,icol_full)
ENDDO
ENDDO
DO il_i = 0,pmaxi ! Loop over columns
DO ik_i = 0,jmaxi
icol_sub = (Jmaxi +1)*il_i + ik_i +1
icol_full = (jdimi +1)*il_i + ik_i +1
CeipjF_kp(irow_sub,icol_sub,ikp) = CeipjF_full(irow_full,icol_full)
ENDDO
ENDDO
ENDDO
ENDDO
! Fill sub array with only usefull polynmials degree CALL closef(fid2)
DO ip_i = 0,Pmaxi ! Loop over rows DEALLOCATE(CeipjF_full)
DO ij_i = 0,Jmaxi DEALLOCATE(CeipjT_full)
irow_sub = (Jmaxi +1)*ip_i + ij_i +1
irow_full = (jdimi +1)*ip_i + ij_i +1 !!!!!!!!!!!!!!! Ion matrices !!!!!!!!!!!!!!
DO il_i = 0,Pmaxi ! Loop over columns ! get the self electron colision matrix
DO ik_i = 0,Jmaxi IF (CO .GT. 0) THEN
icol_sub = (Jmaxi +1)*il_i + ik_i +1 WRITE(mat_filename,'(a,a,a,a,a,a)') TRIM(selfmat_file),&
icol_full = (jdimi +1)*il_i + ik_i +1 'NFLR_',TRIM(ADJUSTL(NFLR_string)),'_kperp_',TRIM(ADJUSTL(kperp_string)),'.h5'
Ciipj__kp (irow_sub,icol_sub,ikp) = Ciipj_full (irow_full,icol_full) ELSE
ENDDO mat_filename = selfmat_file
ENDDO ENDIF
ENDDO
ENDDO ! WRITE(*,*) 'loading : ', mat_filename
CALL openf(mat_filename, fid3, 'r', 'D', mpicomm=comm_p);
CALL allocate_array( Ciipj_full, 1,(pdimi+1)*(jdimi+1), 1,(pdimi+1)*(jdimi+1))
IF ( (pmaxe .EQ. pmaxi) .AND. (jmaxe .EQ. jmaxi) ) THEN ! if same degrees ion and electron matrices are alike so load Ceepj
CALL getarr(fid3, '/Caapj/Ceepj', Ciipj_full) ! get array (moli format)
ELSE
CALL getarr(fid3, '/Caapj/Ciipj', Ciipj_full) ! get array (moli format)
ENDIF
! Fill sub array with only usefull polynmials degree
DO ip_i = 0,Pmaxi ! Loop over rows
DO ij_i = 0,Jmaxi
irow_sub = (Jmaxi +1)*ip_i + ij_i +1
irow_full = (jdimi +1)*ip_i + ij_i +1
DO il_i = 0,Pmaxi ! Loop over columns
DO ik_i = 0,Jmaxi
icol_sub = (Jmaxi +1)*il_i + ik_i +1
icol_full = (jdimi +1)*il_i + ik_i +1
Ciipj__kp (irow_sub,icol_sub,ikp) = Ciipj_full (irow_full,icol_full)
ENDDO
ENDDO
ENDDO
ENDDO
CALL closef(fid3) CALL closef(fid3)
DEALLOCATE(Ciipj_full) DEALLOCATE(Ciipj_full)
! get the Test and Back field electron ion collision matrix
IF (CO .GT. 0) THEN
WRITE(mat_filename,'(a,a,a,a,a,a,a,a)') TRIM(iemat_file),&
'NFLRe_',TRIM(ADJUSTL(NFLR_string)),'_NFLRi_',TRIM(ADJUSTL(NFLR_string)),&
'_kperp_',TRIM(ADJUSTL(kperp_string)),'.h5'
ELSE
mat_filename = iemat_file
ENDIF
! write(*,*) 'loading : ', mat_filename
CALL openf(mat_filename,fid4, 'r', 'D', mpicomm=comm_p);
CALL allocate_array( CiepjT_full, 1,(pdimi+1)*(jdimi+1), 1,(pdimi+1)*(jdimi+1))
CALL allocate_array( CiepjF_full, 1,(pdimi+1)*(jdimi+1), 1,(pdime+1)*(jdime+1))
CALL getarr(fid4, '/Ciepj/CiepjT', CiepjT_full)
CALL getarr(fid4, '/Ciepj/CiepjF', CiepjF_full)
! Fill sub array with only usefull polynmials degree
DO ip_i = 0,Pmaxi ! Loop over rows
DO ij_i = 0,Jmaxi
irow_sub = (Jmaxi +1)*ip_i + ij_i +1
irow_full = (jdimi +1)*ip_i + ij_i +1
DO il_i = 0,Pmaxi ! Loop over columns
DO ik_i = 0,Jmaxi
icol_sub = (Jmaxi +1)*il_i + ik_i +1
icol_full = (jdimi +1)*il_i + ik_i +1
CiepjT_kp(irow_sub,icol_sub,ikp) = CiepjT_full(irow_full,icol_full)
ENDDO
ENDDO
DO il_e = 0,pmaxe ! Loop over columns
DO ik_e = 0,jmaxe
icol_sub = (jmaxe +1)*il_e + ik_e +1
icol_full = (jdime +1)*il_e + ik_e +1
CiepjF_kp(irow_sub,icol_sub,ikp) = CiepjF_full(irow_full,icol_full)
ENDDO
ENDDO
ENDDO
ENDDO
! get the Test and Back field electron ion collision matrix CALL closef(fid4)
IF (CO .GT. 0) THEN DEALLOCATE(CiepjF_full)
WRITE(mat_filename,'(a,a,a,a,a,a,a,a)') TRIM(iemat_file),& DEALLOCATE(CiepjT_full)
'NFLRe_',TRIM(ADJUSTL(NFLR_string)),'_NFLRi_',TRIM(ADJUSTL(NFLR_string)),&
'_kperp_',TRIM(ADJUSTL(kperp_string)),'.h5'
ELSE
mat_filename = iemat_file
ENDIF ENDIF
! write(*,*) 'loading : ', mat_filename
CALL openf(mat_filename,fid4, 'r', 'D', mpicomm=comm_p);
CALL allocate_array( CiepjT_full, 1,(pdimi+1)*(jdimi+1), 1,(pdimi+1)*(jdimi+1))
CALL allocate_array( CiepjF_full, 1,(pdimi+1)*(jdimi+1), 1,(pdime+1)*(jdime+1))
CALL getarr(fid4, '/Ciepj/CiepjT', CiepjT_full)
CALL getarr(fid4, '/Ciepj/CiepjF', CiepjF_full)
! Fill sub array with only usefull polynmials degree
DO ip_i = 0,Pmaxi ! Loop over rows
DO ij_i = 0,Jmaxi
irow_sub = (Jmaxi +1)*ip_i + ij_i +1
irow_full = (jdimi +1)*ip_i + ij_i +1
DO il_i = 0,Pmaxi ! Loop over columns
DO ik_i = 0,Jmaxi
icol_sub = (Jmaxi +1)*il_i + ik_i +1
icol_full = (jdimi +1)*il_i + ik_i +1
CiepjT_kp(irow_sub,icol_sub,ikp) = CiepjT_full(irow_full,icol_full)
ENDDO
ENDDO
DO il_e = 0,pmaxe ! Loop over columns
DO ik_e = 0,jmaxe
icol_sub = (jmaxe +1)*il_e + ik_e +1
icol_full = (jdime +1)*il_e + ik_e +1
CiepjF_kp(irow_sub,icol_sub,ikp) = CiepjF_full(irow_full,icol_full)
ENDDO
ENDDO
ENDDO
ENDDO
CALL closef(fid4)
DEALLOCATE(CiepjF_full)
DEALLOCATE(CiepjT_full)
ENDDO ENDDO
IF (CO .GT. 0) THEN ! Interpolation of the kperp matrix values on kr kz grid IF (CO .GT. 0) THEN ! Interpolation of the kperp matrix values on kr kz grid
......
src/grid_mod.F90
+ 2
− 0
View file @ 607caa60
...@@ -24,6 +24,7 @@ MODULE grid ...@@ -24,6 +24,7 @@ MODULE grid
! For Orszag filter ! For Orszag filter
REAL(dp), PUBLIC, PROTECTED :: two_third_krmax REAL(dp), PUBLIC, PROTECTED :: two_third_krmax
REAL(dp), PUBLIC, PROTECTED :: two_third_kzmax REAL(dp), PUBLIC, PROTECTED :: two_third_kzmax
REAL(dp), PUBLIC, PROTECTED :: two_third_kpmax
! 1D Antialiasing arrays (2/3 rule) ! 1D Antialiasing arrays (2/3 rule)
REAL(dp), DIMENSION(:), ALLOCATABLE, PUBLIC :: AA_r REAL(dp), DIMENSION(:), ALLOCATABLE, PUBLIC :: AA_r
...@@ -251,6 +252,7 @@ CONTAINS ...@@ -251,6 +252,7 @@ CONTAINS
kparray(ikp) = REAL(ikp-1,dp) * deltakr kparray(ikp) = REAL(ikp-1,dp) * deltakr
ENDDO ENDDO
write(*,*) rank_kr, ': ikps = ', ikps, 'ikpe = ',ikpe write(*,*) rank_kr, ': ikps = ', ikps, 'ikpe = ',ikpe
two_third_kpmax = SQRT(two_third_krmax**2+two_third_kzmax**2)
END SUBROUTINE END SUBROUTINE
SUBROUTINE grid_readinputs SUBROUTINE grid_readinputs
......
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