Multidip  1.0
Multi-photon matrix elements
multidip_params Module Reference

Hard-coded parameters of MULTIDIP. More...

Functions/Subroutines

subroutine read_input_namelist (input, order, lusct, lukmt, lubnd, nkset, polar, omega, verbose, rmt_data, first_IP, rasym, raw, erange, mpiio)
 

Variables

integer, parameter ntermsasy = 5
 Number of terms in expansion of the Coulomb-Hankel functions. More...
 
integer, parameter ntermsppt = 3
 Number of terms in expansion of the exponential integral. More...
 
integer, parameter nmaxphotons = bit_size(0)
 The limit on number of photons is nested_cgreen_integ More...
 
integer, parameter max_romb_level = 20
 Maximal nesting level for Romberg integration. More...
 
integer, parameter max_levin_level = 20
 Maximal nesting level for Levin integration. More...
 
integer, parameter cheb_order = 5
 Order of Chebyshev interpolation used in Levin quadrature. More...
 
integer maxtarg = 0
 Maximal number of target states to calculate dipoles for (0 = all) More...
 
logical cache_integrals = .false.
 Cache Coulomb-Green integrals in memory (but disables some threading) More...
 
logical check_dipoles = .true.
 Check that all dipole matrices in molecular_data are nonzero. More...
 
logical closed_interm = .false.
 Consider weakly closed channel in intermediate states (unfinished!) More...
 
logical coulomb_check = .true.
 Skip integrals that cannot be asymptotically approximated well. More...
 
logical print_warnings = .true.
 Print warnings about non-converged integrals. More...
 
logical custom_kmatrices = .false.
 Ignore RSOLVE K-matrices and calculate them from scratch. More...
 
logical extend_istate = .false.
 Continue initial state from the known inner region expansion outwards. More...
 
integer num_integ_algo = 2
 Numerical integration mode (1 = Romberg, 2 = Levin) More...
 
real(wp) epsrel = 1e-6
 Romberg integration relative tolerance for convergence. More...
 
real(wp) coultol = 1e-4
 Coulomb matching relative tolerance (decides whether to use asym. forms) More...
 
real(wp) closed_range = 5.0
 Energy interval (a.u.) before threshold for inclusion of closed channels. More...
 
real(wp), parameter alpha = 1/137.03599907_wp
 Fine structure constant. More...
 
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), dimension(3), parameter compt = ['x', 'y', 'z']
 
integer, dimension(3), parameter carti = [3, 1, 2]
 Position of a Cartesian coordinate in the real spherical basis (y, z, x). More...
 
integer, dimension(3), parameter cartm = [+1, -1, 0]
 Real spherical harmonic m-value corresponding to given a Cartesian coord. More...
 

Detailed Description

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.

Function/Subroutine Documentation

◆ read_input_namelist()

subroutine multidip_params::read_input_namelist ( integer, intent(in)  input,
integer, intent(inout)  order,
integer, dimension(8), intent(inout)  lusct,
integer, dimension(8), intent(inout)  lukmt,
integer, intent(inout)  lubnd,
integer, dimension(8), intent(inout)  nkset,
complex(wp), dimension(3, nmaxphotons), intent(inout)  polar,
real(wp), dimension(nmaxphotons), intent(inout)  omega,
logical, intent(inout)  verbose,
character(len=256), intent(inout)  rmt_data,
real(wp), intent(inout)  first_IP,
real(wp), intent(inout)  rasym,
character(len=256), intent(inout)  raw,
integer, dimension(2), intent(inout)  erange,
logical, intent(inout)  mpiio 
)

Definition at line 72 of file multidip_params.F90.

Here is the caller graph for this function:

Variable Documentation

◆ alpha

real(wp), parameter multidip_params::alpha = 1/137.03599907_wp

Fine structure constant.

Definition at line 55 of file multidip_params.F90.

◆ cache_integrals

logical multidip_params::cache_integrals = .false.

Cache Coulomb-Green integrals in memory (but disables some threading)

Definition at line 43 of file multidip_params.F90.

◆ carti

integer, dimension(3), parameter multidip_params::carti = [3, 1, 2]

Position of a Cartesian coordinate in the real spherical basis (y, z, x).

Definition at line 67 of file multidip_params.F90.

◆ cartm

integer, dimension(3), parameter multidip_params::cartm = [+1, -1, 0]

Real spherical harmonic m-value corresponding to given a Cartesian coord.

Definition at line 68 of file multidip_params.F90.

◆ cheb_order

integer, parameter multidip_params::cheb_order = 5

Order of Chebyshev interpolation used in Levin quadrature.

Definition at line 40 of file multidip_params.F90.

◆ check_dipoles

logical multidip_params::check_dipoles = .true.

Check that all dipole matrices in molecular_data are nonzero.

Definition at line 44 of file multidip_params.F90.

◆ closed_interm

logical multidip_params::closed_interm = .false.

Consider weakly closed channel in intermediate states (unfinished!)

Definition at line 45 of file multidip_params.F90.

◆ closed_range

real(wp) multidip_params::closed_range = 5.0

Energy interval (a.u.) before threshold for inclusion of closed channels.

Definition at line 53 of file multidip_params.F90.

◆ compt

character(len=1), dimension(3), parameter multidip_params::compt = ['x', 'y', 'z']

Definition at line 65 of file multidip_params.F90.

◆ cone

complex(wp), parameter multidip_params::cone = 1

Definition at line 62 of file multidip_params.F90.

◆ coulomb_check

logical multidip_params::coulomb_check = .true.

Skip integrals that cannot be asymptotically approximated well.

Definition at line 46 of file multidip_params.F90.

◆ coultol

real(wp) multidip_params::coultol = 1e-4

Coulomb matching relative tolerance (decides whether to use asym. forms)

Definition at line 52 of file multidip_params.F90.

◆ custom_kmatrices

logical multidip_params::custom_kmatrices = .false.

Ignore RSOLVE K-matrices and calculate them from scratch.

Definition at line 48 of file multidip_params.F90.

◆ czero

complex(wp), parameter multidip_params::czero = 0

Definition at line 61 of file multidip_params.F90.

◆ epsrel

real(wp) multidip_params::epsrel = 1e-6

Romberg integration relative tolerance for convergence.

Definition at line 51 of file multidip_params.F90.

◆ extend_istate

logical multidip_params::extend_istate = .false.

Continue initial state from the known inner region expansion outwards.

Definition at line 49 of file multidip_params.F90.

◆ imu

complex(wp), parameter multidip_params::imu = (0, 1)

Definition at line 63 of file multidip_params.F90.

◆ max_levin_level

integer, parameter multidip_params::max_levin_level = 20

Maximal nesting level for Levin integration.

Definition at line 39 of file multidip_params.F90.

◆ max_romb_level

integer, parameter multidip_params::max_romb_level = 20

Maximal nesting level for Romberg integration.

Definition at line 38 of file multidip_params.F90.

◆ maxtarg

integer multidip_params::maxtarg = 0

Maximal number of target states to calculate dipoles for (0 = all)

Definition at line 42 of file multidip_params.F90.

◆ nmaxphotons

integer, parameter multidip_params::nmaxphotons = bit_size(0)

The limit on number of photons is nested_cgreen_integ

Definition at line 37 of file multidip_params.F90.

◆ ntermsasy

integer, parameter multidip_params::ntermsasy = 5

Number of terms in expansion of the Coulomb-Hankel functions.

Definition at line 35 of file multidip_params.F90.

◆ ntermsppt

integer, parameter multidip_params::ntermsppt = 3

Number of terms in expansion of the exponential integral.

Definition at line 36 of file multidip_params.F90.

◆ num_integ_algo

integer multidip_params::num_integ_algo = 2

Numerical integration mode (1 = Romberg, 2 = Levin)

Definition at line 50 of file multidip_params.F90.

◆ pi

real(wp), parameter multidip_params::pi = 4*atan(1.0_wp)

Definition at line 59 of file multidip_params.F90.

◆ print_warnings

logical multidip_params::print_warnings = .true.

Print warnings about non-converged integrals.

Definition at line 47 of file multidip_params.F90.

◆ rhalf

real(wp), parameter multidip_params::rhalf = 0.5

Definition at line 58 of file multidip_params.F90.

◆ rone

real(wp), parameter multidip_params::rone = 1

Definition at line 57 of file multidip_params.F90.

◆ rzero

real(wp), parameter multidip_params::rzero = 0

Definition at line 56 of file multidip_params.F90.