X_utilities.f90

Go to the documentation of this file.

!------------------------------------------------------------------------------!
!------------------------------------------------------------------------------!
module x_utilities
#include <PB3D_macros.h>
    use str_utilities
    use messages
    use num_vars, only: dp, max_name_ln, iu, max_str_ln
    use grid_vars, only: grid_type
    use x_vars, only: n_mod_x
    
    implicit none
    
    private
    public sec_ind_loc2tot, is_necessary_x, get_sec_x_range, calc_memory_x, &
        &do_x, trim_modes
    
    ! interfaces
    
    interface sec_ind_loc2tot
        module procedure sec_ind_loc2tot_1
        module procedure sec_ind_loc2tot_2
    end interface
    
contains
    function sec_ind_loc2tot_1(id,lim_sec_X) result(res)
        ! input / output
        integer, intent(in) :: id
        integer, intent(in), optional :: lim_sec_x(2)
        integer :: res
        
        ! local variables
        integer :: lim_sec_x_loc(2)                                             ! local lim_sec_X
        
        ! set local lim_sec_X
        lim_sec_x_loc = [1,n_mod_x]
        if (present(lim_sec_x)) lim_sec_x_loc = lim_sec_x
        
        ! set sec_ind_tot
        res = lim_sec_x_loc(1)-1+id
    end function sec_ind_loc2tot_1
    function sec_ind_loc2tot_2(id,jd,lim_sec_X) result(res)
        ! input / output
        integer, intent(in) :: id
        integer, intent(in) :: jd
        integer, intent(in), optional :: lim_sec_x(2,2)
        integer :: res(2)
        
        ! local variables
        integer :: lim_sec_x_loc(2,2)                                           ! local lim_sec_X
        
        ! set local lim_sec_X
        lim_sec_x_loc(:,1) = [1,n_mod_x]
        lim_sec_x_loc(:,2) = [1,n_mod_x]
        if (present(lim_sec_x)) lim_sec_x_loc = lim_sec_x
        
        ! set sec_ind_tot
        res = [lim_sec_x_loc(1,1)-1+id,lim_sec_x_loc(1,2)-1+jd]
    end function sec_ind_loc2tot_2
    
    subroutine get_sec_x_range(sec_X_range_loc,sec_X_range_tot,m,sym,lim_sec_X)
        use x_vars, only: n_mod_x
        use num_utilities, only: c
        
        ! input / output
        integer, intent(inout) :: sec_x_range_loc(2)
        integer, intent(inout) :: sec_x_range_tot(2)
        integer, intent(in) :: m
        logical, intent(in) :: sym
        integer, intent(in), optional :: lim_sec_x(2,2)
        
        ! local variables
        integer :: k_range_loc(2)                                               ! local range of k
        integer :: k                                                            ! counter
        integer :: n_mod                                                        ! local n_mod
        
        ! set number of modes of dimension 1
        if (present(lim_sec_x)) then
            n_mod = lim_sec_x(2,1)-lim_sec_x(1,1)+1
        else
            n_mod = n_mod_x
        endif
        
        ! initialize secondary X range
        k_range_loc = [n_mod+1,0]                                               ! initialize to inverted values, out of bounds
        
        ! find start and end of k range
        find_start: do k = 1,n_mod
            if (is_necessary_x(sym,[k,m],lim_sec_x)) then                       ! first necessary pair [k,m]
                k_range_loc(1) = k
                exit find_start
            end if
        end do find_start
        find_stop: do k = n_mod,k_range_loc(1),-1
            if (is_necessary_x(sym,[k,m],lim_sec_x)) then                       ! last necessary pair [k,m]
                k_range_loc(2) = k
                exit find_stop
            end if
        end do find_stop
        
        ! translate k range and m value to 1-D index
        do k = 1,2
            sec_x_range_loc(k) = c([k_range_loc(k),m],sym,n_mod_x,lim_sec_x)
            sec_x_range_tot(k) = c(sec_ind_loc2tot(k_range_loc(k),m,lim_sec_x),&
                &sym,n_mod_x)
        end do
    end subroutine get_sec_x_range
    
    logical function do_x()
        use num_vars, only: x_jobs_lims, x_job_nr
        
        ! increment perturbation job nr.
        x_job_nr = x_job_nr + 1
        
        if (x_job_nr.le.size(x_jobs_lims,2)) then
            do_x = .true.
        else
            do_x = .false.
        end if
    end function do_x
    
    logical function is_necessary_x(sym,sec_X_id,lim_sec_X) result(res)
        ! input / output
        logical, intent(in) :: sym
        integer, intent(in) :: sec_x_id(2)
        integer, intent(in), optional :: lim_sec_x(2,2)
        
        ! local variables
        integer :: lim_sec_x_loc(2,2)                                           ! local version of lim_sec_X
        
        ! initialize res
        res = .true.
        
        ! modify res depending on symmetry
        if (sym) then
            ! set local lim_sec_X
            lim_sec_x_loc(:,1) = [1,n_mod_x]
            lim_sec_x_loc(:,2) = [1,n_mod_x]
            if (present(lim_sec_x)) lim_sec_x_loc = lim_sec_x
            
            ! set res
            if (lim_sec_x_loc(1,1)+sec_x_id(1).lt.&
                &lim_sec_x_loc(1,2)+sec_x_id(2)) res = .false.
        end if
    end function is_necessary_x
    
    integer function calc_memory_x(ord,arr_size,n_mod,mem_size) result(ierr)
        use iso_c_binding
        
        character(*), parameter :: rout_name = 'calc_memory_X'
        
        ! input / output
        integer, intent(in) :: ord
        integer, intent(in) :: arr_size
        integer, intent(in) :: n_mod
        real(dp), intent(inout) :: mem_size
        
        ! local variables
        integer(C_SIZE_T) :: dp_size                                            ! size of dp
        character(len=max_str_ln) :: err_msg                                    ! error message
        
        ! initialize ierr
        ierr = 0
        
        call lvl_ud(1)
        
        ! get size of complex variable
        dp_size = 2*sizeof(1._dp)                                               ! complex variable
        
        ! calculate memory depending on order
        select case(ord)
            case (1)                                                            ! vectorial data: U, DU
                ! set memory size
                mem_size = 4*arr_size
                mem_size = mem_size*n_mod**ord
                mem_size = mem_size*dp_size
            case (2)                                                            ! tensorial data: PV, KV
                ! set memory size
                mem_size = arr_size
                mem_size = mem_size*(2*n_mod**ord+4*n_mod*(n_mod+1)/2)
                mem_size = mem_size*dp_size
            case default
                ierr = 1
                err_msg = 'Orders > 2 are not implemented'
                chckerr(err_msg)
        end select
        
        ! convert B to MB
        mem_size = mem_size*1.e-6_dp
        
        ! apply 50% safety factor (empirical)
        mem_size = mem_size*1.5_dp
        
        ! test overflow
        if (mem_size.lt.0) then
            ierr = 1
            chckerr('Overflow occured')
        end if
        
        call lvl_ud(-1)
    end function calc_memory_x
 
    integer function trim_modes(mds_i,mds_o,id_lim_i,id_lim_o) result(ierr)
        use x_vars, only: modes_type
        
        character(*), parameter :: rout_name = 'trim_modes'
        
        ! input / output
        type(modes_type), intent(in) :: mds_i
        type(modes_type), intent(in) :: mds_o
        integer, intent(inout) :: id_lim_i(2)
        integer, intent(inout) :: id_lim_o(2)
        
        ! local variables
        integer :: m                                                            ! counter
        character(len=max_str_ln) :: err_msg                                    ! error message
        
        ! initialize ierr
        ierr = 0
        
        ! find limits point where input modes and output modes coincide
        ! (input grid should comprise output grid)
        id_lim_o = [1,size(mds_o%sec,1)]
        id_lim_i = [-1,-1]
        do m = 1,size(mds_i%sec,1)
            if (mds_i%sec(m,1).eq.mds_o%sec(id_lim_o(1),1)) then
                id_lim_i(1) = m
                exit
            end if
        end do
        do m = size(mds_i%sec,1),1,-1
            if (mds_i%sec(m,1).eq.mds_o%sec(id_lim_o(2),1)) then
                id_lim_i(2) = m
                exit
            end if
        end do
        
        ! test
        if (any(id_lim_o.le.0)) then
            ierr = 1
            err_msg = 'Cannot find limits of input modes'
            chckerr(err_msg)
        end if
    end function trim_modes
end module x_utilities