source: trunk/LMDZ.COMMON/libf/misc/sort_mod.F90 @ 2911

Last change on this file since 2911 was 2050, checked in by jvatant, 7 years ago

Major radiative transfer contribution : Add the 'corrk_recombin' option that allows to use
correlated-k for single species instead of pre-mix and enables more flexiblity for variable species.
-> Algorithm inspired from Lacis and Oinas 1991 and Amundsen et al 2016

+ Added 'recombin_corrk.F90' - Important improvements in 'sugas_corrk.90' and 'callcork.F90'

  • Must have the desired variable species as tracers -> TBD : Enable to force composition even if no tracers
  • To have decent CPU time recombining is not done on all gridpoints and wavelenghts but we calculate a gasi/v_recomb variable on the reference corrk-k T,P grid (only for T,P values who match the atmospheric conditions ) which is then processed as a standard pre-mix in optci/v routines, but updated every time tracers on the ref P grid have varied > 1%.


READ CAREFULY :

  • In case of 'corrk_recombin', the variable L_NREFVAR doesn't have the same meaning as before and doesn't stand for the different mixing ratios but the different species.
  • Input corr-k should be found in corrkdir within 'corrk_gcm_IR/VI_XXX.dat' and can contain a 'fixed' specie ( compulsory if you include self-broadening ) that MUST have been created with correct mixing ratios, or a variable specie for which mixing ratio MUST have been set to 1 ( no self-broadening then, assume it's a trace sepecie ) -> You can't neither have CIA of variable species included upstream in the corr-k
File size: 5.1 KB
Line 
1MODULE sort_mod
2
3  ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
4  ! Sorting algorithms
5  ! ~~~~~~~~~~~~~~~~~~~
6  ! Contains quicksort (qsort) and insertion sort (isort)
7  ! -> Sorts an array (real only for now) by increasing values
8  ! -> If wanted, output indexes of permutations if you want to apply the same sorting to other arrays
9  !          ( this is quicker than having other arrays going through the whole process )
10  !
11  ! TODO : Extend interface to integers - Enable decreasing order
12  !
13  ! Author : J. Vatant d'Ollone - 2018
14  ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
15
16  IMPLICIT NONE
17
18  INTERFACE qsort
19     MODULE PROCEDURE qsort_r
20     MODULE PROCEDURE qsort_outp_r
21  END INTERFACE qsort
22
23  INTERFACE isort
24     MODULE PROCEDURE isort_r
25     MODULE PROCEDURE isort_outp_r
26  END INTERFACE isort
27
28CONTAINS
29
30  ! 1-pivot quicksort
31  RECURSIVE SUBROUTINE qsort_r(A,n)
32
33    IMPLICIT NONE
34
35    ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
36    ! I/O
37    !logical,               intent(in)    :: o !! Order T=increasing/F=decreasing
38    integer,               intent(in)    :: n !! Size of array
39    real,    dimension(n), intent(inout) :: A !! Array to sort
40
41    ! Local variables
42    integer :: left, right
43    real    :: random
44    real    :: pivot
45    real    :: temp
46    integer :: marker
47    ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
48
49    if (n > 1) then
50
51       call random_number(random)
52       pivot = A(int(random*real(n-1))+1) ! random pivot (not best performance, but avoids worst-case)
53       left = 0
54       right = n + 1
55
56       do while (left < right)
57          right = right - 1
58          do while (A(right) > pivot)
59             right = right - 1
60          enddo
61          left = left + 1
62          do while (A(left) < pivot)
63             left = left + 1
64          enddo
65          if (left < right) then
66             temp     = A(left)
67             A(left)  = A(right)
68             A(right) = temp
69          endif
70       enddo
71
72       if (left == right) then
73          marker = left + 1
74       else
75          marker = left
76       endif
77
78       call qsort(A(:marker-1), marker-1)
79       call qsort(A(marker:),   n-marker+1)
80
81    endif
82
83  END SUBROUTINE qsort_r
84
85
86  ! 1-pivot quicksort + output array of permutations
87  RECURSIVE SUBROUTINE qsort_outp_r(A,n,P)
88
89    IMPLICIT NONE
90
91    ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
92    ! I/O
93    !logical,               intent(in)    :: o !! Order T=increasing/F=decreasing
94    integer,               intent(in)    :: n !! Size of array
95    integer, dimension(n), intent(inout) :: P !! Output array of permutations
96    real,    dimension(n), intent(inout) :: A !! Array to sort
97
98    ! Local variables
99    integer :: i
100    integer :: left, right
101    real    :: random
102    real    :: pivot
103    real    :: tempA
104    integer :: tempP
105    integer :: marker
106    ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
107
108    if (n > 1) then
109
110       call random_number(random)
111       pivot = A(int(random*real(n-1))+1) ! random pivot (not best performance, but avoids worst-case)
112       left = 0
113       right = n + 1
114
115       do while (left < right)
116          right = right - 1
117          do while (A(right) > pivot)
118             right = right - 1
119          enddo
120          left = left + 1
121          do while (A(left) < pivot)
122             left = left + 1
123          enddo
124          if (left < right) then
125             tempA    = A(left)
126             tempP    = P(left)
127             A(left)  = A(right)
128             P(left)  = P(right)
129             A(right) = tempA
130             P(right) = tempP
131          endif
132       enddo
133
134       if (left == right) then
135          marker = left + 1
136       else
137          marker = left
138       endif
139
140       call qsort(A(:marker-1), marker-1, P(:marker-1))
141       call qsort(A(marker:),   n-marker+1, P(marker:))
142
143    endif
144
145  END SUBROUTINE qsort_outp_r
146
147
148  ! Insertion sort
149  SUBROUTINE isort_r(A,n)
150    IMPLICIT NONE
151
152    INTEGER,               INTENT(in)    :: n !! Size of array
153    REAL,    DIMENSION(n), INTENT(inout) :: A !! Array to sort
154
155    INTEGER ::  i, j
156    REAL    ::  x
157
158    DO i = 2, n
159       x = a(i)
160       j = i - 1
161       DO WHILE (j >= 1)
162          IF (a(j) <= x) EXIT
163          a(j + 1) = a(j)
164          j = j - 1
165       ENDDO
166       a(j + 1) = x
167    ENDDO
168  END SUBROUTINE isort_r
169
170
171  ! Insertion sort + output array of permutations
172  SUBROUTINE isort_outp_r(A,n,P)
173    IMPLICIT NONE
174
175    INTEGER,               INTENT(in)    :: n !! Size of array
176    INTEGER, DIMENSION(n), INTENT(inout) :: P !! Output array of permutations
177    REAL,    DIMENSION(n), INTENT(inout) :: A !! Array to sort
178
179    INTEGER ::  i, j, px
180    REAL    ::  x
181
182    DO i = 2, n
183       x = a(i)
184       px = p(i)
185       j = i - 1
186       DO WHILE (j >= 1)
187          IF (a(j) <= x) EXIT
188          a(j + 1) = a(j)
189          p(j + 1) = p(j)
190          j = j - 1
191       ENDDO
192       a(j + 1) = x
193       p(j + 1) = px
194    ENDDO
195  END SUBROUTINE isort_outp_r
196
197
198END MODULE sort_mod
199
Note: See TracBrowser for help on using the repository browser.