Homec4science
Diffusion COSOlver (improved_sugama) Browse

cosolver/src/load_prec_coeffs_mod.f90
cd1286e2d410improved_sugama

Last Change
Raw File

load_prec_coeffs_mod.f90
View Options

module load_prec_coeffs
! loads the pre-computed coefficients Tljpmf, curlyNpjmn, Pipjm
! Note: - used for GK Coulomb collision operator
use prec_const
use basic
use model
use array
implicit none
PUBLIC
INTEGER, SAVE :: fid_
INTEGER, SAVE :: IDXMAX ! ... maximum indices for Tljpmf
! Tljpmf maximum indices
INTEGER, SAVE :: Tidx_l_max, Tidx_j_max, Tidx_p_max, Tidx_f_max
! CurlyNpjmn maximum indices
INTEGER, SAVE :: Nidx_p_max, Nidx_j_max, Nidx_n_max
CHARACTER(len=*),PARAMETER :: Tljpmf_in = 'log.in'
CHARACTER(len=*),PARAMETER :: CurlyNpjmn_in = 'log.in'
CHARACTER(len=*),PARAMETER :: Pipjm_in = 'Pipjm.in'
CHARACTER(len=*), PARAMETER :: Tljpmf_dir = 'Tljpmf'
CHARACTER(len=*), PARAMETER :: CurlyNpjmn_dir = 'CurlyNpjmn'
contains
!
!___________________________________________________________
!
subroutine load_gk_coulomb_coeffs(dirCab)
! load the Tljpmf and curlyNpjmn coeffs precomputed for gkcoulomb
implicit none
! directory where the coefficients are stored
CHARACTER(len=*), intent(in) :: dirCab ! species depend directory
IF( me .eq. 0) WRITE(*,*) 'Load T_lj^pmf coefficients in ', dirCab, ' ...'
CALL read_Tljpmf(dirCab)
IF( me .eq. 0) WRITE(*,*) '... T_lj^pmf loaded'
IF(me .eq. 0) WRITE(*,*) 'Load curlyNpjm ...'
call load_CurlyNpjmn
IF(me .eq. 0) WRITE(*,*) '... curlyNpjm loaded'
end subroutine load_gk_coulomb_coeffs
!
!___________________________________________________________
!
subroutine read_Tljpmf(dirCab)
!
! read NN Tljpmf coefficients from Tljpmf directory dir
!
implicit none
!
! directory where the coefficients are stored
CHARACTER(len=*), intent(in) :: dirCab ! species depend directory
! locs. vars.
integer :: NN,ENDOF,iNN
CHARACTER(len=256) :: buff_,filename_
INTEGER :: IT1,IT2,IT3,IT4,IT5
INTEGER :: FMAXX
CHARACTER(len=100),DIMENSION(:),ALLOCATABLE :: TableNames
INTEGER :: Ntab,iTab
CHARACTER(len=T5len) :: coeffstr
!
!
! 1. Get the number of NN coefficients
!
FMAXX = max(Jmaxe,Jmaxi) + max(neFLR,niFLR)
IDXMAX = max(FMAXX,JEmaxx,PMMAXX,PMAXI,PMAXE)+1
!
WRITE(filename_, *) trim(adjustl(cosolver_dir))//'/'// Tljpmf_dir // '/' // dirCab // '/'// Tljpmf_in
!
IF( me .eq. 0) WRITE(*, *) 'Read log.in file ', trim(adjustl(filename_))
!
OPEN(fid_,FILE=trim(adjustl(filename_)))
! Count the number of files
NTab = 0
ENDOF = 0
DO WHILE(ENDOF .eq. 0)
READ(fid_,*,IOSTAT =ENDOF) buff_
NTab = NTab +1
ENDDO
!
NTab = NTab -1
!
CLOSE(fid_)
!
! Allocate array of Tljpmf table
IF( allocated(TableNames)) deallocate(TableNames)
ALLOCATE(TableNames(1:NTab))
!
!
! Get the name of the table containing T5 coeff.
!
OPEN(fid_ ,FILE=trim(adjustl(filename_)))
ENDOF = 0
iTab = 1
DO WHILE (ENDOF .eq. 0 )
READ(fid_,*,IOSTAT=ENDOF) buff_
IF (iTab .le. NTab) THEN
TableNames(iTab) = adjustl(trim(buff_))
iTab = iTab +1
ENDIF
ENDDO
CLOSE(fid_)
!
!
Tidx_l_max = 0
Tidx_j_max = 0
Tidx_p_max = 0
Tidx_f_max = 0
!
NN = 0
DO iTab=1,NTab
! name of table to read
write(filename_, *) trim(adjustl(cosolver_dir))//'/'// Tljpmf_dir // '/' // dirCab // '/'// TableNames(iTab)
OPEN(fid_,FILE = trim(adjustl(filename_)))
ENDOF = 0
DO WHILE (ENDOF .eq. 0 )
READ(fid_,*,IOSTAT=ENDOF) IT1,IT2,IT3,IT4,IT5
IF (IT1 .le. IDXMAX .and. IT2 .le. IDXMAX .and. IT3 .le. IDXMAX .and. IT4 .le. IDXMAX .and. IT5 .le. IDXMAX) THEN
IF( IT1 .gt. Tidx_l_max ) Tidx_l_max = IT1
IF( IT2 .gt. Tidx_j_max ) Tidx_j_max = IT2
IF( IT3 .gt. Tidx_p_max ) Tidx_p_max = IT3
IF( IT5 .gt. Tidx_f_max ) Tidx_f_max = IT5
NN = NN+1
ENDIF
ENDDO
CLOSE(fid_)
ENDDO
!
IF( me .eq. 0) WRITE(*,*) '(l,j,p,m,f)< (', Tidx_l_max, ',', Tidx_j_max, ',', Tidx_p_max, ',:,', Tidx_f_max, ')'
!
! 2. Allocate arrays
!
IF( allocated(TljpmfSTR)) deallocate(TljpmfSTR)
allocate(TljpmfSTR(1:NN))
IF( allocated(NNLJPMF)) deallocate(NNLJPMF)
allocate(NNLJPMF(1:NN))
!
! 3. Inquire NN coefficients
!
iNN = 1
!
DO iTab =1, NTab
write(filename_, *) trim(adjustl(cosolver_dir))//'/'// Tljpmf_dir // '/' // dirCab // '/'// TableNames(iTab)
OPEN(FID_,FILE = filename_)
ENDOF = 0
DO WHILE (ENDOF .eq. 0 )
READ(fid_,*,IOSTAT=ENDOF) IT1,IT2,IT3,IT4,IT5,coeffstr
IF (IT1 .le. IDXMAX .and. IT2 .le. IDXMAX .and. IT3 .le. IDXMAX .and. IT4 .le. IDXMAX .and. IT5 .le. IDXMAX ) THEN
!
CALL GET_NNLJPMF(IT1,IT2,IT3,IT4,IT5,NNLJPMF(iNN))
TljpmfSTR(iNN) = coeffstr
!! write(*,*) NNLJPMF(iNN),TljpmfSTR(iNN)
IF(ENDOF .eq. 0) iNN = iNN+1
ENDIF
ENDDO
!
CLOSE(fid_)
ENDDO
!
!
!
! Free memory
DEALLOCATE(TableNames)
!
end subroutine read_Tljpmf
!
!______________________________________________________________________________________________
!
SUBROUTINE GET_NNLJPMF(l,j,p,m,f,NNLJPMF_)
!
! construct character string l//$//j//$//p//$//m//$//f
!
IMPLICIT NONE
!
INTEGER,intent(in) :: l,j,p,m,f
CHARACTER(len=NNLJPMFlen),intent(inout) :: NNLJPMF_
!
CHARACTER*2 :: SIT1,SIT2,SIT3,SIT4,SIT5
!
! Convert int to string
WRITE(SIT1,'(i2)') l
WRITE(SIT2,'(i2)') j
WRITE(SIT3,'(i2)') p
WRITE(SIT4,'(i2)') m
WRITE(SIT5,'(i2)') f
!
NNLJPMF_ = trim(adjustl(SIT1))//'$'//trim(adjustl(SIT2))//'$'//trim(adjustl(SIT3))//'$'//trim(adjustl(SIT4))//'$'//trim(adjustl(SIT5))
!
END SUBROUTINE GET_NNLJPMF
!______________________________________________________________________________________________
!
SUBROUTINE GET_NNPJM(p,j,m,NNPJM_)
!
! construct character string p//$//j//$//m
!
IMPLICIT NONE
!
INTEGER,intent(in) :: p,j,m
CHARACTER(len=NNPJMlen),intent(inout) :: NNPJM_
!
CHARACTER*2 :: SIT1,SIT2,SIT3
!
! Convert int to string
WRITE(SIT1,'(i2)') p
WRITE(SIT2,'(i2)') j
WRITE(SIT3,'(i2)') m
!
NNPJM_ = trim(adjustl(SIT1))//'$'//trim(adjustl(SIT2))//'$'//trim(adjustl(SIT3))
!
END SUBROUTINE GET_NNPJM
!______________________________________________________________________________________________
!
FUNCTION GET_Tljpmf(l,j,p,m,f) RESULT(Tljpmfstr_)
!
! Returns the Tljpmf string from character string array Tljpmfstr
!
IMPLICIT NONE
!
INTEGER,intent(in) :: l,j,p,m,f
CHARACTER(len=T5len) :: Tljpmfstr_
! Local variables
CHARACTER(len=NNLJPMFlen) :: NNLJPMF_
INTEGER :: POS,i,inquireMaxInd
!
!
IF( l .gt. Tidx_l_max .or. j .gt. Tidx_j_max .or. p .gt. Tidx_p_max .or. f .gt. Tidx_f_max) THEN
WRITE(*,* ) 'ERROR Tljpmf: inquire unexisting indices', l,'/' , Tidx_l_max, ',',j,',', Tidx_j_max, ', ', p , '/', Tidx_p_max, ',',f ,'/' ,Tidx_f_max
stop
ENDIF
! IF( l .gt. IDXMAX .or. j .gt. IDXMAX .or. p .gt. IDXMAX .or. m .gt. IDXMAX .or. f .gt. IDXMAX) THEN
! inquireMaxInd = max(l,j,p,m,f)
! write(*,*) 'WARNING !!! NEED TO LOAD MORE T(l,j,p,m,f) .... ! Inquire ', inquireMaxInd
! stop
! ENDIF
POS =0
! Get string id
CALL GET_NNLJPMF(l,j,p,m,f,NNLJPMF_)
! Get Tljpmf position
DO i=1,size(NNLJPMF)
IF( NNLJPMF(i) == NNLJPMF_) THEN
POS = i
GOTO 55
ENDIF
ENDDO
55 CONTINUE
! Retrieve coeff string from array
IF( POS .gt. 0) THEN
Tljpmfstr_ = TljpmfSTR(POS)
ELSE
Tljpmfstr_ = '0.0'
ENDIF
END FUNCTION GET_Tljpmf
!
!____________________________________________________________________________________________
!
subroutine read_Pipjm
!
! read NN Pipjm coefficients
!
implicit none
!
integer :: NN,ENDOF,iNN
CHARACTER(len=100) :: buff_,filename_
INTEGER :: IT1,IT2,IT3
CHARACTER(len=100),DIMENSION(:),ALLOCATABLE :: PipjmTableNames
INTEGER :: Ntab,iTab
CHARACTER(len=T5len) :: coeffstr
!
!
! 1. Get the number of NN coefficients
!
!
OPEN(fid_,FILE=Pipjm_in)
! Count the number of files
NTab = 0
ENDOF = 0
DO WHILE(ENDOF .eq. 0)
READ(fid_,*,IOSTAT =ENDOF) buff_
NTab = NTab +1
ENDDO
!
NTab = NTab -1
!
CLOSE(fid_)
!
! Allocate array of Tljpmf table)
ALLOCATE(PipjmTableNames(1:NTab))
!
!
!
OPEN(fid_,FILE=Pipjm_in)
ENDOF = 0
iTab = 1
DO WHILE (ENDOF .eq. 0 )
READ(fid_,*,IOSTAT=ENDOF) buff_
IF (iTab .le. NTab) THEN
PipjmTableNames(iTab) = adjustl(trim(buff_))
iTab = iTab +1
ENDIF
ENDDO
CLOSE(fid_)
!
! 3. Inquire NN coefficients
!
iNN = 1
!
DO iTab =1, NTab
filename_ = trim(adjustl('./Pipjm/'//PipjmTableNames(iTab)))
OPEN(FID_,FILE = filename_)
ENDOF = 0
DO WHILE (ENDOF .eq. 0 )
READ(fid_,*,IOSTAT=ENDOF) IT1,IT2,IT3,coeffstr
IF (IT1 .le. PMMAXX .and. IT2 .le. JEmaxx .and. IT3 .le. PMMAXX) THEN
!
PipjmSTR(IT1,IT2,IT3) = coeffstr
!
!! write(*,*) PipjmSTR(IT1,IT2,IT3)
IF(ENDOF .eq. 0) iNN = iNN+1
ENDIF
ENDDO
!
CLOSE(fid_)
ENDDO
!
!
!
! Free memory
DEALLOCATE(PipjmTableNames)
!
end subroutine read_Pipjm
!______________________________________________________________________________________________
!
subroutine load_curlyNpjmn
! load GK row Vector Npj from *.csv file to character array
! Note: only important the non-zero coefficients
implicit none
!
! directory where the coefficients are stored
CHARACTER(len=1) :: dirCab ! species depend directory
! local vars.
INTEGER :: ENDOF
INTEGER :: IT1,IT2,IT3,IT4,IT5
INTEGER :: PMAXX_,JMAXX_,MMAXX_,NMAXX_,ROWMAXX_
CHARACTER(len=100),DIMENSION(:),ALLOCATABLE :: TableNames
INTEGER :: Ntab,iTab,NN,iNN
CHARACTER(len=T5len) :: buff_
CHARACTER(len=256) :: filename_
!
! Maximum indices to load
JMAXX_ = JEmaxx ! energy component
MMAXX_ = PMAXX
!! NMAXX_ = NFLRMAXX
NMAXX_ = NPFLRMAXX
ROWMAXX_ = numbb(Pmaxb,Jmaxb)
!
IF( coll_ab .eq. 'ie' .or. coll_ab .eq. 'ii' ) THEN
CALL allocate_array(curlyNpjmn_STR,0,PMMAXX,0,JEmaxx,0,PMMAXX,0,NPFLRMAXX,numbi(0,0),numbi(Pmaxi,Jmaxi))
dirCab = 'i'
ENDIF
IF( coll_ab .eq. 'ei' .or. coll_ab .eq. 'ee' ) THEN
dirCab = 'e'
CALL allocate_array(curlyNpjmn_STR,0,PMMAXX,0,JEmaxx,0,PMMAXX,0,NPFLRMAXX,numbe(0,0),numbe(Pmaxe,Jmaxe))
ENDIF
!
WRITE(filename_, *) trim(adjustl(cosolver_dir))//'/' //CurlyNpjmn_dir // '/' // dirCab // '/'// curlyNpjmn_in
!
IF( me .eq. 0) WRITE(*, *) 'Read log.in file ', trim(adjustl(filename_))
!
OPEN(fid_,FILE = trim(adjustl(filename_)))
!
! 1. Count the number of files and coeff
!
NTab = 0
ENDOF = 0
DO WHILE(ENDOF .eq. 0)
READ(fid_,*,IOSTAT =ENDOF) buff_
NTab = NTab +1
ENDDO
!
NTab = NTab -1
!
CLOSE(fid_)
!
! Allocate array of
IF( allocated(TableNames)) deallocate(TableNames)
ALLOCATE(TableNames(1:NTab))
!
OPEN(fid_,FILE= trim(adjustl(filename_)))
!
ENDOF = 0
iTab = 1
DO WHILE (ENDOF .eq. 0 )
READ(fid_,*,IOSTAT=ENDOF) buff_
IF (iTab .le. NTab) THEN
TableNames(iTab) = adjustl(trim(buff_))
iTab = iTab +1
ENDIF
ENDDO
CLOSE(fid_)
!
Nidx_p_max = 0
Nidx_j_max = 0
Nidx_n_max = 0
!
iNN =1
DO iTab=1,NTab
! name of table to read
write(filename_, *) trim(adjustl(cosolver_dir))//'/' // curlyNpjmn_dir // '/' // dirCab // '/'// TableNames(iTab)
OPEN(fid_,FILE = filename_)
ENDOF = 0
DO WHILE (ENDOF .eq. 0 )
READ(fid_,*,IOSTAT=ENDOF) IT1,IT2,IT3,IT4,IT5,buff_
IF( IT1 .gt. Nidx_p_max ) Nidx_p_max = IT1
IF( IT2 .gt. Nidx_j_max ) Nidx_j_max = IT2
IF( IT4 .gt. Nidx_n_max ) Nidx_n_max = IT4
IF (IT1 .le. PMMAXX .and. IT2 .le. JEmaxx .and. IT3 .le. MMAXX_ .and. IT4 .le. NMAXX_ .and. IT5 .le. ROWMAXX_) THEN
curlyNpjmn_STR(IT1,IT2,IT3,IT4,IT5) = buff_
ENDIF
ENDDO
CLOSE(fid_)
ENDDO
!
IF( PMMAXX .gt. Nidx_p_max .or. JEmaxx .gt. Nidx_j_max .or. NMAXX_ .gt. Nidx_n_max ) THEN
WRITE(*,*) 'Error CurlNpjmn: inquire unexisting indices', PMMAXX , '/', Nidx_p_max, ', ', JEmaxx, '/', Nidx_j_max, ', :,', NMAXX_, '/', Nidx_n_max
STOP
ENDIF
!
! Free memory
DEALLOCATE(TableNames)
!
end subroutine load_curlyNpjmn
!______________________________________________________________________________________________
!
subroutine GET_PJMNCOL(p,j,m,n,col,PJMNCOL_)
!
! construct character string p//$//j//$//m//$//n//$//col
!
implicit none
!
INTEGER,intent(in) :: p,j,m,n,col
CHARACTER(len=PJMNCOLlen),intent(inout) :: PJMNCOL_
!
CHARACTER*2 :: SIT1,SIT2,SIT3,SIT4,SIT5
! Convert int to string
WRITE(SIT1,'(i2)') p
WRITE(SIT2,'(i2)') j
WRITE(SIT3,'(i2)') m
WRITE(SIT4,'(i2)') n
WRITE(SIT5,'(i2)') col
!
PJMNCOL_ = trim(adjustl(SIT1))//'$'//trim(adjustl(SIT2))//'$'//trim(adjustl(SIT3))//'$'//trim(adjustl(SIT4))//'$'//trim(adjustl(SIT5))
!
end subroutine get_PJMNCOL
!______________________________________________________________________________________________
!
end module load_prec_coeffs
c4science · Help