multidip_params.F90

Go to the documentation of this file.

! Copyright 2020
!
! For a comprehensive list of the developers that contributed to these codes
! see the UK-AMOR website.
!
! This file is part of UKRmol-out (UKRmol+ suite).
!
!     UKRmol-out is free software: you can redistribute it and/or modify
!     it under the terms of the GNU General Public License as published by
!     the Free Software Foundation, either version 3 of the License, or
!     (at your option) any later version.
!
!     UKRmol-out is distributed in the hope that it will be useful,
!     but WITHOUT ANY WARRANTY; without even the implied warranty of
!     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
!     GNU General Public License for more details.
!
!     You should have received a copy of the GNU General Public License
!     along with  UKRmol-out (in source/COPYING). Alternatively, you can also visit
!     <https://www.gnu.org/licenses/>.
!
!> \brief   Hard-coded parameters of MULTIDIP
!> \author  J Benda
!> \date    2020 - 2024
!>
!> This module contains mathematical and physical, as well as algorithmical parameters. The latter ones may make it
!> to the input namelist some day. At the moment their casual use is considered advanced.
!>
module multidip_params
 
    use precisn_gbl, only: wp
 
    implicit none
 
    integer,  parameter :: ntermsasy = 5                !< Number of terms in expansion of the Coulomb-Hankel functions
    integer,  parameter :: ntermsppt = 3                !< Number of terms in expansion of the exponential integral
    integer,  parameter :: nmaxphotons = bit_size(0)    !< The limit on number of photons is `nested_cgreen_integ`
    integer,  parameter :: max_romb_level = 20          !< Maximal nesting level for Romberg integration
    integer,  parameter :: max_levin_level = 20         !< Maximal nesting level for Levin integration
    integer,  parameter :: cheb_order = 5               !< Order of Chebyshev interpolation used in Levin quadrature
 
    integer             :: maxtarg = 0                  !< Maximal number of target states to calculate dipoles for (0 = all)
    logical             :: cache_integrals = .false.    !< Cache Coulomb-Green integrals in memory (but disables some threading)
    logical             :: check_dipoles = .true.       !< Check that all dipole matrices in molecular_data are nonzero
    logical             :: closed_interm = .false.      !< Consider weakly closed channel in intermediate states (unfinished!)
    logical             :: coulomb_check = .true.       !< Skip integrals that cannot be asymptotically approximated well
    logical             :: print_warnings = .true.      !< Print warnings about non-converged integrals
    logical             :: custom_kmatrices = .false.   !< Ignore RSOLVE K-matrices and calculate them from scratch
    logical             :: extend_istate = .false.      !< Continue initial state from the known inner region expansion outwards
    logical             :: ion_ion_analytic = .true.    !< Integrate one-dimensional ion-ion integrals using a closed-form formula
    integer             :: num_integ_algo = 2           !< Numerical integration mode (1 = Romberg, 2 = Levin)
    real(wp)            :: epsrel = 1e-6                !< Romberg integration relative tolerance for convergence
    real(wp)            :: coultol = 1e-4               !< Coulomb matching relative tolerance (decides whether to use asym. forms)
    real(wp)            :: closed_range = 5.0           !< Energy interval (a.u.) before threshold for inclusion of closed channels
 
    real(wp), parameter :: alpha = 1/137.03599907_wp    !< Fine structure constant
    real(wp), parameter :: rzero = 0
    real(wp), parameter :: rone = 1
    real(wp), parameter :: rhalf = 0.5
    real(wp), parameter :: pi = 4*atan(1.0_wp)
 
    complex(wp), parameter :: czero = 0
    complex(wp), parameter :: cone = 1
    complex(wp), parameter :: imu = (0, 1)
 
    character(len=1), parameter :: compt(3) = ['x', 'y', 'z']
 
    integer, parameter :: carti(3) = [3, 1, 2]          !< Position of a Cartesian coordinate in the real spherical basis (y, z, x).
    integer, parameter :: cartm(3) = [+1, -1, 0]        !< Real spherical harmonic m-value corresponding to given a Cartesian coord.
 
contains
 
    subroutine read_input_namelist (input, order, lusct, lukmt, lubnd, nkset, polar, omega, verbose, rmt_data, first_IP, rasym, &
                                    raw, erange, mpiio, gpu, lab_polar, lu_pw_dipoles)
 
        integer,     intent(in)    :: input
        integer,     intent(inout) :: order, lusct(8), lukmt(8), lubnd, nkset(8), erange(2), lab_polar(nMaxPhotons), lu_pw_dipoles
        logical,     intent(inout) :: verbose, mpiio, gpu
        real(wp),    intent(inout) :: omega(nMaxPhotons), first_IP, rasym
        complex(wp), intent(inout) :: polar(3, nMaxPhotons)
        character(len=256), intent(inout) :: rmt_data, raw
 
        namelist /mdip/ order, lusct, lukmt, lubnd, nkset, polar, omega, verbose, rmt_data, first_ip, rasym, raw, erange, mpiio, &
                        maxtarg, check_dipoles, closed_interm, coulomb_check, print_warnings, custom_kmatrices, gpu, &
                        num_integ_algo, epsrel, coultol, closed_range, extend_istate, cache_integrals, ion_ion_analytic, &
                        lab_polar, lu_pw_dipoles
 
        order = 1
        lusct = 0
        lukmt = 0
        lubnd = 0
        nkset = 1
        polar = 0
        omega = 0
        rasym = 100
        first_ip = 0
        erange = 0
        verbose = .false.
        rmt_data = 'molecular_data'
        raw = ''
        mpiio = .false.
        gpu = .false.
        lab_polar = 0
        lu_pw_dipoles = 410
 
        read (input, nml = mdip)
 
        if (order < 1) then
            print '(A)', 'Order must be positive'
            stop 1
        end if
 
        if (order > nmaxphotons) then
            print '(A,I0)', 'Order must be <= ', nmaxphotons
            stop 1
        end if
 
        if (raw /= '' .and. raw /= 'sph' .and. raw /= 'xyz' .and. raw /= 'both') then
            print '(A)', 'The parameter "raw" must be empty or one of: "sph", "xyz", "both".'
            stop 1
        end if
 
#ifndef WITH_SCALAPACK
        if (mpiio) then
            print '(A)', 'MPI-IO mode requires compiling MULTIDIP with ScaLAPACK support'
            stop 1
        end if
#endif
 
        polar(:, order + 1:nmaxphotons) = 0
        omega(order + 1:nmaxphotons) = 0
 
    end subroutine read_input_namelist
 
end module multidip_params