Ignore:
Timestamp:
Nov 30, 2018, 12:55:49 AM (6 years ago)
Author:
jvatant
Message:

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:
1 moved

Legend:

Unmodified
Added
Removed

  • trunk/LMDZ.COMMON/libf/misc/sort_mod.F90

    r2049 r2050  
    1 MODULE qsort_mod
     1MODULE sort_mod
    22
    3   ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    4   ! Quick-sort algorithm
    5   ! ~~~~~~~~~~~~~~~~~~~~~
    6   ! -> Sorts an array (real only for now) y by increasing values
    7   ! -> If wanted, applies the same sorting to a companion array x ( useful if you have binned xy data )
     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 )
    810  !
    9   ! TODO : Extend interface to integers
     11  ! TODO : Extend interface to integers - Enable decreasing order
    1012  !
    1113  ! Author : J. Vatant d'Ollone - 2018
    12   ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     14  ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    1315
    1416  IMPLICIT NONE
    1517
    1618  INTERFACE qsort
    17      MODULE PROCEDURE qsort_y_r
    18      MODULE PROCEDURE qsort_xy_r
     19     MODULE PROCEDURE qsort_r
     20     MODULE PROCEDURE qsort_outp_r
    1921  END INTERFACE qsort
     22
     23  INTERFACE isort
     24     MODULE PROCEDURE isort_r
     25     MODULE PROCEDURE isort_outp_r
     26  END INTERFACE isort
    2027
    2128CONTAINS
    2229
    23   RECURSIVE SUBROUTINE qsort_y_r(y,n)
     30  ! 1-pivot quicksort
     31  RECURSIVE SUBROUTINE qsort_r(A,n)
    2432
    2533    IMPLICIT NONE
    … …  
    2735    ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    2836    ! I/O
    29     integer,               intent(in)    :: n !! Size of arrays
    30     real,    dimension(n), intent(inout) :: y !! Values to sort by
     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
    3140
    3241    ! Local variables
    … …  
    3443    real    :: random
    3544    real    :: pivot
    36     real    :: tempy
     45    real    :: temp
    3746    integer :: marker
    3847    ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    … …  
    4150
    4251       call random_number(random)
    43        pivot = y(int(random*real(n-1))+1) ! random pivot (not best performance, but avoids worst-case)
     52       pivot = A(int(random*real(n-1))+1) ! random pivot (not best performance, but avoids worst-case)
    4453       left = 0
    4554       right = n + 1
    … …  
    4756       do while (left < right)
    4857          right = right - 1
    49           do while (y(right) > pivot)
     58          do while (A(right) > pivot)
    5059             right = right - 1
    5160          enddo
    5261          left = left + 1
    53           do while (y(left) < pivot)
     62          do while (A(left) < pivot)
    5463             left = left + 1
    5564          enddo
    5665          if (left < right) then
    57              tempy = y(left)
    58              y(left) = y(right)
    59              y(right) = tempy
     66             temp     = A(left)
     67             A(left)  = A(right)
     68             A(right) = temp
    6069          endif
    6170       enddo
    … …  
    6776       endif
    6877
    69        call qsort(y(:marker-1), marker-1)
    70        call qsort(y(marker:),   n-marker+1)
     78       call qsort(A(:marker-1), marker-1)
     79       call qsort(A(marker:),   n-marker+1)
    7180
    7281    endif
    7382
    74   END SUBROUTINE qsort_y_r
     83  END SUBROUTINE qsort_r
    7584
    7685
    77   RECURSIVE SUBROUTINE qsort_xy_r(y,x,n)
     86  ! 1-pivot quicksort + output array of permutations
     87  RECURSIVE SUBROUTINE qsort_outp_r(A,n,P)
    7888
    7989    IMPLICIT NONE
    … …  
    8191    ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    8292    ! I/O
    83     integer,               intent(in)    :: n !! Size of arrays
    84     real,    dimension(n), intent(inout) :: y !! Values to sort by
    85     real,    dimension(n), intent(inout) :: x !! Companion array that will be sorted same way as y
     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
    8697
    8798    ! Local variables
     99    integer :: i
    88100    integer :: left, right
    89101    real    :: random
    90102    real    :: pivot
    91     real    :: tempy
    92     real    :: tempx
     103    real    :: tempA
     104    integer :: tempP
    93105    integer :: marker
    94     ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 
     106    ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    95107
    96108    if (n > 1) then
    97109
    98110       call random_number(random)
    99        pivot = y(int(random*real(n-1))+1) ! random pivot (not best performance, but avoids worst-case)
     111       pivot = A(int(random*real(n-1))+1) ! random pivot (not best performance, but avoids worst-case)
    100112       left = 0
    101113       right = n + 1
    … …  
    103115       do while (left < right)
    104116          right = right - 1
    105           do while (y(right) > pivot)
     117          do while (A(right) > pivot)
    106118             right = right - 1
    107119          enddo
    108120          left = left + 1
    109           do while (y(left) < pivot)
     121          do while (A(left) < pivot)
    110122             left = left + 1
    111123          enddo
    112124          if (left < right) then
    113              tempy = y(left)
    114              tempx = x(left)
    115              y(left) = y(right)
    116              x(left) = x(right)
    117              y(right) = tempy
    118              x(right) = tempx
     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
    119131          endif
    120132       enddo
    … …  
    126138       endif
    127139
    128        call qsort(y(:marker-1), x(:marker-1), marker-1)
    129        call qsort(y(marker:),   x(marker:),   n-marker+1)
     140       call qsort(A(:marker-1), marker-1, P(:marker-1))
     141       call qsort(A(marker:),   n-marker+1, P(marker:))
    130142
    131143    endif
    132144
    133   END SUBROUTINE qsort_xy_r
     145  END SUBROUTINE qsort_outp_r
    134146
    135147
    136 END MODULE qsort_mod
     148  ! Insertion sort
     149  SUBROUTINE isort_r(A,n)
     150    IMPLICIT NONE
    137151
     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 TracChangeset for help on using the changeset viewer.