Browse Source

Moved work arrays out of fortran. Fixes #47

lon0
Bill Ladwig 7 years ago
parent
commit
71fe678e91
  1. 28
      fortran/rip_cape.f90
  2. 6
      fortran/wrf_fctt.f90
  3. 6
      fortran/wrf_vinterp.f90
  4. 19
      src/wrf/extension.py

28
fortran/rip_cape.f90

@ -272,6 +272,7 @@ END SUBROUTINE DPFCALC @@ -272,6 +272,7 @@ END SUBROUTINE DPFCALC
! NCLFORTSTART
SUBROUTINE DCAPECALC3D(prs,tmk,qvp,ght,ter,sfp,cape,cin,&
prsf, prs_new, tmk_new, qvp_new, ght_new,&
cmsg,mix,mjy,mkzh,ter_follow,&
psafile, errstat, errmsg)
USE wrf_constants, ONLY : CELKEL, G, EZERO, ESLCON1, ESLCON2, &
@ -293,6 +294,14 @@ SUBROUTINE DCAPECALC3D(prs,tmk,qvp,ght,ter,sfp,cape,cin,& @@ -293,6 +294,14 @@ SUBROUTINE DCAPECALC3D(prs,tmk,qvp,ght,ter,sfp,cape,cin,&
REAL(KIND=8), DIMENSION(mix,mjy), INTENT(IN) ::sfp
REAL(KIND=8), DIMENSION(mix,mjy,mkzh), INTENT(OUT) :: cape
REAL(KIND=8), DIMENSION(mix,mjy,mkzh), INTENT(OUT) :: cin
REAL(KIND=8), DIMENSION(mkzh,mix,mjy), INTENT(INOUT) :: prsf
REAL(KIND=8), DIMENSION(mkzh,mix,mjy), INTENT(INOUT) :: prs_new
REAL(KIND=8), DIMENSION(mkzh,mix,mjy), INTENT(INOUT) :: tmk_new
REAL(KIND=8), DIMENSION(mkzh,mix,mjy), INTENT(INOUT) :: qvp_new
REAL(KIND=8), DIMENSION(mkzh,mix,mjy), INTENT(INOUT) :: ght_new
REAL(KIND=8), INTENT(IN) :: cmsg
CHARACTER(LEN=*), INTENT(IN) :: psafile
INTEGER, INTENT(INOUT) :: errstat
@ -308,15 +317,10 @@ SUBROUTINE DCAPECALC3D(prs,tmk,qvp,ght,ter,sfp,cape,cin,& @@ -308,15 +317,10 @@ SUBROUTINE DCAPECALC3D(prs,tmk,qvp,ght,ter,sfp,cape,cin,&
REAL(KIND=8) :: eslift, tmkenv, qvpenv, tonpsadiabat
REAL(KIND=8) :: benamin, dz
REAL(KIND=8), DIMENSION(150) :: buoy, zrel, benaccum
REAL(KIND=8), DIMENSION(mkzh,mix,mjy) :: prsf
REAL(KIND=8), DIMENSION(150) :: psadithte, psadiprs
REAL(KIND=8), DIMENSION(150,150) :: psaditmk
LOGICAL :: elfound
REAL(KIND=8), DIMENSION(mkzh,mix,mjy) :: prs_new
REAL(KIND=8), DIMENSION(mkzh,mix,mjy) :: tmk_new
REAL(KIND=8), DIMENSION(mkzh,mix,mjy) :: qvp_new
REAL(KIND=8), DIMENSION(mkzh,mix,mjy) :: ght_new
! To remove compiler warnings
tmkpari = 0
@ -597,6 +601,7 @@ END SUBROUTINE DCAPECALC3D @@ -597,6 +601,7 @@ END SUBROUTINE DCAPECALC3D
! NCLFORTSTART
SUBROUTINE DCAPECALC2D(prs,tmk,qvp,ght,ter,sfp,cape,cin,&
prsf, prs_new, tmk_new, qvp_new, ght_new,&
cmsg,mix,mjy,mkzh,ter_follow,&
psafile, errstat, errmsg)
USE wrf_constants, ONLY : CELKEL, G, EZERO, ESLCON1, ESLCON2, &
@ -618,6 +623,13 @@ SUBROUTINE DCAPECALC2D(prs,tmk,qvp,ght,ter,sfp,cape,cin,& @@ -618,6 +623,13 @@ SUBROUTINE DCAPECALC2D(prs,tmk,qvp,ght,ter,sfp,cape,cin,&
REAL(KIND=8), DIMENSION(mix,mjy), INTENT(IN) ::sfp
REAL(KIND=8), DIMENSION(mix,mjy,mkzh), INTENT(OUT) :: cape
REAL(KIND=8), DIMENSION(mix,mjy,mkzh), INTENT(OUT) :: cin
REAL(KIND=8), DIMENSION(mkzh,mix,mjy), INTENT(INOUT) :: prsf
REAL(KIND=8), DIMENSION(mkzh,mix,mjy), INTENT(INOUT) :: prs_new
REAL(KIND=8), DIMENSION(mkzh,mix,mjy), INTENT(INOUT) :: tmk_new
REAL(KIND=8), DIMENSION(mkzh,mix,mjy), INTENT(INOUT) :: qvp_new
REAL(KIND=8), DIMENSION(mkzh,mix,mjy), INTENT(INOUT) :: ght_new
REAL(KIND=8), INTENT(IN) :: cmsg
CHARACTER(LEN=*), INTENT(IN) :: psafile
INTEGER, INTENT(INOUT) :: errstat
@ -635,15 +647,11 @@ SUBROUTINE DCAPECALC2D(prs,tmk,qvp,ght,ter,sfp,cape,cin,& @@ -635,15 +647,11 @@ SUBROUTINE DCAPECALC2D(prs,tmk,qvp,ght,ter,sfp,cape,cin,&
REAL(KIND=8) :: eslift, tmkenv, qvpenv, tonpsadiabat
REAL(KIND=8) :: benamin, dz, pup, pdn
REAL(KIND=8), DIMENSION(150) :: buoy, zrel, benaccum
REAL(KIND=8), DIMENSION(mkzh,mix,mjy) :: prsf
REAL(KIND=8), DIMENSION(150) :: psadithte, psadiprs
REAL(KIND=8), DIMENSION(150,150) :: psaditmk
LOGICAL :: elfound
REAL(KIND=8), DIMENSION(mkzh) :: eth_temp
REAL(KIND=8), DIMENSION(mkzh,mix,mjy) :: prs_new
REAL(KIND=8), DIMENSION(mkzh,mix,mjy) :: tmk_new
REAL(KIND=8), DIMENSION(mkzh,mix,mjy) :: qvp_new
REAL(KIND=8), DIMENSION(mkzh,mix,mjy) :: ght_new
! To remove compiler warnings
errstat = 0

6
fortran/wrf_fctt.f90

@ -1,5 +1,5 @@ @@ -1,5 +1,5 @@
!NCLFORTSTART
SUBROUTINE wrfcttcalc(prs, tk, qci, qcw, qvp, ght, ter, ctt, haveqci, nz, ns, ew)
SUBROUTINE wrfcttcalc(prs, tk, qci, qcw, qvp, ght, ter, ctt, pf, haveqci, nz, ns, ew)
USE wrf_constants, ONLY : EPS, USSALR, RD, G, ABSCOEFI, ABSCOEF, CELKEL
IMPLICIT NONE
@ -12,6 +12,8 @@ SUBROUTINE wrfcttcalc(prs, tk, qci, qcw, qvp, ght, ter, ctt, haveqci, nz, ns, ew @@ -12,6 +12,8 @@ SUBROUTINE wrfcttcalc(prs, tk, qci, qcw, qvp, ght, ter, ctt, haveqci, nz, ns, ew
REAL(KIND=8), DIMENSION(ew,ns), INTENT(IN) :: ter
REAL(KIND=8), DIMENSION(ew,ns), INTENT(OUT) :: ctt
REAL(KIND=8), DIMENSION(ew,ns,nz), INTENT(INOUT) :: pf
!NCLEND
! REAL(KIND=8) :: znfac(nz)
@ -20,7 +22,7 @@ SUBROUTINE wrfcttcalc(prs, tk, qci, qcw, qvp, ght, ter, ctt, haveqci, nz, ns, ew @@ -20,7 +22,7 @@ SUBROUTINE wrfcttcalc(prs, tk, qci, qcw, qvp, ght, ter, ctt, haveqci, nz, ns, ew
INTEGER i,j,k,ripk
REAL(KIND=8) :: opdepthu, opdepthd, dp, arg1, fac, prsctt, ratmix
REAL(KIND=8) :: arg2, agl_hgt, vt
REAL(KIND=8), DIMENSION(ew,ns,nz) :: pf
REAL(KIND=8) :: p1, p2
!$OMP PARALLEL

6
fortran/wrf_vinterp.f90

@ -127,7 +127,7 @@ END FUNCTION wrf_intrp_value @@ -127,7 +127,7 @@ END FUNCTION wrf_intrp_value
!NCLFORTSTART
SUBROUTINE wrf_vintrp(datain, dataout, pres, tk, qvp, ght, terrain,&
sfp, smsfp, vcarray, interp_levels, numlevels,&
icase, ew, ns, nz, extrap, vcor, logp, rmsg,&
icase, ew, ns, nz, extrap, vcor, logp, tempout, rmsg,&
errstat, errmsg)
USE wrf_constants, ONLY : ALGERR, SCLHT, EXPON, EXPONI, GAMMA, GAMMAMD, TLCLC1, &
TLCLC2, TLCLC3, TLCLC4, THTECON1, THTECON2, THTECON3, &
@ -146,6 +146,9 @@ SUBROUTINE wrf_vintrp(datain, dataout, pres, tk, qvp, ght, terrain,& @@ -146,6 +146,9 @@ SUBROUTINE wrf_vintrp(datain, dataout, pres, tk, qvp, ght, terrain,&
REAL(KIND=8), DIMENSION(ew,ns,numlevels), INTENT(OUT) :: dataout
REAL(KIND=8), DIMENSION(ew,ns,nz), INTENT(IN) :: vcarray
REAL(KIND=8), DIMENSION(numlevels), INTENT(IN) :: interp_levels
REAL(KIND=8), DIMENSION(ew,ns), INTENT(INOUT) :: tempout
REAL(KIND=8), INTENT(IN) :: rmsg
INTEGER, INTENT(INOUT) :: errstat
CHARACTER(LEN=*), INTENT(INOUT) :: errmsg
@ -156,7 +159,6 @@ SUBROUTINE wrf_vintrp(datain, dataout, pres, tk, qvp, ght, terrain,& @@ -156,7 +159,6 @@ SUBROUTINE wrf_vintrp(datain, dataout, pres, tk, qvp, ght, terrain,&
INTEGER :: i, j, k, kupper !itriv
INTEGER :: ifound, isign !miy,mjx
INTEGER :: log_errcnt, interp_errcnt, interp_errstat
REAL(KIND=8), DIMENSION(ew,ns) :: tempout
REAL(KIND=8) :: rlevel, vlev, diff
REAL(KIND=8) :: tmpvlev
REAL(KIND=8) :: vcp1, vcp0, valp0, valp1

19
src/wrf/extension.py

@ -618,6 +618,14 @@ def _cape(p_hpa, tk, qv, ht, ter, sfp, missing, i3dflag, ter_follow, @@ -618,6 +618,14 @@ def _cape(p_hpa, tk, qv, ht, ter, sfp, missing, i3dflag, ter_follow,
else:
cape_routine = dcapecalc2d
# Work arrays
k_left_shape = (p_hpa.shape[2], p_hpa.shape[0], p_hpa.shape[1])
prsf = np.empty(k_left_shape, np.float64, order="F")
prs_new = np.empty(k_left_shape, np.float64, order="F")
tmk_new = np.empty(k_left_shape, np.float64, order="F")
qvp_new = np.empty(k_left_shape, np.float64, order="F")
ght_new = np.empty(k_left_shape, np.float64, order="F")
# note that p_hpa, tk, qv, and ht have the vertical flipped
result = cape_routine(p_hpa,
tk,
@ -627,6 +635,11 @@ def _cape(p_hpa, tk, qv, ht, ter, sfp, missing, i3dflag, ter_follow, @@ -627,6 +635,11 @@ def _cape(p_hpa, tk, qv, ht, ter, sfp, missing, i3dflag, ter_follow,
sfp,
capeview,
cinview,
prsf,
prs_new,
tmk_new,
qvp_new,
ght_new,
missing,
ter_follow,
psafile,
@ -766,6 +779,8 @@ def _ctt(p_hpa, tk, qice, qcld, qv, ght, ter, haveqci, outview=None): @@ -766,6 +779,8 @@ def _ctt(p_hpa, tk, qice, qcld, qv, ght, ter, haveqci, outview=None):
if outview is None:
outview = np.empty_like(ter)
pf = np.empty(p_hpa.shape[0:3], np.float64, order="F")
result = wrfcttcalc(p_hpa,
tk,
qice,
@ -774,6 +789,7 @@ def _ctt(p_hpa, tk, qice, qcld, qv, ght, ter, haveqci, outview=None): @@ -774,6 +789,7 @@ def _ctt(p_hpa, tk, qice, qcld, qv, ght, ter, haveqci, outview=None):
ght,
ter,
outview,
pf,
haveqci)
return result
@ -859,6 +875,8 @@ def _vintrp(field, pres, tk, qvp, ght, terrain, sfp, smsfp, @@ -859,6 +875,8 @@ def _vintrp(field, pres, tk, qvp, ght, terrain, sfp, smsfp,
outdims = field.shape[0:2] + interp_levels.shape
outview = np.empty(outdims, field.dtype, order="F")
tempout = np.zeros(field.shape[0:2], np.float64, order="F")
errstat = np.array(0)
errmsg = np.zeros(Constants.ERRLEN, "c")
@ -877,6 +895,7 @@ def _vintrp(field, pres, tk, qvp, ght, terrain, sfp, smsfp, @@ -877,6 +895,7 @@ def _vintrp(field, pres, tk, qvp, ght, terrain, sfp, smsfp,
extrap,
vcor,
logp,
tempout,
missing,
errstat,
errmsg)

Loading…
Cancel
Save