Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Separate cloud optics kernels #302

Merged
merged 9 commits into from
Nov 9, 2024
Merged

Separate cloud optics kernels #302

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
9 commits
Select commit Hold shift + click to select a range
39ab878
Cloud optics kernels in separate file, added to APIs
RobertPincus Aug 12, 2024
61d0208
Update cloud kernel names
RobertPincus Aug 12, 2024
463d3d1
Unsaved changes to cloud kernel name
RobertPincus Aug 12, 2024
24bcf73
Fortran API syntax
RobertPincus Aug 13, 2024
bab335f
Syntax
RobertPincus Aug 13, 2024
c9bf841
Update CI with changes from develop
RobertPincus Aug 15, 2024
6894049
Merge branch 'develop' into feature-optics-kernels
RobertPincus Oct 29, 2024
586ef4a
More recent nvfortran compiler for CI
RobertPincus Nov 7, 2024
7a49783
Update CUDA toolkit on Levante
RobertPincus Nov 7, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
5 changes: 3 additions & 2 deletions .gitlab/levante.yml
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -30,7 +30,8 @@ variables:
.nvhpc:
variables:
# Core variables:
FC: /sw/spack-levante/nvhpc-22.5-v4oky3/Linux_x86_64/22.5/compilers/bin/nvfortran
# FC: /sw/spack-levante/nvhpc-22.5-v4oky3/Linux_x86_64/22.5/compilers/bin/nvfortran
FC: /sw/spack-levante/nvhpc-24.9-p7iohv/Linux_x86_64/24.9/compilers/bin/nvfortran
# Convenience variables:
VERSION_FCFLAGS: --version
NFHOME: /sw/spack-levante/netcdf-fortran-4.5.4-syv4qr
Expand Down Expand Up @@ -70,7 +71,7 @@ variables:
- .common-levante
variables:
# Compiler flags used for ICON model:
FCFLAGS: -g -O2 -Mrecursive -Mallocatable=03 -Mstack_arrays -Minfo=accel,inline -acc=gpu,verystrict -gpu=cc80,cuda11.7 -DRTE_USE_${FPMODEL}
FCFLAGS: -g -O2 -Mrecursive -Mallocatable=03 -Mstack_arrays -Minfo=accel,inline -acc=gpu,verystrict -gpu=cc80,cuda11.8 -DRTE_USE_${FPMODEL}
RTE_KERNELS: accel

.nag-cpu:
Expand Down
7 changes: 7 additions & 0 deletions rrtmgp-frontend/Make.depends
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -38,3 +38,10 @@ mo_gas_optics_rrtmgp.o: \
mo_gas_concentrations.o \
mo_gas_optics.o \
mo_gas_optics_rrtmgp_kernels.o mo_gas_optics_rrtmgp.F90

#
# RRTMGP cloud optics
#
mo_cloud_optics_rrtmgp.o: \
$(RTE_FORTRAN_INTERFACE) \
mo_cloud_optics_rrtmgp_kernels.o mo_cloud_optics_rrtmgp.F90
194 changes: 6 additions & 188 deletions rrtmgp-frontend/mo_cloud_optics_rrtmgp.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -28,10 +28,9 @@ module mo_cloud_optics_rrtmgp
ty_optical_props_1scl, &
ty_optical_props_2str, &
ty_optical_props_nstr
use mo_cloud_optics_rrtmgp_kernels, only: &
compute_cld_from_table, compute_cld_from_pade
implicit none
interface pade_eval
module procedure pade_eval_nbnd, pade_eval_1
end interface pade_eval
private
! -----------------------------------------------------------------------------------
type, extends(ty_optical_props), public :: ty_cloud_optics_rrtmgp
Expand Down Expand Up @@ -486,14 +485,14 @@ function cloud_optics(this, &
!
! Liquid
!
call compute_all_from_table(ncol, nlay, nbnd, liqmsk, clwp, reliq, &
call compute_cld_from_table(ncol, nlay, nbnd, liqmsk, clwp, reliq, &
this%liq_nsteps,this%liq_step_size,this%radliq_lwr, &
this%lut_extliq, this%lut_ssaliq, this%lut_asyliq, &
ltau, ltaussa, ltaussag)
!
! Ice
!
call compute_all_from_table(ncol, nlay, nbnd, icemsk, ciwp, reice, &
call compute_cld_from_table(ncol, nlay, nbnd, icemsk, ciwp, reice, &
this%ice_nsteps,this%ice_step_size,this%radice_lwr, &
this%lut_extice(:,:,this%icergh), &
this%lut_ssaice(:,:,this%icergh), &
Expand All @@ -505,13 +504,13 @@ function cloud_optics(this, &
! Hard coded assumptions: order of approximants, three size regimes
!
nsizereg = size(this%pade_extliq,2)
call compute_all_from_pade(ncol, nlay, nbnd, nsizereg, &
call compute_cld_from_pade(ncol, nlay, nbnd, nsizereg, &
liqmsk, clwp, reliq, &
2, 3, this%pade_sizreg_extliq, this%pade_extliq, &
2, 2, this%pade_sizreg_ssaliq, this%pade_ssaliq, &
2, 2, this%pade_sizreg_asyliq, this%pade_asyliq, &
ltau, ltaussa, ltaussag)
call compute_all_from_pade(ncol, nlay, nbnd, nsizereg, &
call compute_cld_from_pade(ncol, nlay, nbnd, nsizereg, &
icemsk, ciwp, reice, &
2, 3, this%pade_sizreg_extice, this%pade_extice(:,:,:,this%icergh), &
2, 2, this%pade_sizreg_ssaice, this%pade_ssaice(:,:,:,this%icergh), &
Expand Down Expand Up @@ -619,185 +618,4 @@ function get_max_radius_ice(this) result(r)

r = this%radice_upr
end function get_max_radius_ice
!--------------------------------------------------------------------------------------------------------------------
!
! Ancillary functions
!
!--------------------------------------------------------------------------------------------------------------------
!
! Linearly interpolate values from a lookup table with "nsteps" evenly-spaced
! elements starting at "offset." The table's second dimension is band.
! Returns 0 where the mask is false.
! We could also try gather/scatter for efficiency
!
subroutine compute_all_from_table(ncol, nlay, nbnd, mask, lwp, re, &
nsteps, step_size, offset, &
tau_table, ssa_table, asy_table, &
tau, taussa, taussag)
integer, intent(in) :: ncol, nlay, nbnd, nsteps
logical(wl), dimension(ncol,nlay), intent(in) :: mask
real(wp), dimension(ncol,nlay), intent(in) :: lwp, re
real(wp), intent(in) :: step_size, offset
real(wp), dimension(nsteps, nbnd), intent(in) :: tau_table, ssa_table, asy_table
real(wp), dimension(ncol,nlay,nbnd) :: tau, taussa, taussag
! ---------------------------
integer :: icol, ilay, ibnd
integer :: index
real(wp) :: fint
real(wp) :: t, ts ! tau, tau*ssa, tau*ssa*g
! ---------------------------
!$acc parallel loop gang vector default(present) collapse(3)
!$omp target teams distribute parallel do simd collapse(3)
do ibnd = 1, nbnd
do ilay = 1,nlay
do icol = 1, ncol
if(mask(icol,ilay)) then
index = min(floor((re(icol,ilay) - offset)/step_size)+1, nsteps-1)
fint = (re(icol,ilay) - offset)/step_size - (index-1)
t = lwp(icol,ilay) * &
(tau_table(index, ibnd) + fint * (tau_table(index+1,ibnd) - tau_table(index,ibnd)))
ts = t * &
(ssa_table(index, ibnd) + fint * (ssa_table(index+1,ibnd) - ssa_table(index,ibnd)))
taussag(icol,ilay,ibnd) = &
ts * &
(asy_table(index, ibnd) + fint * (asy_table(index+1,ibnd) - asy_table(index,ibnd)))
taussa (icol,ilay,ibnd) = ts
tau (icol,ilay,ibnd) = t
else
tau (icol,ilay,ibnd) = 0._wp
taussa (icol,ilay,ibnd) = 0._wp
taussag(icol,ilay,ibnd) = 0._wp
end if
end do
end do
end do
end subroutine compute_all_from_table
!
! Pade functions
!
!---------------------------------------------------------------------------
subroutine compute_all_from_pade(ncol, nlay, nbnd, nsizes, &
mask, lwp, re, &
m_ext, n_ext, re_bounds_ext, coeffs_ext, &
m_ssa, n_ssa, re_bounds_ssa, coeffs_ssa, &
m_asy, n_asy, re_bounds_asy, coeffs_asy, &
tau, taussa, taussag)
integer, intent(in) :: ncol, nlay, nbnd, nsizes
logical(wl), &
dimension(ncol,nlay), intent(in) :: mask
real(wp), dimension(ncol,nlay), intent(in) :: lwp, re
real(wp), dimension(nsizes+1), intent(in) :: re_bounds_ext, re_bounds_ssa, re_bounds_asy
integer, intent(in) :: m_ext, n_ext
real(wp), dimension(nbnd,nsizes,0:m_ext+n_ext), &
intent(in) :: coeffs_ext
integer, intent(in) :: m_ssa, n_ssa
real(wp), dimension(nbnd,nsizes,0:m_ssa+n_ssa), &
intent(in) :: coeffs_ssa
integer, intent(in) :: m_asy, n_asy
real(wp), dimension(nbnd,nsizes,0:m_asy+n_asy), &
intent(in) :: coeffs_asy
real(wp), dimension(ncol,nlay,nbnd) :: tau, taussa, taussag
! ---------------------------
integer :: icol, ilay, ibnd, irad
real(wp) :: t, ts

!$acc parallel loop gang vector default(present) collapse(3)
!$omp target teams distribute parallel do simd collapse(3)
do ibnd = 1, nbnd
do ilay = 1, nlay
do icol = 1, ncol
if(mask(icol,ilay)) then
!
! Finds index into size regime table
! This works only if there are precisely three size regimes (four bounds) and it's
! previously guaranteed that size_bounds(1) <= size <= size_bounds(4)
!
irad = min(floor((re(icol,ilay) - re_bounds_ext(2))/re_bounds_ext(3))+2, 3)
t = lwp(icol,ilay) * &
pade_eval(ibnd, nbnd, nsizes, m_ext, n_ext, irad, re(icol,ilay), coeffs_ext)

irad = min(floor((re(icol,ilay) - re_bounds_ssa(2))/re_bounds_ssa(3))+2, 3)
! Pade approximants for co-albedo can sometimes be negative
ts = t * (1._wp - max(0._wp, &
pade_eval(ibnd, nbnd, nsizes, m_ssa, n_ssa, irad, re(icol,ilay), coeffs_ssa)))

irad = min(floor((re(icol,ilay) - re_bounds_asy(2))/re_bounds_asy(3))+2, 3)
taussag(icol,ilay,ibnd) = &
ts * &
pade_eval(ibnd, nbnd, nsizes, m_asy, n_asy, irad, re(icol,ilay), coeffs_asy)

taussa (icol,ilay,ibnd) = ts
tau (icol,ilay,ibnd) = t
else
tau (icol,ilay,ibnd) = 0._wp
taussa (icol,ilay,ibnd) = 0._wp
taussag(icol,ilay,ibnd) = 0._wp
end if
end do
end do
end do

end subroutine compute_all_from_pade
!---------------------------------------------------------------------------
!
! Evaluate Pade approximant of order [m/n]
!
function pade_eval_nbnd(nbnd, nrads, m, n, irad, re, pade_coeffs)
integer, intent(in) :: nbnd, nrads, m, n, irad
real(wp), dimension(nbnd, nrads, 0:m+n), &
intent(in) :: pade_coeffs
real(wp), intent(in) :: re
real(wp), dimension(nbnd) :: pade_eval_nbnd

integer :: iband
real(wp) :: numer, denom
integer :: i

do iband = 1, nbnd
denom = pade_coeffs(iband,irad,n+m)
do i = n-1+m, 1+m, -1
denom = pade_coeffs(iband,irad,i)+re*denom
end do
denom = 1._wp +re*denom

numer = pade_coeffs(iband,irad,m)
do i = m-1, 1, -1
numer = pade_coeffs(iband,irad,i)+re*numer
end do
numer = pade_coeffs(iband,irad,0) +re*numer

pade_eval_nbnd(iband) = numer/denom
end do
end function pade_eval_nbnd
!---------------------------------------------------------------------------
!
! Evaluate Pade approximant of order [m/n]
!
function pade_eval_1(iband, nbnd, nrads, m, n, irad, re, pade_coeffs)
!$acc routine seq
!$omp declare target
!
integer, intent(in) :: iband, nbnd, nrads, m, n, irad
real(wp), dimension(nbnd, nrads, 0:m+n), &
intent(in) :: pade_coeffs
real(wp), intent(in) :: re
real(wp) :: pade_eval_1

real(wp) :: numer, denom
integer :: i

denom = pade_coeffs(iband,irad,n+m)
do i = n-1+m, 1+m, -1
denom = pade_coeffs(iband,irad,i)+re*denom
end do
denom = 1._wp +re*denom

numer = pade_coeffs(iband,irad,m)
do i = m-1, 1, -1
numer = pade_coeffs(iband,irad,i)+re*numer
end do
numer = pade_coeffs(iband,irad,0) +re*numer

pade_eval_1 = numer/denom
end function pade_eval_1
end module mo_cloud_optics_rrtmgp
6 changes: 5 additions & 1 deletion rrtmgp-kernels/Make.depends
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,6 +1,10 @@
RRTMGP_FORTRAN_KERNELS = mo_gas_optics_rrtmgp_kernels.o
RRTMGP_FORTRAN_KERNELS = mo_gas_optics_rrtmgp_kernels.o mo_cloud_optics_rrtmgp_kernels.o

#
# Gas optics
#
mo_gas_optics_rrtmgp_kernels.o: $(RTE_FORTRAN_KERNELS) mo_gas_optics_rrtmgp_kernels.F90
#
# Cloud optics
#
mo_cloud_optics_rrtmgp_kernels.o: $(RTE_FORTRAN_KERNELS) mo_cloud_optics_rrtmgp_kernels.F90
58 changes: 58 additions & 0 deletions rrtmgp-kernels/api/mo_cloud_optics_rrtmgp_kernels.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,58 @@
module mo_cloud_optics_rrtmgp_kernels
use mo_rte_kind, only : wp, wl
implicit none
private
public :: compute_cld_from_table, compute_cld_from_pade
interface
!---------------------------------------------------------------------------
!
! Linearly interpolate values from a lookup table with "nsteps" evenly-spaced
! elements starting at "offset." The table's second dimension is band.
! Returns 0 where the mask is false.
! We could also try gather/scatter for efficiency
!
subroutine compute_cld_from_table(ncol, nlay, nbnd, mask, lwp, re, &
nsteps, step_size, offset, &
tau_table, ssa_table, asy_table, &
tau, taussa, taussag) bind(C, name="rrtmgp_compute_cld_from_table")
use mo_rte_kind, only : wp, wl
integer, intent(in) :: ncol, nlay, nbnd, nsteps
logical(wl), dimension(ncol,nlay), intent(in) :: mask
real(wp), dimension(ncol,nlay), intent(in) :: lwp, re
real(wp), intent(in) :: step_size, offset
real(wp), dimension(nsteps, nbnd), intent(in) :: tau_table, ssa_table, asy_table
real(wp), dimension(ncol,nlay,nbnd) :: tau, taussa, taussag
end subroutine compute_cld_from_table

!---------------------------------------------------------------------------
!
! Pade functions
!
!---------------------------------------------------------------------------
subroutine compute_cld_from_pade(ncol, nlay, nbnd, nsizes, &
mask, lwp, re, &
m_ext, n_ext, re_bounds_ext, coeffs_ext, &
m_ssa, n_ssa, re_bounds_ssa, coeffs_ssa, &
m_asy, n_asy, re_bounds_asy, coeffs_asy, &
tau, taussa, taussag) bind(C, name="rrtmgp_compute_cld_from_pade")
use mo_rte_kind, only : wp, wl
integer, intent(in) :: ncol, nlay, nbnd, nsizes
logical(wl), &
dimension(ncol,nlay), intent(in) :: mask
real(wp), dimension(ncol,nlay), intent(in) :: lwp, re
real(wp), dimension(nsizes+1), intent(in) :: re_bounds_ext, re_bounds_ssa, re_bounds_asy
integer, intent(in) :: m_ext, n_ext
real(wp), dimension(nbnd,nsizes,0:m_ext+n_ext), &
intent(in) :: coeffs_ext
integer, intent(in) :: m_ssa, n_ssa
real(wp), dimension(nbnd,nsizes,0:m_ssa+n_ssa), &
intent(in) :: coeffs_ssa
integer, intent(in) :: m_asy, n_asy
real(wp), dimension(nbnd,nsizes,0:m_asy+n_asy), &
intent(in) :: coeffs_asy
real(wp), dimension(ncol,nlay,nbnd) :: tau, taussa, taussag
end subroutine compute_cld_from_pade
end interface
!---------------------------------------------------------------------------
end module mo_cloud_optics_rrtmgp_kernels

52 changes: 52 additions & 0 deletions rrtmgp-kernels/api/rrtmgp_kernels.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,7 @@ This header files defines the C bindings for the kernels used in RRTMGP

extern "C"
{
/* Gas optics kernels */
void rrtmgp_interpolation(
const int& ncol, const int& nlay,
const int& ngas, const int& nflav, const int& neta,
Expand Down Expand Up @@ -120,4 +121,55 @@ extern "C"
Float* lev_src, // [out] (ncol,nlay+1,ngpt)
Float* sfc_src_jac // [out] (ncol, ngpt)
);

/* Cloud optics kernels */
void rrtmgp_compute_tau_rayleigh(
const int& ncol, const int& nlay, const int& nband, const int& ngpt,
const int& ngas, const int& nflav, const int& neta, const int& npres, const int& ntemp,
const int* gpoint_flavor, // (2,ngpt)
const int* band_lims_gpt, // (2,nbnd)
const Float* krayl, // (ntemp,neta,ngpt,2)
const int& idx_h2o,
const Float* col_dry, // (ncol,nlay)
const Float* col_gas, // (ncol,nlay,ngas+1)
const Float* fminor, // (2,2,ncol,nlay,nflav)
const int* jeta, // (2, ncol,nlay,nflav)
const Bool* tropo, // (ncol,nlay)
const int* jtemp, // (ncol,nlay)
Float* tau_rayleigh // [inout] (ncol,nlay.ngpt)
);

void rrtmgp_compute_cld_from_table(
const int& ncol, int& nlay, int& nbnd, int& nsteps,
const Bool* mask, // (ncol,nlay)
const Float* lwp, // (ncol,nlay)
const Float* re, // (ncol,nlay)
const Float& step_size,
const Float& offset,
const Float* tau_table, // (nsteps, nbnd)
const Float* ssa_table, // (nsteps, nbnd)
const Float* asy_table, // (nsteps, nbnd)
Float* tau, // (ncol,nlay,nbnd)
Float* taussa, // (ncol,nlay,nbnd)
Float* taussag // (ncol,nlay,nbnd)
);

void rrtmgp_compute_cld_from_pade(
const int& ncol, int& nlay, int& nbnd, int& nsizes,
const Bool* mask, // (ncol,nlay)
const Float* lwp, // (ncol,nlay)
const Float* re, // (ncol,nlay)
const Float* re_bounds_ext, // (nsizes+1)
const Float* re_bounds_ssa, // (nsizes+1)
const Float* re_bounds_asy, // (nsizes+1)
const int& m_ext, int& n_ext,
const Float* coeffs_ext, // (nbnd,nsizes,0:m_ext+n_ext)
const int& m_ssa, int& n_ssa,
const Float* coeffs_ssa, // (nbnd,nsizes,0:m_ssa+n_ssa)
const int& m_asy, int& n_asy,
const Float* coeffs_asy, // (nbnd,nsizes,0:m_asy+n_asy)
Float* tau, // (ncol,nlay,nbnd)
Float* taussa, // (ncol,nlay,nbnd)
Float* taussag // (ncol,nlay,nbnd)
);
}
Loading
Loading