diff --git a/CMakeLists.txt b/CMakeLists.txt index ded2d23..f05fed3 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -1,7 +1,7 @@ cmake_minimum_required(VERSION "3.14") cmake_policy(SET CMP0148 OLD) # To suppress a warning emerging from scikit-build -project(Clima LANGUAGES Fortran C VERSION "0.4.3") +project(Clima LANGUAGES Fortran C VERSION "0.4.4") set(CMAKE_Fortran_MODULE_DIRECTORY "${CMAKE_BINARY_DIR}/modules") diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index 20f5d73..fe1e54d 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -49,7 +49,6 @@ target_link_libraries(clima h5fortran futils finterp - mrgrnk dop853 minpack ${LAPACK_LIBRARIES} diff --git a/src/dependencies/CMakeLists.txt b/src/dependencies/CMakeLists.txt index f948e1d..0471ba0 100644 --- a/src/dependencies/CMakeLists.txt +++ b/src/dependencies/CMakeLists.txt @@ -22,9 +22,9 @@ CPMAddPackage( # futils CPMAddPackage( NAME futils - VERSION 0.1.8 + VERSION 0.1.9 GITHUB_REPOSITORY "Nicholaswogan/futils" - GIT_TAG "v0.1.8" + GIT_TAG "v0.1.9" EXCLUDE_FROM_ALL ON ) @@ -63,6 +63,3 @@ CPMAddPackage( add_library(minpack ${minpack_SOURCE_DIR}/src/minpack.f90 ) - -# mrgrnk is from ORDERPACK 2.0: http://www.fortran-2000.com/rank/#1.1 -add_library(mrgrnk mrgrnk_mod.f90) diff --git a/src/dependencies/mrgrnk_mod.f90 b/src/dependencies/mrgrnk_mod.f90 deleted file mode 100644 index e0a457d..0000000 --- a/src/dependencies/mrgrnk_mod.f90 +++ /dev/null @@ -1,609 +0,0 @@ -Module mrgrnk_mod -Integer, Parameter :: kdp = selected_real_kind(15) -public :: mrgrnk -private :: kdp -private :: R_mrgrnk, I_mrgrnk, D_mrgrnk -interface mrgrnk - module procedure D_mrgrnk, R_mrgrnk, I_mrgrnk -end interface mrgrnk -contains - -Subroutine D_mrgrnk (XDONT, IRNGT) -! __________________________________________________________ -! MRGRNK = Merge-sort ranking of an array -! For performance reasons, the first 2 passes are taken -! out of the standard loop, and use dedicated coding. -! __________________________________________________________ -! __________________________________________________________ - Real (kind=kdp), Dimension (:), Intent (In) :: XDONT - Integer, Dimension (:), Intent (Out) :: IRNGT -! __________________________________________________________ - Real (kind=kdp) :: XVALA, XVALB -! - Integer, Dimension (SIZE(IRNGT)) :: JWRKT - Integer :: LMTNA, LMTNC, IRNG1, IRNG2 - Integer :: NVAL, IIND, IWRKD, IWRK, IWRKF, JINDA, IINDA, IINDB -! - NVAL = Min (SIZE(XDONT), SIZE(IRNGT)) - Select Case (NVAL) - Case (:0) - Return - Case (1) - IRNGT (1) = 1 - Return - Case Default - Continue - End Select -! -! Fill-in the index array, creating ordered couples -! - Do IIND = 2, NVAL, 2 - If (XDONT(IIND-1) <= XDONT(IIND)) Then - IRNGT (IIND-1) = IIND - 1 - IRNGT (IIND) = IIND - Else - IRNGT (IIND-1) = IIND - IRNGT (IIND) = IIND - 1 - End If - End Do - If (Modulo(NVAL, 2) /= 0) Then - IRNGT (NVAL) = NVAL - End If -! -! We will now have ordered subsets A - B - A - B - ... -! and merge A and B couples into C - C - ... -! - LMTNA = 2 - LMTNC = 4 -! -! First iteration. The length of the ordered subsets goes from 2 to 4 -! - Do - If (NVAL <= 2) Exit -! -! Loop on merges of A and B into C -! - Do IWRKD = 0, NVAL - 1, 4 - If ((IWRKD+4) > NVAL) Then - If ((IWRKD+2) >= NVAL) Exit -! -! 1 2 3 -! - If (XDONT(IRNGT(IWRKD+2)) <= XDONT(IRNGT(IWRKD+3))) Exit -! -! 1 3 2 -! - If (XDONT(IRNGT(IWRKD+1)) <= XDONT(IRNGT(IWRKD+3))) Then - IRNG2 = IRNGT (IWRKD+2) - IRNGT (IWRKD+2) = IRNGT (IWRKD+3) - IRNGT (IWRKD+3) = IRNG2 -! -! 3 1 2 -! - Else - IRNG1 = IRNGT (IWRKD+1) - IRNGT (IWRKD+1) = IRNGT (IWRKD+3) - IRNGT (IWRKD+3) = IRNGT (IWRKD+2) - IRNGT (IWRKD+2) = IRNG1 - End If - Exit - End If -! -! 1 2 3 4 -! - If (XDONT(IRNGT(IWRKD+2)) <= XDONT(IRNGT(IWRKD+3))) Cycle -! -! 1 3 x x -! - If (XDONT(IRNGT(IWRKD+1)) <= XDONT(IRNGT(IWRKD+3))) Then - IRNG2 = IRNGT (IWRKD+2) - IRNGT (IWRKD+2) = IRNGT (IWRKD+3) - If (XDONT(IRNG2) <= XDONT(IRNGT(IWRKD+4))) Then -! 1 3 2 4 - IRNGT (IWRKD+3) = IRNG2 - Else -! 1 3 4 2 - IRNGT (IWRKD+3) = IRNGT (IWRKD+4) - IRNGT (IWRKD+4) = IRNG2 - End If -! -! 3 x x x -! - Else - IRNG1 = IRNGT (IWRKD+1) - IRNG2 = IRNGT (IWRKD+2) - IRNGT (IWRKD+1) = IRNGT (IWRKD+3) - If (XDONT(IRNG1) <= XDONT(IRNGT(IWRKD+4))) Then - IRNGT (IWRKD+2) = IRNG1 - If (XDONT(IRNG2) <= XDONT(IRNGT(IWRKD+4))) Then -! 3 1 2 4 - IRNGT (IWRKD+3) = IRNG2 - Else -! 3 1 4 2 - IRNGT (IWRKD+3) = IRNGT (IWRKD+4) - IRNGT (IWRKD+4) = IRNG2 - End If - Else -! 3 4 1 2 - IRNGT (IWRKD+2) = IRNGT (IWRKD+4) - IRNGT (IWRKD+3) = IRNG1 - IRNGT (IWRKD+4) = IRNG2 - End If - End If - End Do -! -! The Cs become As and Bs -! - LMTNA = 4 - Exit - End Do -! -! Iteration loop. Each time, the length of the ordered subsets -! is doubled. -! - Do - If (LMTNA >= NVAL) Exit - IWRKF = 0 - LMTNC = 2 * LMTNC -! -! Loop on merges of A and B into C -! - Do - IWRK = IWRKF - IWRKD = IWRKF + 1 - JINDA = IWRKF + LMTNA - IWRKF = IWRKF + LMTNC - If (IWRKF >= NVAL) Then - If (JINDA >= NVAL) Exit - IWRKF = NVAL - End If - IINDA = 1 - IINDB = JINDA + 1 -! -! Shortcut for the case when the max of A is smaller -! than the min of B. This line may be activated when the -! initial set is already close to sorted. -! - IF (XDONT(IRNGT(JINDA)) <= XDONT(IRNGT(IINDB))) CYCLE -! -! One steps in the C subset, that we build in the final rank array -! -! Make a copy of the rank array for the merge iteration -! - JWRKT (1:LMTNA) = IRNGT (IWRKD:JINDA) -! - XVALA = XDONT (JWRKT(IINDA)) - XVALB = XDONT (IRNGT(IINDB)) -! - Do - IWRK = IWRK + 1 -! -! We still have unprocessed values in both A and B -! - If (XVALA > XVALB) Then - IRNGT (IWRK) = IRNGT (IINDB) - IINDB = IINDB + 1 - If (IINDB > IWRKF) Then -! Only A still with unprocessed values - IRNGT (IWRK+1:IWRKF) = JWRKT (IINDA:LMTNA) - Exit - End If - XVALB = XDONT (IRNGT(IINDB)) - Else - IRNGT (IWRK) = JWRKT (IINDA) - IINDA = IINDA + 1 - If (IINDA > LMTNA) Exit! Only B still with unprocessed values - XVALA = XDONT (JWRKT(IINDA)) - End If -! - End Do - End Do -! -! The Cs become As and Bs -! - LMTNA = 2 * LMTNA - End Do -! - Return -! -End Subroutine D_mrgrnk - -Subroutine R_mrgrnk (XDONT, IRNGT) -! __________________________________________________________ -! MRGRNK = Merge-sort ranking of an array -! For performance reasons, the first 2 passes are taken -! out of the standard loop, and use dedicated coding. -! __________________________________________________________ -! _________________________________________________________ - Real, Dimension (:), Intent (In) :: XDONT - Integer, Dimension (:), Intent (Out) :: IRNGT -! __________________________________________________________ - Real :: XVALA, XVALB -! - Integer, Dimension (SIZE(IRNGT)) :: JWRKT - Integer :: LMTNA, LMTNC, IRNG1, IRNG2 - Integer :: NVAL, IIND, IWRKD, IWRK, IWRKF, JINDA, IINDA, IINDB -! - NVAL = Min (SIZE(XDONT), SIZE(IRNGT)) - Select Case (NVAL) - Case (:0) - Return - Case (1) - IRNGT (1) = 1 - Return - Case Default - Continue - End Select -! -! Fill-in the index array, creating ordered couples -! - Do IIND = 2, NVAL, 2 - If (XDONT(IIND-1) <= XDONT(IIND)) Then - IRNGT (IIND-1) = IIND - 1 - IRNGT (IIND) = IIND - Else - IRNGT (IIND-1) = IIND - IRNGT (IIND) = IIND - 1 - End If - End Do - If (Modulo(NVAL, 2) /= 0) Then - IRNGT (NVAL) = NVAL - End If -! -! We will now have ordered subsets A - B - A - B - ... -! and merge A and B couples into C - C - ... -! - LMTNA = 2 - LMTNC = 4 -! -! First iteration. The length of the ordered subsets goes from 2 to 4 -! - Do - If (NVAL <= 2) Exit -! -! Loop on merges of A and B into C -! - Do IWRKD = 0, NVAL - 1, 4 - If ((IWRKD+4) > NVAL) Then - If ((IWRKD+2) >= NVAL) Exit -! -! 1 2 3 -! - If (XDONT(IRNGT(IWRKD+2)) <= XDONT(IRNGT(IWRKD+3))) Exit -! -! 1 3 2 -! - If (XDONT(IRNGT(IWRKD+1)) <= XDONT(IRNGT(IWRKD+3))) Then - IRNG2 = IRNGT (IWRKD+2) - IRNGT (IWRKD+2) = IRNGT (IWRKD+3) - IRNGT (IWRKD+3) = IRNG2 -! -! 3 1 2 -! - Else - IRNG1 = IRNGT (IWRKD+1) - IRNGT (IWRKD+1) = IRNGT (IWRKD+3) - IRNGT (IWRKD+3) = IRNGT (IWRKD+2) - IRNGT (IWRKD+2) = IRNG1 - End If - Exit - End If -! -! 1 2 3 4 -! - If (XDONT(IRNGT(IWRKD+2)) <= XDONT(IRNGT(IWRKD+3))) Cycle -! -! 1 3 x x -! - If (XDONT(IRNGT(IWRKD+1)) <= XDONT(IRNGT(IWRKD+3))) Then - IRNG2 = IRNGT (IWRKD+2) - IRNGT (IWRKD+2) = IRNGT (IWRKD+3) - If (XDONT(IRNG2) <= XDONT(IRNGT(IWRKD+4))) Then -! 1 3 2 4 - IRNGT (IWRKD+3) = IRNG2 - Else -! 1 3 4 2 - IRNGT (IWRKD+3) = IRNGT (IWRKD+4) - IRNGT (IWRKD+4) = IRNG2 - End If -! -! 3 x x x -! - Else - IRNG1 = IRNGT (IWRKD+1) - IRNG2 = IRNGT (IWRKD+2) - IRNGT (IWRKD+1) = IRNGT (IWRKD+3) - If (XDONT(IRNG1) <= XDONT(IRNGT(IWRKD+4))) Then - IRNGT (IWRKD+2) = IRNG1 - If (XDONT(IRNG2) <= XDONT(IRNGT(IWRKD+4))) Then -! 3 1 2 4 - IRNGT (IWRKD+3) = IRNG2 - Else -! 3 1 4 2 - IRNGT (IWRKD+3) = IRNGT (IWRKD+4) - IRNGT (IWRKD+4) = IRNG2 - End If - Else -! 3 4 1 2 - IRNGT (IWRKD+2) = IRNGT (IWRKD+4) - IRNGT (IWRKD+3) = IRNG1 - IRNGT (IWRKD+4) = IRNG2 - End If - End If - End Do -! -! The Cs become As and Bs -! - LMTNA = 4 - Exit - End Do -! -! Iteration loop. Each time, the length of the ordered subsets -! is doubled. -! - Do - If (LMTNA >= NVAL) Exit - IWRKF = 0 - LMTNC = 2 * LMTNC -! -! Loop on merges of A and B into C -! - Do - IWRK = IWRKF - IWRKD = IWRKF + 1 - JINDA = IWRKF + LMTNA - IWRKF = IWRKF + LMTNC - If (IWRKF >= NVAL) Then - If (JINDA >= NVAL) Exit - IWRKF = NVAL - End If - IINDA = 1 - IINDB = JINDA + 1 -! -! Shortcut for the case when the max of A is smaller -! than the min of B. This line may be activated when the -! initial set is already close to sorted. -! - IF (XDONT(IRNGT(JINDA)) <= XDONT(IRNGT(IINDB))) CYCLE -! -! One steps in the C subset, that we build in the final rank array -! -! Make a copy of the rank array for the merge iteration -! - JWRKT (1:LMTNA) = IRNGT (IWRKD:JINDA) -! - XVALA = XDONT (JWRKT(IINDA)) - XVALB = XDONT (IRNGT(IINDB)) -! - Do - IWRK = IWRK + 1 -! -! We still have unprocessed values in both A and B -! - If (XVALA > XVALB) Then - IRNGT (IWRK) = IRNGT (IINDB) - IINDB = IINDB + 1 - If (IINDB > IWRKF) Then -! Only A still with unprocessed values - IRNGT (IWRK+1:IWRKF) = JWRKT (IINDA:LMTNA) - Exit - End If - XVALB = XDONT (IRNGT(IINDB)) - Else - IRNGT (IWRK) = JWRKT (IINDA) - IINDA = IINDA + 1 - If (IINDA > LMTNA) Exit! Only B still with unprocessed values - XVALA = XDONT (JWRKT(IINDA)) - End If -! - End Do - End Do -! -! The Cs become As and Bs -! - LMTNA = 2 * LMTNA - End Do -! - Return -! -End Subroutine R_mrgrnk -Subroutine I_mrgrnk (XDONT, IRNGT) -! __________________________________________________________ -! MRGRNK = Merge-sort ranking of an array -! For performance reasons, the first 2 passes are taken -! out of the standard loop, and use dedicated coding. -! __________________________________________________________ -! __________________________________________________________ - Integer, Dimension (:), Intent (In) :: XDONT - Integer, Dimension (:), Intent (Out) :: IRNGT -! __________________________________________________________ - Integer :: XVALA, XVALB -! - Integer, Dimension (SIZE(IRNGT)) :: JWRKT - Integer :: LMTNA, LMTNC, IRNG1, IRNG2 - Integer :: NVAL, IIND, IWRKD, IWRK, IWRKF, JINDA, IINDA, IINDB -! - NVAL = Min (SIZE(XDONT), SIZE(IRNGT)) - Select Case (NVAL) - Case (:0) - Return - Case (1) - IRNGT (1) = 1 - Return - Case Default - Continue - End Select -! -! Fill-in the index array, creating ordered couples -! - Do IIND = 2, NVAL, 2 - If (XDONT(IIND-1) <= XDONT(IIND)) Then - IRNGT (IIND-1) = IIND - 1 - IRNGT (IIND) = IIND - Else - IRNGT (IIND-1) = IIND - IRNGT (IIND) = IIND - 1 - End If - End Do - If (Modulo(NVAL, 2) /= 0) Then - IRNGT (NVAL) = NVAL - End If -! -! We will now have ordered subsets A - B - A - B - ... -! and merge A and B couples into C - C - ... -! - LMTNA = 2 - LMTNC = 4 -! -! First iteration. The length of the ordered subsets goes from 2 to 4 -! - Do - If (NVAL <= 2) Exit -! -! Loop on merges of A and B into C -! - Do IWRKD = 0, NVAL - 1, 4 - If ((IWRKD+4) > NVAL) Then - If ((IWRKD+2) >= NVAL) Exit -! -! 1 2 3 -! - If (XDONT(IRNGT(IWRKD+2)) <= XDONT(IRNGT(IWRKD+3))) Exit -! -! 1 3 2 -! - If (XDONT(IRNGT(IWRKD+1)) <= XDONT(IRNGT(IWRKD+3))) Then - IRNG2 = IRNGT (IWRKD+2) - IRNGT (IWRKD+2) = IRNGT (IWRKD+3) - IRNGT (IWRKD+3) = IRNG2 -! -! 3 1 2 -! - Else - IRNG1 = IRNGT (IWRKD+1) - IRNGT (IWRKD+1) = IRNGT (IWRKD+3) - IRNGT (IWRKD+3) = IRNGT (IWRKD+2) - IRNGT (IWRKD+2) = IRNG1 - End If - Exit - End If -! -! 1 2 3 4 -! - If (XDONT(IRNGT(IWRKD+2)) <= XDONT(IRNGT(IWRKD+3))) Cycle -! -! 1 3 x x -! - If (XDONT(IRNGT(IWRKD+1)) <= XDONT(IRNGT(IWRKD+3))) Then - IRNG2 = IRNGT (IWRKD+2) - IRNGT (IWRKD+2) = IRNGT (IWRKD+3) - If (XDONT(IRNG2) <= XDONT(IRNGT(IWRKD+4))) Then -! 1 3 2 4 - IRNGT (IWRKD+3) = IRNG2 - Else -! 1 3 4 2 - IRNGT (IWRKD+3) = IRNGT (IWRKD+4) - IRNGT (IWRKD+4) = IRNG2 - End If -! -! 3 x x x -! - Else - IRNG1 = IRNGT (IWRKD+1) - IRNG2 = IRNGT (IWRKD+2) - IRNGT (IWRKD+1) = IRNGT (IWRKD+3) - If (XDONT(IRNG1) <= XDONT(IRNGT(IWRKD+4))) Then - IRNGT (IWRKD+2) = IRNG1 - If (XDONT(IRNG2) <= XDONT(IRNGT(IWRKD+4))) Then -! 3 1 2 4 - IRNGT (IWRKD+3) = IRNG2 - Else -! 3 1 4 2 - IRNGT (IWRKD+3) = IRNGT (IWRKD+4) - IRNGT (IWRKD+4) = IRNG2 - End If - Else -! 3 4 1 2 - IRNGT (IWRKD+2) = IRNGT (IWRKD+4) - IRNGT (IWRKD+3) = IRNG1 - IRNGT (IWRKD+4) = IRNG2 - End If - End If - End Do -! -! The Cs become As and Bs -! - LMTNA = 4 - Exit - End Do -! -! Iteration loop. Each time, the length of the ordered subsets -! is doubled. -! - Do - If (LMTNA >= NVAL) Exit - IWRKF = 0 - LMTNC = 2 * LMTNC -! -! Loop on merges of A and B into C -! - Do - IWRK = IWRKF - IWRKD = IWRKF + 1 - JINDA = IWRKF + LMTNA - IWRKF = IWRKF + LMTNC - If (IWRKF >= NVAL) Then - If (JINDA >= NVAL) Exit - IWRKF = NVAL - End If - IINDA = 1 - IINDB = JINDA + 1 -! -! Shortcut for the case when the max of A is smaller -! than the min of B. This line may be activated when the -! initial set is already close to sorted. -! - IF (XDONT(IRNGT(JINDA)) <= XDONT(IRNGT(IINDB))) CYCLE -! -! One steps in the C subset, that we build in the final rank array -! -! Make a copy of the rank array for the merge iteration -! - JWRKT (1:LMTNA) = IRNGT (IWRKD:JINDA) -! - XVALA = XDONT (JWRKT(IINDA)) - XVALB = XDONT (IRNGT(IINDB)) -! - Do - IWRK = IWRK + 1 -! -! We still have unprocessed values in both A and B -! - If (XVALA > XVALB) Then - IRNGT (IWRK) = IRNGT (IINDB) - IINDB = IINDB + 1 - If (IINDB > IWRKF) Then -! Only A still with unprocessed values - IRNGT (IWRK+1:IWRKF) = JWRKT (IINDA:LMTNA) - Exit - End If - XVALB = XDONT (IRNGT(IINDB)) - Else - IRNGT (IWRK) = JWRKT (IINDA) - IINDA = IINDA + 1 - If (IINDA > LMTNA) Exit! Only B still with unprocessed values - XVALA = XDONT (JWRKT(IINDA)) - End If -! - End Do - End Do -! -! The Cs become As and Bs -! - LMTNA = 2 * LMTNA - End Do -! - Return -! -End Subroutine I_mrgrnk -end module diff --git a/src/radtran/clima_radtran_radiate.f90 b/src/radtran/clima_radtran_radiate.f90 index 2d5abf3..172f54d 100644 --- a/src/radtran/clima_radtran_radiate.f90 +++ b/src/radtran/clima_radtran_radiate.f90 @@ -440,7 +440,7 @@ subroutine k_aee(op, kset, zenith_u, zenith_weights, surface_albedo, surface_emi subroutine k_rorr(op, kset, zenith_u, zenith_weights, surface_albedo, surface_emissivity, cols, rw, rz) use futils, only: rebin - use mrgrnk_mod, only: mrgrnk + use futils_mrgrnk, only: mrgrnk use clima_radtran_types, only: OpticalProperties, Ksettings use clima_radtran_types, only: RadiateZWrk, RadiateXSWrk