Line data Source code
1 : ! ***************************************************************************************************
2 : ! Copyright (C) 2020-2023 Green-X library
3 : ! This file is distributed under the terms of the APACHE2 License.
4 : !
5 : ! ***************************************************************************************************
6 : !> \brief This module contains auxiliary procedures and data structures for the main minimax routines
7 : ! ***************************************************************************************************
8 : module minimax_utils
9 : #if defined HAVE_CONFIG_H
10 : #include "config.h"
11 : #endif
12 : #include "abi_common.h"
13 : #include "gx_common.h"
14 : use defs_basis, only: dp
15 : use m_errors
16 : !use kinds, only: dp
17 : implicit none
18 :
19 : private
20 :
21 : type, public :: er_aw_aux
22 : ! Sorted array of the energy ranges
23 : real(kind=dp), dimension(:), allocatable :: energy_range
24 : ! Matrices with coefficients and weights per energy region
25 : real(kind=dp), dimension(:, :), allocatable :: aw_erange_matrix
26 : contains
27 : procedure :: get_coeff_weight => coeffs_and_weights
28 : end type er_aw_aux
29 :
30 : !> Transformation types
31 : integer, parameter, public :: cosine_tw = 1
32 : integer, parameter, public :: cosine_wt = 2
33 : integer, parameter, public :: sine_tw = 3
34 :
35 : contains
36 :
37 : !> \brief Find first element in unsorted array that is strictly greater than a given value
38 : !> This algorithm is O(n), difficult to do better with unsorted arrays
39 : !! @param[in] lenght - lenght of sorted array
40 : !! @param[in] einter - sorted array of the energy intervals
41 : !! @param[in] eval - the energy value
42 4 : function find_erange(length, einter, eval) result(idx)
43 : integer, intent(in) :: length
44 : real(kind=dp), dimension(length), intent(in) :: einter
45 : real(kind=dp), intent(in) :: eval
46 : integer :: idx
47 :
48 : ! Auxiliary variables
49 : integer :: jdx
50 : real(kind=dp) :: tmp_min_max
51 :
52 : ! Begin work
53 4 : tmp_min_max = huge(0.0_dp)
54 4 : idx = length + 1
55 :
56 113 : do jdx = 1, length
57 113 : if (eval < einter(jdx) .and. einter(jdx) < tmp_min_max) then
58 109 : idx = jdx
59 109 : tmp_min_max = einter(jdx)
60 : end if
61 : end do
62 :
63 4 : end function find_erange
64 :
65 : !> \brief Selects the energy region and scales weights and coefficients
66 : !! @param[in] grid_size - the grid size
67 : !! @param[in] bup - length of the energy region array
68 : !! @param[in] e_range - the selected energy range
69 : !! @param[inout] e_ratio - an heuristic correction factor
70 : !! @param[inout] ac_we - vector containing coefficients and weights
71 4 : subroutine coeffs_and_weights(this, grid_size, bup, e_range, ac_we, e_ratio)
72 : class(er_aw_aux), intent(in) :: this
73 : integer, intent(in) :: grid_size
74 : integer, intent(in) :: bup
75 : real(kind=dp), intent(in) :: e_range
76 : real(kind=dp), dimension(:), intent(inout) :: ac_we
77 : real(kind=dp), intent(inout) :: e_ratio
78 :
79 : ! Internal variables
80 : integer :: ien
81 :
82 : ! Select energy region
83 4 : ien = find_erange(bup, this%energy_range, e_range)
84 :
85 : ! Scale grids for large sizes when erange falls in the first energy range
86 4 : if (ien == 1 .and. grid_size > 20) then
87 0 : e_ratio = this%energy_range(1) / e_range
88 0 : if (e_ratio > 1.5_dp) then
89 0 : e_ratio = e_ratio / 1.5_dp
90 : endif
91 : end if
92 :
93 124 : ac_we(:) = this%aw_erange_matrix(:, ien)
94 :
95 4 : end subroutine coeffs_and_weights
96 :
97 0 : end module minimax_utils
|