From ca50a6b9f77cb964e885f336bc4b7f283fc15f70 Mon Sep 17 00:00:00 2001 From: Peter Norris Date: Tue, 28 Dec 2021 14:11:32 -0500 Subject: [PATCH 1/9] pmn: add TAUTTX efficiently --- .../GEOSsolar_GridComp/GEOS_SolarGridComp.F90 | 194 +++++++++++------- 1 file changed, 122 insertions(+), 72 deletions(-) diff --git a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 index 5e7c129bb..f0a8a1274 100644 --- a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 +++ b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 @@ -988,12 +988,18 @@ subroutine SetServices ( GC, RC ) VLOCATION = MAPL_VLocationNone, __RC__) call MAPL_AddExportSpec(GC, & - LONG_NAME = 'in_cloud_optical_thickness_of_all_clouds', & + LONG_NAME = 'in_cloud_optical_thickness_of_all_clouds__deprecated', & UNITS = '1' , & SHORT_NAME = 'TAUTT', & DIMS = MAPL_DimsHorzOnly, & VLOCATION = MAPL_VLocationNone, __RC__) -!? PMN ... implement the approx fix as another variable ... this one is broken + + call MAPL_AddExportSpec(GC, & + LONG_NAME = 'in_cloud_optical_thickness_of_all_clouds__improved', & + UNITS = '1' , & + SHORT_NAME = 'TAUTTX', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__) call MAPL_AddExportSpec(GC, & LONG_NAME = 'in_cloud_optical_thickness_for_ice_clouds', & @@ -4314,13 +4320,18 @@ subroutine UPDATE_EXPORT(IM,JM,LM, RC) real, pointer, dimension(:,:,:,:) :: TAUCLD, HYDROMETS, REFF real, pointer, dimension(:,:,:) :: TAUI,TAUW,TAUR,TAUS + ! for efficiency + real, allocatable, dimension(:,:) :: aCLDL,aCLDM,aCLDH + real, allocatable, dimension(:,:) :: aTAUL,aTAUM,aTAUH + real, allocatable, dimension(:,:) :: aCLDT + real, dimension(LM ) :: DUM1D real, dimension(LM,4) :: DUM2D real, pointer, dimension(:,:) :: TDUST,TSALT,TSO4,TBC,TOC - real, pointer, dimension(:,:) :: CLDH,CLDM,CLDL,CLDT, & - TAUH,TAUM,TAUL,TAUT, & - CLDTMP, CLDPRS + real, pointer, dimension(:,:) :: CLDH,CLDM,CLDL,CLDT, & + TAUH,TAUM,TAUL,TAUT,TAUTX, & + CLDTMP,CLDPRS type (ESMF_FieldBundle) :: BUNDLE type (ESMF_Field) :: FIELD @@ -4447,55 +4458,49 @@ subroutine UPDATE_EXPORT(IM,JM,LM, RC) call MAPL_GetPointer(EXPORT , TAUM, 'TAUMD', __RC__) call MAPL_GetPointer(EXPORT , TAUH, 'TAUHI', __RC__) call MAPL_GetPointer(EXPORT , TAUT, 'TAUTT', __RC__) + call MAPL_GetPointer(EXPORT , TAUTX, 'TAUTTX', __RC__) call MAPL_GetPointer(EXPORT , CLDTMP, 'CLDTMP', __RC__) call MAPL_GetPointer(EXPORT , CLDPRS, 'CLDPRS', __RC__) - if(associated(FCLD)) FCLD = CLIN + if (associated(FCLD)) FCLD = CLIN - if(associated(CLDH)) then - CLDH = 0. + if (associated(CLDH) .or. associated(CLDT) .or. associated(TAUTX)) then + allocate(aCLDH(IM,JM),__STAT__) + aCLDH = 0. do l=1,LCLDMH-1 - CLDH = max(CLDH,CLIN(:,:,L)) + aCLDH = max(aCLDH,CLIN(:,:,L)) end do + if (associated(CLDH)) CLDH = aCLDH end if - if(associated(CLDM)) then - CLDM = 0. + if (associated(CLDM) .or. associated(CLDT) .or. associated(TAUTX)) then + allocate(aCLDM(IM,JM),__STAT__) + aCLDM = 0. do l=LCLDMH,LCLDLM-1 - CLDM = max(CLDM,CLIN(:,:,L)) + aCLDM = max(aCLDM,CLIN(:,:,L)) end do + if (associated(CLDM)) CLDM = aCLDM end if - if(associated(CLDL)) then - CLDL = 0. + if (associated(CLDL) .or. associated(CLDT) .or. associated(TAUTX)) then + allocate(aCLDL(IM,JM),__STAT__) + aCLDL = 0. do l=LCLDLM,LM - CLDL = max(CLDL,CLIN(:,:,L)) + aCLDL = max(aCLDL,CLIN(:,:,L)) end do + if (associated(CLDL)) CLDL = aCLDL end if - if(associated(CLDT)) then - CLD = 0. - do l=1,LCLDMH-1 - CLD = max(CLD,CLIN(:,:,L)) - end do - CLDT = (1-CLD) - CLD = 0. - do l= LCLDMH,LCLDLM-1 - CLD = max(CLD,CLIN(:,:,L)) - end do - CLDT = CLDT*(1-CLD) - CLD = 0. - do l=LCLDLM,LM - CLD = max(CLD,CLIN(:,:,L)) - end do - CLDT = 1.0 - CLDT*(1-CLD) + if (associated(CLDT) .or. associated(TAUTX)) then + allocate(aCLDT(IM,JM),__STAT__) + aCLDT = 1. - (1-aCLDH)*(1-aCLDM)*(1-aCLDL) + if (associated(CLDT)) CLDT = aCLDT end if - if(associated(TAUI ).or.associated(TAUW ).or. & - associated(TAUR ).or.associated(TAUS ).or. & - associated(TAUL ).or.associated(TAUM ).or. & - associated(TAUH ).or.associated(TAUT ).or. & - associated(CLDTMP).or.associated(CLDPRS)) then + if (associated(TAUI) .or. associated(TAUW) .or. associated(TAUR) .or. associated(TAUS).or. & + associated(TAUL) .or. associated(TAUM) .or. associated(TAUH) .or. & + associated(TAUT) .or. associated(TAUTX) .or. & + associated(CLDTMP) .or. associated(CLDPRS)) then allocate( TAUCLD(IM,JM,LM,4), __STAT__) allocate(HYDROMETS(IM,JM,LM,4), __STAT__) @@ -4510,17 +4515,17 @@ subroutine UPDATE_EXPORT(IM,JM,LM, RC) ! 3 Falling Liquid (Rain) ! 4 Falling Ice (Rain) - REFF(:,:,:,1) = RRI * 1.0e6 ! REFF must be in microns - REFF(:,:,:,2) = RRL * 1.0e6 - REFF(:,:,:,3) = RRR * 1.0e6 - REFF(:,:,:,4) = RRS * 1.0e6 + REFF(:,:,:,1) = RRI * 1.e6 ! REFF must be in microns + REFF(:,:,:,2) = RRL * 1.e6 + REFF(:,:,:,3) = RRR * 1.e6 + REFF(:,:,:,4) = RRS * 1.e6 HYDROMETS(:,:,:,1) = RQI HYDROMETS(:,:,:,2) = RQL HYDROMETS(:,:,:,3) = RQR HYDROMETS(:,:,:,4) = RQS - TAUCLD = 0.0 + TAUCLD = 0. ! Due to the generic use of this routine, it currently works on one column at a time, ! thus the need for the array sections below. @@ -4528,60 +4533,100 @@ subroutine UPDATE_EXPORT(IM,JM,LM, RC) ! NOTE: Dummy arrays are passed into outputs 1 and 3 because these are currently only ! used in sorad.F90. - DO I = 1, IM - DO J = 1, JM - CALL GETVISTAU(LM,ZTH(I,J),DP(I,J,:),CLIN(I,J,:),REFF(I,J,:,:),HYDROMETS(I,J,:,:),LCLDMH,LCLDLM,& - DUM2D(:,:),TAUCLD(I,J,:,:),DUM1D(:)) + DO I = 1,IM + DO J = 1,JM + CALL GETVISTAU( & + LM,ZTH(I,J),DP(I,J,:),& + CLIN(I,J,:),REFF(I,J,:,:),HYDROMETS(I,J,:,:),& + LCLDMH,LCLDLM,& + DUM2D(:,:),TAUCLD(I,J,:,:),DUM1D(:)) END DO END DO - if(associated(TAUI)) TAUI = TAUCLD(:,:,:,1) - if(associated(TAUW)) TAUW = TAUCLD(:,:,:,2) - if(associated(TAUR)) TAUR = TAUCLD(:,:,:,3) - if(associated(TAUS)) TAUS = TAUCLD(:,:,:,4) + if (associated(TAUI)) TAUI = TAUCLD(:,:,:,1) + if (associated(TAUW)) TAUW = TAUCLD(:,:,:,2) + if (associated(TAUR)) TAUR = TAUCLD(:,:,:,3) + if (associated(TAUS)) TAUS = TAUCLD(:,:,:,4) + ! use the total hydrometor optical thickness for the general opticl thicknesses below TAUCLD(:,:,:,1) = TAUCLD(:,:,:,1) + TAUCLD(:,:,:,2) + TAUCLD(:,:,:,3) + TAUCLD(:,:,:,4) - if(associated(TAUH)) then - TAUH = 0. + ! TAU[HML] are correct because GETVISTAU produces in-cloud optical thicknesses for + ! 'effective clouds' extended-out and diluted to the maximum cloud fraction in each + ! pressure super-band [LMT]. + + if (associated(TAUH) .or. associated(TAUT) .or. associated(TAUTX)) then + allocate(aTAUH(IM,JM),__STAT__) + aTAUH = 0. do l=1,LCLDMH-1 - TAUH = TAUH + TAUCLD(:,:,L,1) + aTAUH = aTAUH + TAUCLD(:,:,L,1) end do + if (associated(TAUH)) TAUH = aTAUH end if - if(associated(TAUM)) then - TAUM = 0. + if (associated(TAUM) .or. associated(TAUT) .or. associated(TAUTX)) then + allocate(aTAUM(IM,JM),__STAT__) + aTAUM = 0. do l=LCLDMH,LCLDLM-1 - TAUM = TAUM + TAUCLD(:,:,L,1) + aTAUM = aTAUM + TAUCLD(:,:,L,1) end do + if (associated(TAUM)) TAUM = aTAUM end if - if(associated(TAUL)) then - TAUL = 0. + if (associated(TAUL) .or. associated(TAUT) .or. associated(TAUTX)) then + allocate(aTAUL(IM,JM),__STAT__) + aTAUL = 0. do l=LCLDLM,LM - TAUL = TAUL + TAUCLD(:,:,L,1) + aTAUL = aTAUL + TAUCLD(:,:,L,1) end do + if (associated(TAUL)) TAUL = aTAUL end if - if(associated(TAUT)) then - TAUT = 0. - do l=1,LM - TAUT = TAUT + TAUCLD(:,:,L,1) - end do + ! TAUT however is broken because the three super-bands are randomly overlapped + ! and with different effective cloud fractions. It has been broken but used for + ! a long time. It should be considered deprecated. TAUTX below is an improved + ! version. + + if (associated(TAUT)) TAUT = aTAUH + aTAUM + aTAUL + + ! As noted above, one cannot simply add TAUL, TAUM and TAUH to get a column + ! in-cloud optical thickness, because the actual column value depends on the + ! overlap of these bands. This overlap is here assumed random. We can express + ! the approximate column in-cloud optical thickness in terms of the sum over + ! the 2**3 - 1 combinations with some cloud in at least one of the 3 bands, + ! each with their respective fractions. For random overlap, CLDL*CLDM*CLDH of + ! the gridcolumn would have a column TAU of TAUL+TAUM+TAUH, CLDL*CLDM*(1-CLDH) + ! would have a column TAU of TAUL+TAUM, etc. Then, for an in-cloud column TAU, + ! the sum of the 7 must be normalized by the random column cloud fraction + ! CLDT = 1 – (1-CLDL)*(1-CLDM)*(1-CLDH). + ! Not surprisingly this gives + ! TAUTX = (TAUL*CLDL + TAUM*CLDM + TAUH*CLDH) / CLDT, + ! because we assume we can linearly average optical thickness among the comb- + ! inations. This assumption is questionable, since cloud radiative properties + ! are non-linear in optical thickness. This is why TAUTX is approximate. But + ! its the best we SIMPLY can do. + + if (associated(TAUTX)) then + TAUTX = 0. + where (aCLDT > 0.) TAUTX = (aTAUL*aCLDL + aTAUM*aCLDM + aTAUH*aCLDH) / CLDT end if - if(associated(CLDTMP).or.associated(CLDPRS)) then - call MAPL_GetResource(MAPL, TAUCRIT , 'TAUCRIT:', DEFAULT=0.10, __RC__) + if (allocated(aTAUH)) deallocate(aTAUH,__STAT__) + if (allocated(aTAUM)) deallocate(aTAUM,__STAT__) + if (allocated(aTAUL)) deallocate(aTAUL,__STAT__) + + if (associated(CLDTMP) .or. associated(CLDPRS)) then + call MAPL_GetResource(MAPL,TAUCRIT,'TAUCRIT:',DEFAULT=0.10,__RC__) - if(associated(CLDTMP)) CLDTMP = MAPL_UNDEF - if(associated(CLDPRS)) CLDPRS = MAPL_UNDEF + if (associated(CLDTMP)) CLDTMP = MAPL_UNDEF + if (associated(CLDPRS)) CLDPRS = MAPL_UNDEF - do l=LM,1,-1 - if(associated(CLDTMP)) then - where(TAUCLD(:,:,L,1)>TAUCRIT) CLDTMP = T(:,:,L) + do L=LM,1,-1 + if (associated(CLDTMP)) then + where (TAUCLD(:,:,L,1) > TAUCRIT) CLDTMP = T(:,:,L) end if - if(associated(CLDPRS)) then - where(TAUCLD(:,:,L,1)>TAUCRIT) CLDPRS = PLL(:,:,L-1) + if (associated(CLDPRS)) then + where (TAUCLD(:,:,L,1) > TAUCRIT) CLDPRS = PLL(:,:,L-1) end if end do end if @@ -4593,6 +4638,11 @@ subroutine UPDATE_EXPORT(IM,JM,LM, RC) end if + if (allocated(aCLDH)) deallocate(aCLDH,__STAT__) + if (allocated(aCLDM)) deallocate(aCLDM,__STAT__) + if (allocated(aCLDL)) deallocate(aCLDL,__STAT__) + if (allocated(aCLDT)) deallocate(aCLDT,__STAT__) + ! Fill Albedos !------------- From 91d01a82c975b02c96bf53bc3ce5ac9f5c85e236 Mon Sep 17 00:00:00 2001 From: Peter Norris Date: Sun, 9 Jan 2022 15:20:42 -0500 Subject: [PATCH 2/9] pmn: add TAUxxPAR diagnostics plus cosmetics --- .../GEOSirrad_GridComp/GEOS_IrradGridComp.F90 | 2 +- .../rrtmg_lw/gcm_model/src/rrtmg_lw_rad.F90 | 8 +- .../GEOSsolar_GridComp/GEOS_SolarGridComp.F90 | 113 ++++++++++++++--- .../gcm_model/src/rrtmg_sw_cldprmc.F90 | 18 ++- .../rrtmg_sw/gcm_model/src/rrtmg_sw_rad.F90 | 43 +++++-- .../gcm_model/src/rrtmg_sw_spcvmc.F90 | 119 +++++++++++++++++- 6 files changed, 265 insertions(+), 38 deletions(-) diff --git a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSirrad_GridComp/GEOS_IrradGridComp.F90 b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSirrad_GridComp/GEOS_IrradGridComp.F90 index 48b52cc49..f201db8a2 100644 --- a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSirrad_GridComp/GEOS_IrradGridComp.F90 +++ b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSirrad_GridComp/GEOS_IrradGridComp.F90 @@ -3345,7 +3345,7 @@ subroutine Lw_Driver(IM,JM,LM,LATS,LONS,CoresPerNode,RC) do I = 1,IM IJ = IJ + 1 - ! convert super-band clearCounts to cloud fractions + ! convert super-layer clearCounts to cloud fractions if(associated(CLDTTLW)) then CLDTTLW(I,J) = 1.0 - CLEARCOUNTS(IJ,1)/float(NGPTLW) endif diff --git a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSirrad_GridComp/RRTMG/rrtmg_lw/gcm_model/src/rrtmg_lw_rad.F90 b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSirrad_GridComp/RRTMG/rrtmg_lw/gcm_model/src/rrtmg_lw_rad.F90 index 90b6816df..26d208feb 100644 --- a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSirrad_GridComp/RRTMG/rrtmg_lw/gcm_model/src/rrtmg_lw_rad.F90 +++ b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSirrad_GridComp/RRTMG/rrtmg_lw/gcm_model/src/rrtmg_lw_rad.F90 @@ -181,7 +181,7 @@ subroutine rrtmg_lw( & ! ----- Output ----- - ! subcolumn clear counts for Tot|High|Mid|Low bands + ! subcolumn clear counts for Tot|High|Mid|Low super-layers integer, intent(out) :: clearCounts(ncol,4) real, intent(out) :: uflx (ncol,nlay+1) ! Total sky LW upward flux [W/m2] @@ -426,7 +426,7 @@ subroutine rrtmg_lw_part( & ! ----- Output ----- - ! subcolumn clear counts for Tot|High|Mid|Low bands + ! subcolumn clear counts for Tot|High|Mid|Low super-layers integer, intent(out) :: clearCounts(ncol,4) real, intent(out) :: uflx (ncol,nlay+1) ! Total sky LW upward flux [W/m2] @@ -484,7 +484,7 @@ subroutine rrtmg_lw_part( & real :: ciwpmc (nlay,ngptlw,pncol) ! cloud ice water path [g/m2] real :: clwpmc (nlay,ngptlw,pncol) ! cloud liq water path [g/m2] real :: taucmc (nlay,ngptlw,pncol) ! cloud optical depth - integer :: p_clearCounts (4,pncol) ! for super-band cld fractions + integer :: p_clearCounts (4,pncol) ! for super-layer cld fractions ! cloudy for ANY subcol/gpoint of column? logical :: cloudy (nlay,pncol) @@ -554,7 +554,7 @@ subroutine rrtmg_lw_part( & cldymc, ciwpmc, clwpmc, & seed_order=[1,2,3,4]) - ! for super-band cloud fractions + ! for super-layer cloud fractions call clearCounts_threeBand( & pncol, pncol, ngptlw, nlay, cloudLM, cloudMH, cldymc, & diff --git a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 index f0a8a1274..2d6e6607b 100644 --- a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 +++ b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 @@ -1029,6 +1029,34 @@ subroutine SetServices ( GC, RC ) DIMS = MAPL_DimsHorzVert, & VLOCATION = MAPL_VLocationCenter, __RC__) + call MAPL_AddExportSpec(GC, & + LONG_NAME = 'in_cloud_optical_thickness_of_low_clouds_RRTMG_PAR', & + UNITS = '1' , & + SHORT_NAME = 'TAULOPAR', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__) + + call MAPL_AddExportSpec(GC, & + LONG_NAME = 'in_cloud_optical_thickness_of_middle_clouds_RRTMG_PAR', & + UNITS = '1' , & + SHORT_NAME = 'TAUMDPAR', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__) + + call MAPL_AddExportSpec(GC, & + LONG_NAME = 'in_cloud_optical_thickness_of_high_clouds_RRTMG_PAR', & + UNITS = '1' , & + SHORT_NAME = 'TAUHIPAR', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__) + + call MAPL_AddExportSpec(GC, & + LONG_NAME = 'in_cloud_optical_thickness_of_all_clouds_RRTMG_PAR', & + UNITS = '1' , & + SHORT_NAME = 'TAUTTPAR', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__) + call MAPL_AddExportSpec(GC, & LONG_NAME = 'surface_net_downward_shortwave_flux_assuming_clear_sky',& UNITS = 'W m-2', & @@ -1982,17 +2010,24 @@ subroutine SORADCORE(IM,JM,LM,NO_AERO,CURRTIME,LoadBalance,RC) real, dimension(IM,JM) :: ZTH, SLR logical, dimension(IM,JM) :: daytime - ! super-band RRTMG cloud fraction exports + ! super-layer RRTMG cloud fraction exports real, pointer, dimension(:,:) :: CLDTTSW real, pointer, dimension(:,:) :: CLDHISW real, pointer, dimension(:,:) :: CLDMDSW real, pointer, dimension(:,:) :: CLDLOSW + ! super-layer RRTMG PAR optical thickness exports + real, pointer, dimension(:,:) :: TAUTTPAR + real, pointer, dimension(:,:) :: TAUHIPAR + real, pointer, dimension(:,:) :: TAUMDPAR + real, pointer, dimension(:,:) :: TAULOPAR + ! DAYTIME ONLY COPY OF VARIABLES real, pointer, dimension(: ) :: ALBNR,ALBNF,ALBVR,ALBVF,ZT,SLR1D, & - UVRR,UVRF,PARR,PARF,NIRR,NIRF, & - Ig1D,Jg1D,CLDTS,CLDHS,CLDMS,CLDLS + UVRR,UVRF,PARR,PARF,NIRR,NIRF,Ig1D,Jg1D, & + CLDTS,CLDHS,CLDMS,CLDLS, & + TAUTP,TAUHP,TAUMP,TAULP real, pointer, dimension(:,:,:) :: FCLD,TAUI,TAUW,CLIN,RRL,RRI,RQI,RQL,RQR real, pointer, dimension(:,:,:) :: DP, PLL real, pointer, dimension(:,:,:) :: RAERO @@ -2197,12 +2232,18 @@ subroutine SORADCORE(IM,JM,LM,NO_AERO,CURRTIME,LoadBalance,RC) end do end do - ! super-band RRTMG cloud fraction exports + ! super-layer RRTMG cloud fraction exports call MAPL_GetPointer(EXPORT, CLDTTSW, 'CLDTTSW', __RC__) call MAPL_GetPointer(EXPORT, CLDHISW, 'CLDHISW', __RC__) call MAPL_GetPointer(EXPORT, CLDMDSW, 'CLDMDSW', __RC__) call MAPL_GetPointer(EXPORT, CLDLOSW, 'CLDLOSW', __RC__) + ! super-layer RRTMG PAR optical thickness exports + call MAPL_GetPointer(EXPORT, TAUTTPAR, 'TAUTTPAR', __RC__) + call MAPL_GetPointer(EXPORT, TAUHIPAR, 'TAUHIPAR', __RC__) + call MAPL_GetPointer(EXPORT, TAUMDPAR, 'TAUMDPAR', __RC__) + call MAPL_GetPointer(EXPORT, TAULOPAR, 'TAULOPAR', __RC__) + call MAPL_TimerOff(MAPL,"-MISC") ! Load balancing by packing the lit points and sharing work with night regions @@ -2257,11 +2298,12 @@ subroutine SORADCORE(IM,JM,LM,NO_AERO,CURRTIME,LoadBalance,RC) ! component needs the LATS, SLR and ZTH from MAPL and the global ! gridcolumn indicies Ig and Jg. The outputs are the INTERNAL ! variables being refreshed plus four cloud fraction diagnostics -! (CLDTTSW, CLDHISW, CLDMDSW, CLDLOSW). +! (CLDTTSW, CLDHISW, CLDMDSW, CLDLOSW) & four optical thickness +! diagnostics (TAUTTPAR, TAUHIPAR, TAUMDPAR, TAULOPAR). !-------------------------------------------------------------- NumInp = size(ImportSpec) + 5 - NumOut = size(InternalSpec) + 4 + NumOut = size(InternalSpec) + 8 allocate(SlicesInp(NumInp), NamesInp(NumInp), & SlicesOut(NumOut), NamesOut(NumOut), __STAT__) @@ -2528,21 +2570,29 @@ subroutine SORADCORE(IM,JM,LM,NO_AERO,CURRTIME,LoadBalance,RC) OUTPUT_VARS_1: do k=1,NumOut - if (k < NumOut-3) then + if (k < NumOut-7) then ! internal outputs call MAPL_VarSpecGet(InternalSpec(k), & DIMS=dims, SHORT_NAME=NamesOut(k), __RC__) else ! cloud fraction outputs dims = MAPL_DIMSHORZONLY - if (k == NumOut-3) then + if (k == NumOut-7) then NamesOut(k) = "CLDTTSW" - else if (k == NumOut-2) then + else if (k == NumOut-6) then NamesOut(k) = "CLDHISW" - else if (k == NumOut-1) then + else if (k == NumOut-5) then NamesOut(k) = "CLDMDSW" - else + else if (k == NumOut-4) then NamesOut(k) = "CLDLOSW" + else if (k == NumOut-3) then + NamesOut(k) = "TAUTTPAR" + else if (k == NumOut-2) then + NamesOut(k) = "TAUHIPAR" + else if (k == NumOut-1) then + NamesOut(k) = "TAUMDPAR" + else + NamesOut(k) = "TAULOPAR" end if end if @@ -2640,6 +2690,14 @@ subroutine SORADCORE(IM,JM,LM,NO_AERO,CURRTIME,LoadBalance,RC) CLDMS => ptr2(1:Num2do,1) case('CLDLOSW') CLDLS => ptr2(1:Num2do,1) + case('TAUTTPAR') + TAUTP => ptr2(1:Num2do,1) + case('TAUHIPAR') + TAUHP => ptr2(1:Num2do,1) + case('TAUMDPAR') + TAUMP => ptr2(1:Num2do,1) + case('TAULOPAR') + TAULP => ptr2(1:Num2do,1) end select enddo OUTPUT_VARS_2 @@ -3444,7 +3502,7 @@ subroutine SORADCORE(IM,JM,LM,NO_AERO,CURRTIME,LoadBalance,RC) allocate(O3_R (size(Q,1),size(Q,2)),__STAT__) allocate(CO2_R (size(Q,1),size(Q,2)),__STAT__) allocate(CH4_R (size(Q,1),size(Q,2)),__STAT__) - ! super-band cloud fractions + ! super-layer cloud fractions allocate(CLEARCOUNTS (size(Q,1),4),__STAT__) ! output fluxes allocate(SWUFLX (size(Q,1),size(Q,2)+1),__STAT__) @@ -3699,6 +3757,7 @@ subroutine SORADCORE(IM,JM,LM,NO_AERO,CURRTIME,LoadBalance,RC) LM-LCLDLM+1, LM-LCLDMH+1, NORMFLX, & CLEARCOUNTS, SWUFLX, SWDFLX, SWUFLXC, SWDFLXC, & NIRR_R, NIRF_R, PARR_R, PARF_R, UVRR_R, UVRF_R,& + TAUTP, TAUHP, TAUMP, TAULP, & BNDSOLVAR, INDSOLVAR, SOLCYCFRAC) call MAPL_TimerOff(MAPL,"--RRTMG_RUN") @@ -3717,7 +3776,7 @@ subroutine SORADCORE(IM,JM,LM,NO_AERO,CURRTIME,LoadBalance,RC) ! required outputs ! ---------------- - ! convert super-band clearCounts to cloud fractions + ! convert super-layer clearCounts to cloud fractions CLDTS(:) = 1. - CLEARCOUNTS(:,1)/float(NGPTSW) CLDHS(:) = 1. - CLEARCOUNTS(:,2)/float(NGPTSW) CLDMS(:) = 1. - CLEARCOUNTS(:,3)/float(NGPTSW) @@ -3822,13 +3881,13 @@ subroutine SORADCORE(IM,JM,LM,NO_AERO,CURRTIME,LoadBalance,RC) ! internal 3D outputs call ESMFL_StateGetPointerToData(INTERNAL, ptr3, NamesOut(k), __RC__) call ReOrder(BufOut(L1),ptr3,daytime,NumMax,HorzDims,size(ptr3,3),UNPACKIT) - else if (k < NumOut-3) then + else if (k < NumOut-7) then ! internal 2D outputs call ESMFL_StateGetPointerToData(INTERNAL, ptr2, NamesOut(k), __RC__) call ReOrder(BufOut(L1),ptr2,daytime,NumMax,HorzDims,1,UNPACKIT) else ! cloud fraction outputs (2D) - if (NamesOut(k) == "CLDTTSW") then + if (NamesOut(k) == "CLDTTSW") then if (associated(CLDTTSW)) then call ReOrder(BufOut(L1),CLDTTSW,daytime,NumMax,HorzDims,1,UNPACKIT) WHERE (.not.daytime) CLDTTSW = MAPL_UNDEF @@ -3848,6 +3907,26 @@ subroutine SORADCORE(IM,JM,LM,NO_AERO,CURRTIME,LoadBalance,RC) call ReOrder(BufOut(L1),CLDLOSW,daytime,NumMax,HorzDims,1,UNPACKIT) WHERE (.not.daytime) CLDLOSW = MAPL_UNDEF end if + else if (NamesOut(k) == "TAUTTPAR") then + if (associated(TAUTTPAR)) then + call ReOrder(BufOut(L1),TAUTTPAR,daytime,NumMax,HorzDims,1,UNPACKIT) + WHERE (.not.daytime) TAUTTPAR = MAPL_UNDEF + end if + else if (NamesOut(k) == "TAUHIPAR") then + if (associated(TAUHIPAR)) then + call ReOrder(BufOut(L1),TAUHIPAR,daytime,NumMax,HorzDims,1,UNPACKIT) + WHERE (.not.daytime) TAUHIPAR = MAPL_UNDEF + end if + else if (NamesOut(k) == "TAUMDPAR") then + if (associated(TAUMDPAR)) then + call ReOrder(BufOut(L1),TAUMDPAR,daytime,NumMax,HorzDims,1,UNPACKIT) + WHERE (.not.daytime) TAUMDPAR = MAPL_UNDEF + end if + else if (NamesOut(k) == "TAULOPAR") then + if (associated(TAULOPAR)) then + call ReOrder(BufOut(L1),TAULOPAR,daytime,NumMax,HorzDims,1,UNPACKIT) + WHERE (.not.daytime) TAULOPAR = MAPL_UNDEF + end if end if end if @@ -4553,7 +4632,7 @@ subroutine UPDATE_EXPORT(IM,JM,LM, RC) ! TAU[HML] are correct because GETVISTAU produces in-cloud optical thicknesses for ! 'effective clouds' extended-out and diluted to the maximum cloud fraction in each - ! pressure super-band [LMT]. + ! pressure super-layers [LMH]. if (associated(TAUH) .or. associated(TAUT) .or. associated(TAUTX)) then allocate(aTAUH(IM,JM),__STAT__) @@ -4582,7 +4661,7 @@ subroutine UPDATE_EXPORT(IM,JM,LM, RC) if (associated(TAUL)) TAUL = aTAUL end if - ! TAUT however is broken because the three super-bands are randomly overlapped + ! TAUT however is broken because the three super-layers are randomly overlapped ! and with different effective cloud fractions. It has been broken but used for ! a long time. It should be considered deprecated. TAUTX below is an improved ! version. diff --git a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/RRTMG/rrtmg_sw/gcm_model/src/rrtmg_sw_cldprmc.F90 b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/RRTMG/rrtmg_sw/gcm_model/src/rrtmg_sw_cldprmc.F90 index 693bab6e4..5a64ee597 100644 --- a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/RRTMG/rrtmg_sw/gcm_model/src/rrtmg_sw_cldprmc.F90 +++ b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/RRTMG/rrtmg_sw/gcm_model/src/rrtmg_sw_cldprmc.F90 @@ -52,10 +52,10 @@ subroutine cldprmc_sw(pncol, ncol, nlay, iceflag, liqflag, & integer, intent(in) :: liqflag ! see definitions logical, intent(in) :: cldymc (nlay,ngptsw,pncol) ! cloudy or not? [mcica] - real, intent(in) :: ciwpmc (nlay,ngptsw,pncol) ! cloud ice water path [mcica] - real, intent(in) :: clwpmc (nlay,ngptsw,pncol) ! cloud liq water path [mcica] - real, intent(in) :: relqmc (nlay,pncol) ! cloud liq ptle eff radius (um) - real, intent(in) :: reicmc (nlay,pncol) ! cloud ice ptle eff radius (um) + real, intent(in) :: ciwpmc (nlay,ngptsw,pncol) ! cloud ice water path [mcica] [g/m2] + real, intent(in) :: clwpmc (nlay,ngptsw,pncol) ! cloud liq water path [mcica] [g/m2] + real, intent(in) :: relqmc (nlay,pncol) ! cloud liq ptle eff radius [um] + real, intent(in) :: reicmc (nlay,pncol) ! cloud ice ptle eff radius [um] ! Specific definition of reicmc depends on iceflag: ! iceflag = 1: ice effective radius, r_ec, (Ebert and Curry, 1992), @@ -93,6 +93,16 @@ subroutine cldprmc_sw(pncol, ncol, nlay, iceflag, liqflag, & extcoice, gice, ssacoice, forwice, & extcoliq, gliq, ssacoliq, forwliq + ! Notes by PMN: + ! The optical properties per unit ice (liquid) amount (e.g., extcoice) + ! are currently per band only but are redundantly re-calculated for + ! each g-point in the band that has non-zero mcica ice (liquid) amount. + ! This is inefficient if two or more band g-points have ice (liquid). + ! I have resisted the urge to remove this inefficiency because a future + ! improvement may well McICA-generate an ice (liquid) effective radius + ! that covaries with the condensate amount. In this more general case, + ! extcoice, etc., will vary with g-point within a band. + ! --------------------------------------------------------- ! Calculation of absorption coefficients due to ice clouds. ! --------------------------------------------------------- diff --git a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/RRTMG/rrtmg_sw/gcm_model/src/rrtmg_sw_rad.F90 b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/RRTMG/rrtmg_sw/gcm_model/src/rrtmg_sw_rad.F90 index fc6ea6c0f..c902f70a1 100644 --- a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/RRTMG/rrtmg_sw/gcm_model/src/rrtmg_sw_rad.F90 +++ b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/RRTMG/rrtmg_sw/gcm_model/src/rrtmg_sw_rad.F90 @@ -76,6 +76,7 @@ subroutine rrtmg_sw ( & cloudLM, cloudMH, normFlx, & clearCounts, swuflx, swdflx, swuflxc, swdflxc, & nirr, nirf, parr, parf, uvrr, uvrf, & + tautp, tauhp, taump, taulp, & ! optional inputs bndscl, indsolvar, solcycfrac) @@ -232,7 +233,7 @@ subroutine rrtmg_sw ( & ! ----- Outputs ----- - ! subcolumn clear counts for Tot|High|Mid|Low bands + ! subcolumn clear counts for Tot|High|Mid|Low super-layers integer, intent(out) :: clearCounts(ncol,4) real, intent(out) :: swuflx (ncol,nlay+1) ! Total sky SW up flux (W/m2) @@ -248,6 +249,9 @@ subroutine rrtmg_sw ( & real, intent(out) :: uvrr (ncol) ! UV direct down SW flux (W/m2) real, intent(out) :: uvrf (ncol) ! UV diffuse down SW flux (W/m2) + ! In-cloud PAR optical thickness for Tot|High|Mid|Low super-layers + real, intent(out), dimension (ncol) :: tautp, tauhp, taump, taulp + ! ----- Locals ----- integer :: pncol @@ -375,6 +379,7 @@ subroutine rrtmg_sw ( & cloudLM, cloudMH, normFlx, & clearCounts, swuflx, swdflx, swuflxc, swdflxc, & nirr, nirf, parr, parf, uvrr, uvrf, & + tautp, tauhp, taump, taulp, & ! optional inputs bndscl, indsolvar, solcycfrac) @@ -397,6 +402,7 @@ subroutine rrtmg_sw_sub ( & cloudLM, cloudMH, normFlx, & clearCounts, swuflx, swdflx, swuflxc, swdflxc, & nirr, nirf, parr, parf, uvrr, uvrf, & + tautp, tauhp, taump, taulp, & ! optional inputs bndscl, indsolvar, solcycfrac) @@ -473,7 +479,7 @@ subroutine rrtmg_sw_sub ( & real, intent(in) :: galdir (gncol) ! Near-IR surface albedo: direct rad real, intent(in) :: galdif (gncol) ! Near-IR surface albedo: diffuse rad - ! super-band cloud fraction boundaries + ! super-layer cloud fraction boundaries integer, intent(in) :: cloudLM ! Low-mid integer, intent(in) :: cloudMH ! Mid-high @@ -481,7 +487,7 @@ subroutine rrtmg_sw_sub ( & ! ----- Outputs ----- - ! subcolumn clear counts for Tot|High|Mid|Low bands + ! subcolumn clear counts for Tot|High|Mid|Low super-layers integer, intent(out) :: clearCounts(gncol,4) real, intent(out) :: swuflx (gncol,nlay+1) ! Total sky SW up flux (W/m2) @@ -497,6 +503,9 @@ subroutine rrtmg_sw_sub ( & real, intent(out) :: uvrr (gncol) ! UV direct down SW flux (w/m2) real, intent(out) :: uvrf (gncol) ! UV diffuse down SW flux (w/m2) + ! In-cloud PAR optical thickness for Tot|High|Mid|Low super-layers + real, intent(out), dimension (gncol) :: tautp, tauhp, taump, taulp + ! ----- Locals ----- ! Control @@ -563,9 +572,9 @@ subroutine rrtmg_sw_sub ( & real :: zm (nlay,pncol) ! mid-layer hgt for cld overlap [m] logical :: cldymcl (nlay,ngptsw,pncol) ! cloud or not? [mcica] - real :: ciwpmcl (nlay,ngptsw,pncol) ! in-cloud ice water path [mcica] - real :: clwpmcl (nlay,ngptsw,pncol) ! in-cloud liquid water path [mcica] - integer :: p_clearCounts (4,pncol) ! for super-band cld fractions + real :: ciwpmcl (nlay,ngptsw,pncol) ! in-cloud ice water path [mcica] [g/m2] + real :: clwpmcl (nlay,ngptsw,pncol) ! in-cloud liq water path [mcica] [g/m2] + integer :: p_clearCounts (4,pncol) ! for super-layer cld fractions real :: taucmc (nlay,ngptsw,pncol) ! in-cloud optical depth [mcica] real :: taormc (nlay,ngptsw,pncol) ! unscaled in-cloud optl depth [mcica] @@ -597,6 +606,9 @@ subroutine rrtmg_sw_sub ( & real, dimension (pncol) :: & znirr, znirf, zparr, zparf, zuvrr, zuvrf + + ! in-cloud PAR optical thicknesses + real, dimension (pncol) :: ztautp, ztauhp, ztaump, ztaulp ! Solar variability multipliers ! ----------------------------- @@ -1129,7 +1141,7 @@ subroutine rrtmg_sw_sub ( & cldymcl, ciwpmcl, clwpmcl, & seed_order=[4,3,2,1]) - ! for super-band cloud fractions + ! for super-layer cloud fractions call clearCounts_threeBand( & pncol, ncol, ngptsw, nlay, cloudLM, cloudMH, cldymcl, & p_clearCounts) @@ -1163,10 +1175,12 @@ subroutine rrtmg_sw_sub ( & laytrop, jp, jt, jt1, & colch4, colco2, colh2o, colmol, colo2, colo3, & fac00, fac01, fac10, fac11, & + cloudLM, cloudMH, & selffac, selffrac, indself, forfac, forfrac, indfor, & zbbfd, zbbfu, zbbcd, zbbcu, zuvfd, zuvcd, znifd, znicd, & zbbfddir, zbbcddir, zuvfddir, zuvcddir, znifddir, znicddir,& - znirr,znirf,zparr,zparf,zuvrr,zuvrf) + znirr, znirf, zparr, zparf, zuvrr, zuvrf, & + ztautp, ztauhp, ztaump, ztaulp) ! Copy out up and down, clear and total sky fluxes to output arrays. ! Vertical indexing goes from bottom to top; reverse here for GCM if necessary. @@ -1176,7 +1190,7 @@ subroutine rrtmg_sw_sub ( & do icol = 1,ncol gicol = gicol_clr(icol + cols - 1) - ! super-band clear counts + ! super-layer clear counts do n = 1,4 clearCounts (gicol,n) = ngptsw end do @@ -1189,6 +1203,12 @@ subroutine rrtmg_sw_sub ( & swdflx (gicol,ilev) = zbbfd(ilev,icol) enddo + ! super-layer optical thicknesses + tautp(gicol) = 0. + tauhp(gicol) = 0. + taump(gicol) = 0. + taulp(gicol) = 0. + enddo ! surface broadband fluxes @@ -1215,7 +1235,10 @@ subroutine rrtmg_sw_sub ( & swuflx (gicol,ilev) = zbbfu(ilev,icol) swdflx (gicol,ilev) = zbbfd(ilev,icol) enddo - + tautp(gicol) = ztautp(icol) + tauhp(gicol) = ztauhp(icol) + taump(gicol) = ztaump(icol) + taulp(gicol) = ztaulp(icol) enddo do icol = 1,ncol diff --git a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/RRTMG/rrtmg_sw/gcm_model/src/rrtmg_sw_spcvmc.F90 b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/RRTMG/rrtmg_sw/gcm_model/src/rrtmg_sw_spcvmc.F90 index 6b4fea3fc..23df82d80 100644 --- a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/RRTMG/rrtmg_sw/gcm_model/src/rrtmg_sw_spcvmc.F90 +++ b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/RRTMG/rrtmg_sw/gcm_model/src/rrtmg_sw_spcvmc.F90 @@ -36,10 +36,12 @@ subroutine spcvmc_sw ( & laytrop, jp, jt, jt1, & colch4, colco2, colh2o, colmol, colo2, colo3, & fac00, fac01, fac10, fac11, & + cloudLM, cloudMH, & selffac, selffrac, indself, forfac, forfrac, indfor, & pbbfd, pbbfu, pbbcd, pbbcu, puvfd, puvcd, pnifd, pnicd, & pbbfddir, pbbcddir, puvfddir, puvcddir, pnifddir, pnicddir, & - znirr, znirf, zparr, zparf, zuvrr, zuvrf) + znirr, znirf, zparr, zparf, zuvrr, zuvrf, & + ztautp, ztauhp, ztaump, ztaulp) ! --------------------------------------------------------------------------- ! ! Purpose: Contains spectral loop to compute the shortwave radiative fluxes, @@ -129,6 +131,10 @@ subroutine spcvmc_sw ( & real, intent(in), dimension (nlay,pncol) & :: fac00, fac01, fac10, fac11 + ! pressure super-layer interface levels for optical thicknesses + integer, intent(in) :: cloudLM ! Low-mid + integer, intent(in) :: cloudMH ! Mid-high + ! ------- Output ------- real, intent(out) :: pbbcd (nlay+1,pncol) @@ -148,16 +154,23 @@ subroutine spcvmc_sw ( & real, intent(out) :: pnicddir (nlay+1,pncol) real, intent(out) :: pnifddir (nlay+1,pncol) + ! surface band fluxes (direct and TOTAL) real, intent(out), dimension(pncol) :: & znirr, znirf, zparr, zparf, zuvrr, zuvrf + ! in-cloud PAR optical thicknesses + real, intent(out), dimension (pncol) :: & + ztautp, ztauhp, ztaump, ztaulp + ! ------- Local ------- integer :: icol integer :: jk, ikl integer :: iw, jb, ibm - real :: zf, zwf, zincflx + real :: zf, zwf, zincflx, wgt + real :: stautp, stauhp, staump, staulp + real :: wtautp, wtauhp, wtaump, wtaulp real :: zgco (nlay,ngptsw,pncol) real :: zomco (nlay,ngptsw,pncol) @@ -475,6 +488,108 @@ subroutine spcvmc_sw ( & end do enddo + ! diagnostic in-cloud optical thicknesses in PAR super-band + ! (weighted across and within bands by TOA incident flux) + ! ------------------------------------------------------- + ztautp = 0. + ztauhp = 0. + ztaump = 0. + ztaulp = 0. + + ! can only be non-zero for potentially cloudy columns + if (cc == 2) then + + do icol = 1,ncol + + wtautp = 0. + wtauhp = 0. + wtaump = 0. + wtaulp = 0. + + do iw = 1,ngptsw + jb = ngb(iw) + ibm = jb - 15 + + ! band weights + if (ibm >= 10 .and. ibm <= 11) then + ! Photosynthetically active radiation (PAR) + wgt = 1.0 + else if (ibm == 9) then + ! half PAR, half NIR + wgt = 0.5 + else + ! no contribution to PAR + cycle + end if + + ! TOA flux weighting (adjustment for zenith angle + ! not needed since normalized for each icol) + if (isolvar < 0) then + zincflx = adjflux(jb) * zsflxzen(iw,icol) + else + zincflx = adjflux(jb) * ssi (iw,icol) + endif + wgt = wgt * zincflx + + ! low pressure layer + staulp = sum(ptaormc(1:cloudLM,iw,icol),dim=1) + if (staulp > 0.) then + ztaulp(icol) = ztaulp(icol) + wgt * staulp + wtaulp = wtaulp + wgt + end if + + ! mid pressure layer + staump = sum(ptaormc(cloudLM+1:cloudMH,iw,icol),dim=1) + if (staump > 0.) then + ztaump(icol) = ztaump(icol) + wgt * staump + wtaump = wtaump + wgt + end if + + ! high pressure layer + stauhp = sum(ptaormc(cloudMH+1:nlay,iw,icol),dim=1) + if (stauhp > 0.) then + ztauhp(icol) = ztauhp(icol) + wgt * stauhp + wtauhp = wtauhp + wgt + end if + + ! whole subcolumn + stautp = staulp + staump + stauhp + if (stautp > 0.) then + ztautp(icol) = ztautp(icol) + wgt * stautp + wtautp = wtautp + wgt + end if + + end do ! iw + + ! normalize + if (wtautp > 0.) then + ztautp(icol) = ztautp(icol) / wtautp + else + ztautp(icol) = 0. + end if + + if (wtauhp > 0.) then + ztauhp(icol) = ztauhp(icol) / wtauhp + else + ztauhp(icol) = 0. + end if + + if (wtaump > 0.) then + ztaump(icol) = ztaump(icol) / wtaump + else + ztaump(icol) = 0. + end if + + if (wtaulp > 0.) then + ztaulp(icol) = ztaulp(icol) / wtaulp + else + ztaulp(icol) = 0. + end if + + end do ! icol + + end if ! cc==2 + end subroutine spcvmc_sw ! -------------------------------------------------------------------- From 78eac2f86db51cc335fb3ba6c93532ad26590f31 Mon Sep 17 00:00:00 2001 From: Peter Norris Date: Tue, 11 Jan 2022 18:03:49 -0500 Subject: [PATCH 3/9] pmn: add REFRESH-frequency COSZSW as seen by SolarGC REFRESH --- .../GEOSsolar_GridComp/GEOS_SolarGridComp.F90 | 16 ++++++++++++++++ 1 file changed, 16 insertions(+) diff --git a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 index 2d6e6607b..8d8659c3f 100644 --- a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 +++ b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 @@ -1255,6 +1255,13 @@ subroutine SetServices ( GC, RC ) DIMS = MAPL_DimsHorzOnly, & VLOCATION = MAPL_VLocationNone, __RC__) + call MAPL_AddExportSpec(GC, & + LONG_NAME = 'cosine_of_the_solar_zenith_angle_of_Solar_REFRESH', & + UNITS = '1', & + SHORT_NAME = 'COSZSW', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__) + call MAPL_AddExportSpec(GC, & SHORT_NAME = 'CLDTMP', & LONG_NAME = 'cloud_top_temperature', & @@ -2022,6 +2029,9 @@ subroutine SORADCORE(IM,JM,LM,NO_AERO,CURRTIME,LoadBalance,RC) real, pointer, dimension(:,:) :: TAUMDPAR real, pointer, dimension(:,:) :: TAULOPAR + ! cosine solar zenith angle used by REFRESH + real, pointer, dimension(:,:) :: COSZSW + ! DAYTIME ONLY COPY OF VARIABLES real, pointer, dimension(: ) :: ALBNR,ALBNF,ALBVR,ALBVF,ZT,SLR1D, & @@ -2232,6 +2242,9 @@ subroutine SORADCORE(IM,JM,LM,NO_AERO,CURRTIME,LoadBalance,RC) end do end do + ! cosine solar zenith angle used by REFRESH + call MAPL_GetPointer(EXPORT, COSZSW, 'COSZSW', __RC__) + ! super-layer RRTMG cloud fraction exports call MAPL_GetPointer(EXPORT, CLDTTSW, 'CLDTTSW', __RC__) call MAPL_GetPointer(EXPORT, CLDHISW, 'CLDHISW', __RC__) @@ -2270,6 +2283,9 @@ subroutine SORADCORE(IM,JM,LM,NO_AERO,CURRTIME,LoadBalance,RC) NumLit = size(ZTH) end if + ! write out the cosine solar zenith angle actually used by REFRESH + if (associated(COSZSW)) COSZSW = ZTH + ! Create a balancing strategy. This is a collective call on the communicator ! of the current VM. The original, unbalanced local work consists of (OrgLen) ! NumLit soundings, which may be zero. The local work after implementing the From 342651b59fe7316aacbe7218a695cabd9a6a01bc Mon Sep 17 00:00:00 2001 From: Peter Norris Date: Sun, 16 Jan 2022 17:33:49 -0500 Subject: [PATCH 4/9] pmn: vertical ordering agnostic subcol gen and add CLDxxSWHB diagnostics --- .../GEOS_RadiationShared/cloud_subcol_gen.F90 | 187 ++++++++-- .../GEOSirrad_GridComp/GEOS_IrradGridComp.F90 | 100 ++--- .../GEOSsolar_GridComp/GEOS_SolarGridComp.F90 | 343 +++++++++++++++--- 3 files changed, 491 insertions(+), 139 deletions(-) diff --git a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOS_RadiationShared/cloud_subcol_gen.F90 b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOS_RadiationShared/cloud_subcol_gen.F90 index 5de2f0197..9bbfca36e 100644 --- a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOS_RadiationShared/cloud_subcol_gen.F90 +++ b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOS_RadiationShared/cloud_subcol_gen.F90 @@ -18,6 +18,18 @@ module cloud_subcol_gen ! cloud output explicit by changing to a logical cldy_stoch in {true,false}, ! which may permit some exterior efficiency improvements. + ! PMN 2022/01: Original version only worked for layer ordering from 1=bottom + ! upward. New version auto-detects vertical ordering and works for either. + ! See details under generate_stochastic_clouds() and clearCounts_threeBand(). + ! NB: The POPULATION statistics of the subcolumn ensemble should not depend + ! on the vertical ordering. But if you want an exact SAMPLE replication of + ! a generation (such as the bottom-up generation used in RRTMG) then you + ! should ensure the same input ordering (in the RRTMG case, by reversing the + ! vertical ordering from the GEOS-5 top-down ordering BEFORE input to these + ! routines). Otherwise (i.e. if POPULATION statistics equivalence is OK) + ! then save time by avoiding any explicit array re-ordering and let the + ! auto-detect feature do its job. + use cloud_condensate_inhomogeneity, only: & condensate_inhomogeneous, zcw_lookup use iso_fortran_env, only : error_unit @@ -156,7 +168,7 @@ subroutine generate_stochastic_clouds( & ! of the initial seeding on each gridcolumn's PRNG stream. As such, it provides a ! way of generating a different set of random numbers for say LW and SW at the same ! timetep, when the layer pressures used for seed generation are the same. - ! For example, LW can be run with [1,2,3,4] and SW with [4,2,1,3]. + ! For example, LW can be run with [1,2,3,4] and SW with [4,3,2,1]. ! NOTE: Using the same seed_order for LW and SW still does not mean using the "same ! subcolumn cloud ensemble" because the number of subcolumns (g-points) in LW and SW ! differ. Also, of course, McICA .ne. ICA, meaning that each subcolumn in McICA gets @@ -167,6 +179,14 @@ subroutine generate_stochastic_clouds( & ! seeds for LW and SW will create an unwanted correlation of some sort, they can use ! a different seed_order for each. ! + ! Vertical ordering: + ! The ordering is determined by play(1,1) > play(nlay,1) for surface=1 upward, and vica + ! versa. Either vertical ordering is permitted, but the generation always proceeds from + ! layer 1 to layer nlay. That shouldn't matter in a statistical sense. The seed order + ! does NOT need to be changed depending on the vertical ordering. The four layers nearest + ! the surface are always used (in the seed_order sense). Of course, the vertical ordering + ! of all profile inputs (zmid,play,cldfrac,ciwp,clwp) should be the same, and the output + ! vertical ordering (of [cldy|ciwp|clwp]_stoch) is the same as the input ordering. !--------------------------------------------------------------------------------------- integer, intent(in) :: dncol ! Dimensioned number of gridcols @@ -174,6 +194,8 @@ subroutine generate_stochastic_clouds( & integer, intent(in) :: nsubcol ! #Subcols to generate / gridcol integer, intent(in) :: nlay ! Number of model layers + ! Note: only absolute layer separations, abs(zmid(k1,:)-zmid(k2,:)), matter + ! so "zero-reference" of zmid is not important real, intent(in) :: zmid (nlay,dncol) ! Height of midpoints [m] real, intent(in) :: alat (dncol) ! Latitude of gridcolumn [radians] integer, intent(in) :: doy ! Day of year @@ -206,6 +228,9 @@ subroutine generate_stochastic_clouds( & ! ----- Locals ----- + ! vertical ordering + logical :: surface_at_one + ! correlation length for cloud presence and condensate real, dimension(dncol) :: adl, rdl @@ -235,6 +260,9 @@ subroutine generate_stochastic_clouds( & ! indices integer :: icol, isubcol, ilay, n + ! detect vertical ordering + surface_at_one = (play(1,1) > play(nlay,1)) + ! save for speed cond_inhomo = condensate_inhomogeneous() @@ -284,12 +312,13 @@ subroutine generate_stochastic_clouds( & do icol = 1,ncol ! exponential inter-layer correlations ... + ! [alpha|rcorr](k) is correlation between layers k and k-1 do ilay = 2,nlay - alpha(ilay) = exp( -(zmid(ilay,icol)-zmid(ilay-1,icol)) / adl(icol) ) + alpha(ilay) = exp( -abs(zmid(ilay,icol)-zmid(ilay-1,icol)) / adl(icol) ) enddo if (cond_inhomo) then do ilay = 2,nlay - rcorr(ilay) = exp( -(zmid(ilay,icol)-zmid(ilay-1,icol)) / rdl(icol) ) + rcorr(ilay) = exp( -abs(zmid(ilay,icol)-zmid(ilay-1,icol)) / rdl(icol) ) enddo endif @@ -345,7 +374,12 @@ subroutine generate_stochastic_clouds( & ! pseed. ! isolate lower pressures for seeding - pseed = play(1:4,icol) * 100. ! [Pa] + if (surface_at_one) then + pseed = play(1:4,icol) + else + pseed = play(nlay:nlay-3:-1,icol) + end if + pseed = pseed * 100. ! [Pa] ! PMN 2021-11-12: Try [daPa=decaPascals] instead. ! pseed = play(1:4,icol) * 10. ! [daPa] @@ -557,9 +591,24 @@ subroutine clearCounts_threeBand( & ! ---------------------------------------------------------------- ! Count clear subcolumns per gridcolumn for whole column ! and for three pressure bands. - ! layers [1, cloudLM] are in low pressure band - ! layers [cloudLM+1, cloudMH] are in mid pressure band - ! layers [cloudMH+1, nlay ] are in high pressure band + ! Verical ordering is auto-detected via: + ! cloudLM < cloudMH => surface at level 1 + ! cloudLM > cloudMH => surface at level nlay + ! and concomitant layer definitions: + ! layers [1, cloudLM] are in low pressure band + ! layers [cloudLM+1, cloudMH] are in mid pressure band + ! layers [cloudMH+1, nlay ] are in high pressure band + ! for surface at level 1, + ! with cloudLM = nlay - LCLDLM + 1 + ! and cloudMH = nlay - LCLDMH + 1, + ! and + ! layers [1, cloudMH-1] are in high pressure band + ! layers [cloudMH, cloudLM-1] are in mid pressure band + ! layers [cloudLM, nlay ] are in low pressure band + ! for surface at level nlay (TOA at level 1) + ! with cloudLM = LCLDLM + ! and cloudMH = LCLDMH, + ! all in terms of GEOS-5's LCLDLM and LCLDMH (TOA at level 1). ! ---------------------------------------------------------------- integer, intent(in) :: dncol ! Dimensioned number of gridcols @@ -597,41 +646,101 @@ subroutine clearCounts_threeBand( & if (.not. cloud_found) & clearCnts(1,icol) = clearCnts(1,icol) + 1 - ! high pressure band - cloud_found = .false. - do ilay = cloudMH+1, nlay - if (cldy_stoch(ilay,isubcol,icol)) then - cloud_found = .true. - exit - endif - enddo - if (.not. cloud_found) & - clearCnts(2,icol) = clearCnts(2,icol) + 1 + enddo ! subcolumns + enddo ! gridcolumns + + if (cloudLM < cloudMH) then - ! mid pressure band - cloud_found = .false. - do ilay = cloudLM+1, cloudMH - if (cldy_stoch(ilay,isubcol,icol)) then - cloud_found = .true. - exit - endif - enddo - if (.not. cloud_found) & - clearCnts(3,icol) = clearCnts(3,icol) + 1 + ! surface at level 1 ordering ... - ! low pressure band - cloud_found = .false. - do ilay = 1, cloudLM - if (cldy_stoch(ilay,isubcol,icol)) then - cloud_found = .true. - exit - endif - enddo - if (.not. cloud_found) & - clearCnts(4,icol) = clearCnts(4,icol) + 1 + do icol = 1,ncol + do isubcol = 1,nsubcol - enddo ! subcolumns - enddo ! gridcolumns + ! low pressure band + cloud_found = .false. + do ilay = 1, cloudLM + if (cldy_stoch(ilay,isubcol,icol)) then + cloud_found = .true. + exit + endif + enddo + if (.not. cloud_found) & + clearCnts(4,icol) = clearCnts(4,icol) + 1 + + ! mid pressure band + cloud_found = .false. + do ilay = cloudLM+1, cloudMH + if (cldy_stoch(ilay,isubcol,icol)) then + cloud_found = .true. + exit + endif + enddo + if (.not. cloud_found) & + clearCnts(3,icol) = clearCnts(3,icol) + 1 + + ! high pressure band + cloud_found = .false. + do ilay = cloudMH+1, nlay + if (cldy_stoch(ilay,isubcol,icol)) then + cloud_found = .true. + exit + endif + enddo + if (.not. cloud_found) & + clearCnts(2,icol) = clearCnts(2,icol) + 1 + + enddo ! subcolumns + enddo ! gridcolumns + + elseif (cloudLM > cloudMH) then + + ! TOA at level 1 ordering ... + + do icol = 1,ncol + do isubcol = 1,nsubcol + + ! high pressure band + cloud_found = .false. + do ilay = 1, cloudMH-1 + if (cldy_stoch(ilay,isubcol,icol)) then + cloud_found = .true. + exit + endif + enddo + if (.not. cloud_found) & + clearCnts(2,icol) = clearCnts(2,icol) + 1 + + ! mid pressure band + cloud_found = .false. + do ilay = cloudMH, cloudLM-1 + if (cldy_stoch(ilay,isubcol,icol)) then + cloud_found = .true. + exit + endif + enddo + if (.not. cloud_found) & + clearCnts(3,icol) = clearCnts(3,icol) + 1 + + ! low pressure band + cloud_found = .false. + do ilay = cloudLM, nlay + if (cldy_stoch(ilay,isubcol,icol)) then + cloud_found = .true. + exit + endif + enddo + if (.not. cloud_found) & + clearCnts(4,icol) = clearCnts(4,icol) + 1 + + enddo ! subcolumns + enddo ! gridcolumns + + else + + write(error_unit,*) 'file:', __FILE__, ', line:', __LINE__ + error stop 'invalid pressure super-layers!' + + end if end subroutine clearCounts_threeBand diff --git a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSirrad_GridComp/GEOS_IrradGridComp.F90 b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSirrad_GridComp/GEOS_IrradGridComp.F90 index f201db8a2..6b8d83ce0 100644 --- a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSirrad_GridComp/GEOS_IrradGridComp.F90 +++ b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSirrad_GridComp/GEOS_IrradGridComp.F90 @@ -783,33 +783,33 @@ subroutine SetServices ( GC, RC ) DIMS = MAPL_DimsHorzOnly, & VLOCATION = MAPL_VLocationNone, __RC__ ) - call MAPL_AddExportSpec(GC, & - SHORT_NAME = 'CLDTTLW', & - LONG_NAME = 'total_cloud_area_fraction_rrtmg_lw', & - UNITS = '1', & - DIMS = MAPL_DimsHorzOnly, & - VLOCATION = MAPL_VLocationNone, __RC__ ) - - call MAPL_AddExportSpec(GC, & - SHORT_NAME = 'CLDHILW', & - LONG_NAME = 'high-level_cloud_area_fraction_rrtmg_lw', & - UNITS = '1', & - DIMS = MAPL_DimsHorzOnly, & - VLOCATION = MAPL_VLocationNone, __RC__ ) - - call MAPL_AddExportSpec(GC, & - SHORT_NAME = 'CLDMDLW', & - LONG_NAME = 'mid-level_cloud_area_fraction_rrtmg_lw', & - UNITS = '1', & - DIMS = MAPL_DimsHorzOnly, & - VLOCATION = MAPL_VLocationNone, __RC__ ) - - call MAPL_AddExportSpec(GC, & - SHORT_NAME = 'CLDLOLW', & - LONG_NAME = 'low_level_cloud_area_fraction_rrtmg_lw', & - UNITS = '1', & - DIMS = MAPL_DimsHorzOnly, & - VLOCATION = MAPL_VLocationNone, __RC__ ) + call MAPL_AddExportSpec(GC, & + SHORT_NAME = 'CLDTTLW', & + LONG_NAME = 'total_cloud_area_fraction_rrtmg_lw_REFRESH', & + UNITS = '1', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__ ) + + call MAPL_AddExportSpec(GC, & + SHORT_NAME = 'CLDHILW', & + LONG_NAME = 'high-level_cloud_area_fraction_rrtmg_lw_REFRESH', & + UNITS = '1', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__ ) + + call MAPL_AddExportSpec(GC, & + SHORT_NAME = 'CLDMDLW', & + LONG_NAME = 'mid-level_cloud_area_fraction_rrtmg_lw_REFRESH', & + UNITS = '1', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__ ) + + call MAPL_AddExportSpec(GC, & + SHORT_NAME = 'CLDLOLW', & + LONG_NAME = 'low_level_cloud_area_fraction_rrtmg_lw_REFRESH', & + UNITS = '1', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__ ) ! Irrad does not have a "real" internal state. To update the net_longwave_flux ! due to the change of surface temperature every time step, we keep @@ -1701,27 +1701,39 @@ subroutine Lw_Driver(IM,JM,LM,LATS,LONS,CoresPerNode,RC) REFF(:,:,:,KRAIN ) = RR * 1.0e6 REFF(:,:,:,KSNOW ) = RS * 1.0e6 -! Determine the model level seperating high and middle clouds -!------------------------------------------------------------ +! Determine the model level separating high-middle and low-middle clouds +!----------------------------------------------------------------------- - LCLDMH = 1 - do K = 1, LM - if( PREF(K) >= PRS_MID_HIGH ) then - LCLDMH = K - exit - end if - end do + _ASSERT(PRS_MID_HIGH > PREF(1) , 'mid-high pressure band boundary too high!') + _ASSERT(PRS_LOW_MID > PRS_MID_HIGH, 'pressure band misordering!') + _ASSERT(PRS_LOW_MID < PREF(LM) , 'low-mid pressure band boundary too low!') -! Determine the model level seperating low and middle clouds -!----------------------------------------------------------- + ! find mid-high interface level + k = 1 + do while ( PREF(k) < PRS_MID_HIGH ) + k=k+1 + end do + LCLDMH = k + ! Guaranteed that LCLDMH > 1 (by first ASSERT above) + ! and that PREF(LCLDMH) >= PRS_MID_HIGH (by while loop) - LCLDLM = LM - do K = 1, LM - if( PREF(K) >= PRS_LOW_MID ) then - LCLDLM = K - exit - end if + ! find low-mid interface level + do while ( PREF(k) < PRS_LOW_MID ) + k=k+1 end do + LCLDLM = k + ! Guaranteed that LCLDLM <= LM (by third assert above) + ! and that PREF(LCLDLM) >= PRS_LOW_MID (by while loop) + + ! But it's still possible that LCLDLM == LCLDMH if the + ! interface pressures are too close. We now ASSERT to + ! prevent this. + _ASSERT(LCLDMH < LCLDLM, 'PRS_LOW_MID and PRS_MID_HIGH are too close!') + + ! now we have 1 < LCLDMH < LCLDLM <= LM and can use: + ! layers [1, LCLDMH-1] are in high pressure band + ! layers [LCLDMH, LCLDLM-1] are in mid pressure band + ! layers [LCLDLM, LM ] are in low pressure band call MAPL_GetPointer(EXPORT, CLDTTLW, 'CLDTTLW', RC=STATUS); VERIFY_(STATUS) call MAPL_GetPointer(EXPORT, CLDHILW, 'CLDHILW', RC=STATUS); VERIFY_(STATUS) diff --git a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 index 8d8659c3f..bccbbab9a 100644 --- a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 +++ b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 @@ -938,33 +938,61 @@ subroutine SetServices ( GC, RC ) DIMS = MAPL_DimsHorzOnly, & VLOCATION = MAPL_VLocationNone, __RC__) - call MAPL_AddExportSpec(GC, & - SHORT_NAME = 'CLDTTSW', & - LONG_NAME = 'total_cloud_area_fraction_rrtmg_sw', & - UNITS = '1', & - DIMS = MAPL_DimsHorzOnly, & - VLOCATION = MAPL_VLocationNone, __RC__ ) - - call MAPL_AddExportSpec(GC, & - SHORT_NAME = 'CLDHISW', & - LONG_NAME = 'high-level_cloud_area_fraction_rrtmg_sw', & - UNITS = '1', & - DIMS = MAPL_DimsHorzOnly, & - VLOCATION = MAPL_VLocationNone, __RC__ ) - - call MAPL_AddExportSpec(GC, & - SHORT_NAME = 'CLDMDSW', & - LONG_NAME = 'mid-level_cloud_area_fraction_rrtmg_sw', & - UNITS = '1', & - DIMS = MAPL_DimsHorzOnly, & - VLOCATION = MAPL_VLocationNone, __RC__ ) - - call MAPL_AddExportSpec(GC, & - SHORT_NAME = 'CLDLOSW', & - LONG_NAME = 'low-level_cloud_area_fraction_rrtmg_sw', & - UNITS = '1', & - DIMS = MAPL_DimsHorzOnly, & - VLOCATION = MAPL_VLocationNone, __RC__ ) + call MAPL_AddExportSpec(GC, & + SHORT_NAME = 'CLDTTSW', & + LONG_NAME = 'total_cloud_area_fraction_rrtmg_sw_REFRESH', & + UNITS = '1', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__ ) + + call MAPL_AddExportSpec(GC, & + SHORT_NAME = 'CLDHISW', & + LONG_NAME = 'high-level_cloud_area_fraction_rrtmg_sw_REFRESH', & + UNITS = '1', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__ ) + + call MAPL_AddExportSpec(GC, & + SHORT_NAME = 'CLDMDSW', & + LONG_NAME = 'mid-level_cloud_area_fraction_rrtmg_sw_REFRESH', & + UNITS = '1', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__ ) + + call MAPL_AddExportSpec(GC, & + SHORT_NAME = 'CLDLOSW', & + LONG_NAME = 'low-level_cloud_area_fraction_rrtmg_sw_REFRESH', & + UNITS = '1', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__ ) + + call MAPL_AddExportSpec(GC, & + SHORT_NAME = 'CLDTTSWHB', & + LONG_NAME = 'total_cloud_area_fraction_rrtmg_sw_HEARTBEAT', & + UNITS = '1', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__ ) + + call MAPL_AddExportSpec(GC, & + SHORT_NAME = 'CLDHISWHB', & + LONG_NAME = 'high-level_cloud_area_fraction_rrtmg_sw_HEARTBEAT', & + UNITS = '1', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__ ) + + call MAPL_AddExportSpec(GC, & + SHORT_NAME = 'CLDMDSWHB', & + LONG_NAME = 'mid-level_cloud_area_fraction_rrtmg_sw_HEARTBEAT', & + UNITS = '1', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__ ) + + call MAPL_AddExportSpec(GC, & + SHORT_NAME = 'CLDLOSWHB', & + LONG_NAME = 'low-level_cloud_area_fraction_rrtmg_sw_HEARTBEAT', & + UNITS = '1', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__ ) call MAPL_AddExportSpec(GC, & LONG_NAME = 'in_cloud_optical_thickness_of_low_clouds', & @@ -1029,33 +1057,33 @@ subroutine SetServices ( GC, RC ) DIMS = MAPL_DimsHorzVert, & VLOCATION = MAPL_VLocationCenter, __RC__) - call MAPL_AddExportSpec(GC, & - LONG_NAME = 'in_cloud_optical_thickness_of_low_clouds_RRTMG_PAR', & - UNITS = '1' , & - SHORT_NAME = 'TAULOPAR', & - DIMS = MAPL_DimsHorzOnly, & - VLOCATION = MAPL_VLocationNone, __RC__) - - call MAPL_AddExportSpec(GC, & - LONG_NAME = 'in_cloud_optical_thickness_of_middle_clouds_RRTMG_PAR', & - UNITS = '1' , & - SHORT_NAME = 'TAUMDPAR', & - DIMS = MAPL_DimsHorzOnly, & - VLOCATION = MAPL_VLocationNone, __RC__) - - call MAPL_AddExportSpec(GC, & - LONG_NAME = 'in_cloud_optical_thickness_of_high_clouds_RRTMG_PAR', & - UNITS = '1' , & - SHORT_NAME = 'TAUHIPAR', & - DIMS = MAPL_DimsHorzOnly, & - VLOCATION = MAPL_VLocationNone, __RC__) - - call MAPL_AddExportSpec(GC, & - LONG_NAME = 'in_cloud_optical_thickness_of_all_clouds_RRTMG_PAR', & - UNITS = '1' , & - SHORT_NAME = 'TAUTTPAR', & - DIMS = MAPL_DimsHorzOnly, & - VLOCATION = MAPL_VLocationNone, __RC__) + call MAPL_AddExportSpec(GC, & + LONG_NAME = 'in_cloud_optical_thickness_of_low_clouds_RRTMG_PAR_REFRESH', & + UNITS = '1' , & + SHORT_NAME = 'TAULOPAR', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__) + + call MAPL_AddExportSpec(GC, & + LONG_NAME = 'in_cloud_optical_thickness_of_middle_clouds_RRTMG_PAR_REFRESH', & + UNITS = '1' , & + SHORT_NAME = 'TAUMDPAR', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__) + + call MAPL_AddExportSpec(GC, & + LONG_NAME = 'in_cloud_optical_thickness_of_high_clouds_RRTMG_PAR_REFRESH', & + UNITS = '1' , & + SHORT_NAME = 'TAUHIPAR', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__) + + call MAPL_AddExportSpec(GC, & + LONG_NAME = 'in_cloud_optical_thickness_of_all_clouds_RRTMG_PAR_REFRESH', & + UNITS = '1' , & + SHORT_NAME = 'TAUTTPAR', & + DIMS = MAPL_DimsHorzOnly, & + VLOCATION = MAPL_VLocationNone, __RC__) call MAPL_AddExportSpec(GC, & LONG_NAME = 'surface_net_downward_shortwave_flux_assuming_clear_sky',& @@ -1698,22 +1726,39 @@ subroutine RUN ( GC, IMPORT, EXPORT, CLOCK, RC ) call MAPL_GetPointer(IMPORT, PREF, 'PREF', __RC__) -! Determine the model level seperating high-middle and low-middle clouds +! Determine the model level separating high-middle and low-middle clouds !----------------------------------------------------------------------- _ASSERT(PRS_MID_HIGH > PREF(1) , 'mid-high pressure band boundary too high!') _ASSERT(PRS_LOW_MID > PRS_MID_HIGH, 'pressure band misordering!') _ASSERT(PRS_LOW_MID < PREF(LM) , 'low-mid pressure band boundary too low!') + ! find mid-high interface level k = 1 do while ( PREF(k) < PRS_MID_HIGH ) k=k+1 end do LCLDMH = k - do while ( PREF(k) < PRS_LOW_MID ) + ! Guaranteed that LCLDMH > 1 (by first ASSERT above) + ! and that PREF(LCLDMH) >= PRS_MID_HIGH (by while loop) + + ! find low-mid interface level + do while ( PREF(k) < PRS_LOW_MID ) k=k+1 end do LCLDLM = k + ! Guaranteed that LCLDLM <= LM (by third assert above) + ! and that PREF(LCLDLM) >= PRS_LOW_MID (by while loop) + + ! But it's still possible that LCLDLM == LCLDMH if the + ! interface pressures are too close. We now ASSERT to + ! prevent this. + _ASSERT(LCLDMH < LCLDLM, 'PRS_LOW_MID and PRS_MID_HIGH are too close!') + + ! now we have 1 < LCLDMH < LCLDLM <= LM and can use: + ! layers [1, LCLDMH-1] are in high pressure band + ! layers [LCLDMH, LCLDLM-1] are in mid pressure band + ! layers [LCLDLM, LM ] are in low pressure band call MAPL_TimerOff(MAPL,"PRELIMS",__RC__) @@ -4428,6 +4473,22 @@ subroutine UPDATE_EXPORT(IM,JM,LM, RC) TAUH,TAUM,TAUL,TAUT,TAUTX, & CLDTMP,CLDPRS + ! super-layer RRTMG cloud fraction exports on heartbeat + real, pointer, dimension(:,:) :: CLDTTSWHB + real, pointer, dimension(:,:) :: CLDHISWHB + real, pointer, dimension(:,:) :: CLDMDSWHB + real, pointer, dimension(:,:) :: CLDLOSWHB + + ! locals supporting CLD??SWHB + integer :: rpart, pncol, ncld + real :: plmid(LM), tlev(LM-1), cfac(LM) + integer, allocatable, dimension(:) :: icld, jcld + real, allocatable, dimension(:) :: alat + real, allocatable, dimension(:,:) :: zmid, play, cldfrac, ciwp, clwp + logical, allocatable, dimension(:,:,:) :: cldymcl + real, allocatable, dimension(:,:,:) :: ciwpmcl, clwpmcl + integer, allocatable, dimension(:,:) :: clearCounts + type (ESMF_FieldBundle) :: BUNDLE type (ESMF_Field) :: FIELD type (ESMF_TimeInterval) :: DELT @@ -4557,6 +4618,11 @@ subroutine UPDATE_EXPORT(IM,JM,LM, RC) call MAPL_GetPointer(EXPORT , CLDTMP, 'CLDTMP', __RC__) call MAPL_GetPointer(EXPORT , CLDPRS, 'CLDPRS', __RC__) + call MAPL_GetPointer(EXPORT , CLDLOSWHB, 'CLDLOSWHB', __RC__) + call MAPL_GetPointer(EXPORT , CLDMDSWHB, 'CLDMDSWHB', __RC__) + call MAPL_GetPointer(EXPORT , CLDHISWHB, 'CLDHISWHB', __RC__) + call MAPL_GetPointer(EXPORT , CLDTTSWHB, 'CLDTTSWHB', __RC__) + if (associated(FCLD)) FCLD = CLIN if (associated(CLDH) .or. associated(CLDT) .or. associated(TAUTX)) then @@ -4592,6 +4658,171 @@ subroutine UPDATE_EXPORT(IM,JM,LM, RC) if (associated(CLDT)) CLDT = aCLDT end if + ! CLD??SWHB: + ! Special heartbeat versions of RRTMG generated cloud fractions ... + ! These are expensive because they require a call to the cloud generator, + ! which normally is only done inside the IRRAD and SOLAR REFRESHes (and + ! for the SOLAR case only on the sunlit portion of the globe). We provide + ! these here as a means of validation, but they should not be regularly + ! exported since they will slow down SOLAR considerably, and since the + ! equivalent REFRESH-generated versions (without the "HB" suffix), which + ! are updated only at the REFRESH frequency (~hourly), should be fine in + ! most cases (especially for longer term averages). + ! NB: Filled on all globe unlike the RESFRESH version CLD??SW. + + if (associated(CLDLOSWHB) .or. associated(CLDMDSWHB) .or. & + associated(CLDHISWHB) .or. associated(CLDTTSWHB)) then + + ! default to clear columns which do not need subcolumn generation + if (associated(CLDLOSWHB)) CLDLOSWHB = 0. + if (associated(CLDMDSWHB)) CLDMDSWHB = 0. + if (associated(CLDHISWHB)) CLDHISWHB = 0. + if (associated(CLDTTSWHB)) CLDTTSWHB = 0. + + ! partition size pncol for cloudy columns to conserve memory & improve efficiency + call MAPL_GetResource(MAPL,rpart,'RRTMGSW_PARTITION_SIZE:',DEFAULT=0,__RC__) + if (rpart > 0) then + pncol = rpart + else + pncol = 2 + end if + + ! space for partition: + ! The partition stores up cloudy gridcolumns to process in batch + allocate(icld ( pncol),__STAT__) + allocate(jcld ( pncol),__STAT__) + allocate(zmid (LM, pncol),__STAT__) + allocate(alat ( pncol),__STAT__) + allocate(play (LM, pncol),__STAT__) + allocate(cldfrac(LM, pncol),__STAT__) + allocate(ciwp (LM, pncol),__STAT__) + allocate(clwp (LM, pncol),__STAT__) + allocate(cldymcl(LM,ngptsw,pncol),__STAT__) + allocate(ciwpmcl(LM,ngptsw,pncol),__STAT__) + allocate(clwpmcl(LM,ngptsw,pncol),__STAT__) + allocate(clearCounts(4, pncol),__STAT__) + + ! start with empty partition + ncld = 0 + +! after this test ... make sure DOY used consistently by refresh and update in model +! I guess its only now being used in update, but was its use in refresh really consistent? +! this may be a non-zero-diff bug fix later ... co2 by DOY = hb but DOY for RRTMG should be for REFRESH style time + + ! loop over domain + do j = 1,JM + do i = 1,IM + + ! load up cloudy columns to partition + if (any(CLIN(i,j,:) > 0.) then + + ! cloudy column + ncld = ncld + 1 + icld (ncld) = i + jcld (ncld) = j + alat (ncld) = LATS(i,j) + + ! Note: unlike RRTMG we do not reverse the levels. This is a + ! technicality and will not alter the POPULATION stats of the + ! generation and saves time (see notes under cloud_subcol_gen). + ! If an exact replication of RRTMGSW is required, can reverse + ! vertical ordering here ... but an exact replication will + ! also require saving the exact cldfrac used by the REFRESH + ! into the internal state for use here as well. + + plmid = 0.5 * (PLL(i,j,0:LM-1) + PLL(i,j,1:LM)) + play (:,ncld) = plmid / 100. ! hPa + cldfrac(:,ncld) = CLIN(i,j,:) + + ! cloud water paths converted from g/g to g/m^2 + cfac = 1.02 * 100 * (PLL(i,j,1:LM)-PLL(i,j,0:LM-1)) + ciwp(:,ncld) = cfac * RQI(i,j,:) + clwp(:,ncld) = cfac * RQL(i,j,:) + + ! interior interface temperatures + ! * "0-based" but extema at 0 and LM not needed for zmid; + ! * RRTMG call code uses layer delP (DPR) but any multiple + ! of it, specifically cfac, is equivalent. + + tlev = (T(i,j,1:LM-1) * cfac(2:LM) + T(i,j,2:LM) * cfac(1:LM-1)) & + / ( cfac(2:LM) + cfac(1:LM-1)) + + ! Calculate the LAYER (mid-point) heights. + ! The interlayer distances are needed for the calculations + ! of inter-layer correlation for cloud overlapping. Only + ! *relative* distances matter, so wolog set zmid(LM) = 0. + zmid(LM,ncld) = 0. + do k = LM-1, 1, -1 + ! dz ~ RT/g x dp/p by hysrostatic eqn and ideal gas eqn. + ! The jump from LAYER k+1 to k is centered on LEVEL k + ! since the LEVEL indices are zero-based + zmid(k,ncld) = zmid(k+1,ncld) + MAPL_RGAS * tlev(k) / MAPL_GRAV & + * (plmid(k+1) - plmid(k)) / PLL(i,j,k) + end do + + end if ! cloudy column + + ! nothing to process yet? + if (ncld == 0) cycle + + ! process the partition? + if (ncld == pncol & ! partition is full so process + .or. i == IM .and. j == JM & ! partition is partially full + & but no more gridcolumns + & so process what have. + ) then + + ! McICA subcolumn generation + call generate_stochastic_clouds( & + pncol, ncld, ngptsw, LM, & + zmid, alat, doy, & + play, cldfrac, ciwp, clwp, 1.e-20, & + cldymcl, ciwpmcl, clwpmcl, & + seed_order=[4,3,2,1]) + + ! for super-layer cloud fractions + call clearCounts_threeBand( & + pncol, ncld, ngptsw, LM, LCLDLM, LCLDMH, cldymcl, & + clearCounts) + + ! convert super-layer clearCounts to cloud fractions + if (associated(CLDTTSWHB)) then + do n = 1,ncld + CLDTTSWHB(icld(n),jcld(n)) = 1. - clearCounts(n,1)/float(ngptsw) + end do + end if + if (associated(CLDHISWHB)) then + do n = 1,ncld + CLDHISWHB(icld(n),jcld(n)) = 1. - clearCounts(n,2)/float(ngptsw) + end do + end if + if (associated(CLDMDSWHB)) then + do n = 1,ncld + CLDMDSWHB(icld(n),jcld(n)) = 1. - clearCounts(n,3)/float(ngptsw) + end do + end if + if (associated(CLDLOSWHB)) then + do n = 1,ncld + CLDLOSWHB(icld(n),jcld(n)) = 1. - clearCounts(n,4)/float(ngptsw) + end do + end if + + ! restart partition + ncld = 0 + + end if ! process partition + + end do ! i + end do ! j + + ! clean up + deallocate(icld,jcld,__STAT__) + deallocate(zmid,alat,play,cldfrac,ciwp,clwp,__STAT__) + deallocate(cldymcl,ciwpmcl,clwpmcl,__STAT__) + deallocate(clearCounts,__STAT__) + + end if + if (associated(TAUI) .or. associated(TAUW) .or. associated(TAUR) .or. associated(TAUS).or. & associated(TAUL) .or. associated(TAUM) .or. associated(TAUH) .or. & associated(TAUT) .or. associated(TAUTX) .or. & From 9c8fc599744a8922847b4b74c75717764a41d860 Mon Sep 17 00:00:00 2001 From: Peter Norris Date: Sun, 16 Jan 2022 17:45:46 -0500 Subject: [PATCH 5/9] pmn: minor bugfix --- .../GEOSsolar_GridComp/GEOS_SolarGridComp.F90 | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 index bccbbab9a..ea7f33d45 100644 --- a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 +++ b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 @@ -4714,7 +4714,7 @@ subroutine UPDATE_EXPORT(IM,JM,LM, RC) do i = 1,IM ! load up cloudy columns to partition - if (any(CLIN(i,j,:) > 0.) then + if (any(CLIN(i,j,:) > 0.)) then ! cloudy column ncld = ncld + 1 @@ -4812,8 +4812,8 @@ subroutine UPDATE_EXPORT(IM,JM,LM, RC) end if ! process partition - end do ! i - end do ! j + end do ! i + end do ! j ! clean up deallocate(icld,jcld,__STAT__) @@ -4821,7 +4821,7 @@ subroutine UPDATE_EXPORT(IM,JM,LM, RC) deallocate(cldymcl,ciwpmcl,clwpmcl,__STAT__) deallocate(clearCounts,__STAT__) - end if + end if ! CLD??SWHB if (associated(TAUI) .or. associated(TAUW) .or. associated(TAUR) .or. associated(TAUS).or. & associated(TAUL) .or. associated(TAUM) .or. associated(TAUH) .or. & From 35c3f7ac31ed4a155b300446024e862e8ada2ec4 Mon Sep 17 00:00:00 2001 From: Peter Norris Date: Sun, 16 Jan 2022 17:58:12 -0500 Subject: [PATCH 6/9] pmn: minor bugfix --- .../GEOSsolar_GridComp/GEOS_SolarGridComp.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 index ea7f33d45..ff7ef4990 100644 --- a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 +++ b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 @@ -4768,8 +4768,8 @@ subroutine UPDATE_EXPORT(IM,JM,LM, RC) ! process the partition? if (ncld == pncol & ! partition is full so process .or. i == IM .and. j == JM & ! partition is partially full - & but no more gridcolumns - & so process what have. + & ! but no more gridcolumns + & ! so process what have. ) then ! McICA subcolumn generation From 4272a40fd025d99fbe375191507a13278a14c397 Mon Sep 17 00:00:00 2001 From: Peter Norris Date: Sun, 16 Jan 2022 18:11:43 -0500 Subject: [PATCH 7/9] pmn: minor bugfix --- .../GEOSsolar_GridComp/GEOS_SolarGridComp.F90 | 10 ++++------ 1 file changed, 4 insertions(+), 6 deletions(-) diff --git a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 index ff7ef4990..a3711fa22 100644 --- a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 +++ b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 @@ -4765,12 +4765,10 @@ subroutine UPDATE_EXPORT(IM,JM,LM, RC) ! nothing to process yet? if (ncld == 0) cycle - ! process the partition? - if (ncld == pncol & ! partition is full so process - .or. i == IM .and. j == JM & ! partition is partially full - & ! but no more gridcolumns - & ! so process what have. - ) then + ! process the partition if its full or if its partially + ! full but there are no more gridcolumns left. + + if (ncld == pncol .or. i == IM .and. j == JM) then ! McICA subcolumn generation call generate_stochastic_clouds( & From abc099902117eb3dc3b20b5de490a14949c9f66f Mon Sep 17 00:00:00 2001 From: Peter Norris Date: Mon, 17 Jan 2022 08:05:29 -0500 Subject: [PATCH 8/9] pmn: fixed clearCounts miss-indexing bug for CLDxxSWHB --- .../GEOSsolar_GridComp/GEOS_SolarGridComp.F90 | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 index a3711fa22..c961ff972 100644 --- a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 +++ b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 @@ -4786,22 +4786,22 @@ subroutine UPDATE_EXPORT(IM,JM,LM, RC) ! convert super-layer clearCounts to cloud fractions if (associated(CLDTTSWHB)) then do n = 1,ncld - CLDTTSWHB(icld(n),jcld(n)) = 1. - clearCounts(n,1)/float(ngptsw) + CLDTTSWHB(icld(n),jcld(n)) = 1. - clearCounts(1,n)/float(ngptsw) end do end if if (associated(CLDHISWHB)) then do n = 1,ncld - CLDHISWHB(icld(n),jcld(n)) = 1. - clearCounts(n,2)/float(ngptsw) + CLDHISWHB(icld(n),jcld(n)) = 1. - clearCounts(2,n)/float(ngptsw) end do end if if (associated(CLDMDSWHB)) then do n = 1,ncld - CLDMDSWHB(icld(n),jcld(n)) = 1. - clearCounts(n,3)/float(ngptsw) + CLDMDSWHB(icld(n),jcld(n)) = 1. - clearCounts(3,n)/float(ngptsw) end do end if if (associated(CLDLOSWHB)) then do n = 1,ncld - CLDLOSWHB(icld(n),jcld(n)) = 1. - clearCounts(n,4)/float(ngptsw) + CLDLOSWHB(icld(n),jcld(n)) = 1. - clearCounts(4,n)/float(ngptsw) end do end if From bcb7eecc69cccf85fb2fae920d9cc7cfeaf2c957 Mon Sep 17 00:00:00 2001 From: Peter Norris Date: Wed, 19 Jan 2022 15:35:59 -0500 Subject: [PATCH 9/9] pmn: include new diagnostic usage comments and warnings --- .../GEOSirrad_GridComp/GEOS_IrradGridComp.F90 | 8 +++++ .../GEOSsolar_GridComp/GEOS_SolarGridComp.F90 | 34 +++++++++++++++++++ 2 files changed, 42 insertions(+) diff --git a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSirrad_GridComp/GEOS_IrradGridComp.F90 b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSirrad_GridComp/GEOS_IrradGridComp.F90 index 6b8d83ce0..801140179 100644 --- a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSirrad_GridComp/GEOS_IrradGridComp.F90 +++ b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSirrad_GridComp/GEOS_IrradGridComp.F90 @@ -783,6 +783,14 @@ subroutine SetServices ( GC, RC ) DIMS = MAPL_DimsHorzOnly, & VLOCATION = MAPL_VLocationNone, __RC__ ) +! Note: the four CLDxxLW diagnostics below represent super-layer cloud +! fractions based on the subcolumn cloud generation called in RRTMG LW. +! They are global fields but generated only at the LW REFRESH frequency, +! NOT at the heartbeat. As such, they are useful for diagnostic comparisons +! with CLDTT above and with the full CLDxx set from the Solar GC. But they +! should NOT be used to subsample fields that are produced on the model +! heartbeat (e.g. subsampling for cloud presence). + call MAPL_AddExportSpec(GC, & SHORT_NAME = 'CLDTTLW', & LONG_NAME = 'total_cloud_area_fraction_rrtmg_lw_REFRESH', & diff --git a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 index c961ff972..def6e6b7c 100644 --- a/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 +++ b/GEOSagcm_GridComp/GEOSphysics_GridComp/GEOSradiation_GridComp/GEOSsolar_GridComp/GEOS_SolarGridComp.F90 @@ -938,6 +938,20 @@ subroutine SetServices ( GC, RC ) DIMS = MAPL_DimsHorzOnly, & VLOCATION = MAPL_VLocationNone, __RC__) +! Note: the four CLDxxSW diagnostics below represent super-layer cloud +! fractions based on the subcolumn cloud generation called in RRTMG SW. +! They are sunlit only fields and generated only at the SW REFRESH frequency, +! NOT at the heartbeat. As such, they are useful for diagnostic comparisons +! with with the full CLDxx set above. But they should NOT be used to subsample +! fields that are produced on the model heartbeat (e.g. subsampling for cloud +! presence). Note, also, that when comparing CLDxxSW with CLDxx, it is better +! to subsample both with COSZSW >= cmin, (e.g., 0.25). This COSZSW is a +! REFRESH-frequency version of MCOSZ and, as such, is most appropriate for +! subsampling REFRESH-frequency fields like CLDxxSW. Of course, you can also +! subsample CLDxx with COSZSW since CLDxx are global. By sampling both +! CLDxxSW and CLDxx with COSZSW you get a fair apples-to-apples comparison +! between the two. + call MAPL_AddExportSpec(GC, & SHORT_NAME = 'CLDTTSW', & LONG_NAME = 'total_cloud_area_fraction_rrtmg_sw_REFRESH', & @@ -966,6 +980,22 @@ subroutine SetServices ( GC, RC ) DIMS = MAPL_DimsHorzOnly, & VLOCATION = MAPL_VLocationNone, __RC__ ) +! Note: the four CLDxxSWHB diagnostics below represent super-layer cloud +! fractions based on essentially the same subcolumn cloud generation used +! by RRTMG SW but called from within the SOLAR UPDATE at the HEARTBEAT. +! They are GLOBAL (not just sunlit) fields and generated on the heartbeat. +! But because subcolumn cloud generation is EXPENSIVE, asking for any of +! these exports will DOUBLE the cost of running the SOLAR GC. As such, +! they are for SPECIAL VALIDATION PURPOSES ONLY. No cost is incurred if +! they are not exported. Note, also, that they are NOT EXACTLY heartbeat +! versions of CLDxxSW, since they sample the heartbeat cloud fractions +! not the less frequent snapshots used at REFRESH-frequency, and also since +! the generation inside UPDATE is on non-flipped vertical fields. This latter +! difference should be statistically insignificant. A re-coding to use vert- +! ically flipped fields as per RRTMG SW is possible but will be slightly +! slower, and was deemed not necessary since the first cloud fraction +! frequency difference will likely dominate. + call MAPL_AddExportSpec(GC, & SHORT_NAME = 'CLDTTSWHB', & LONG_NAME = 'total_cloud_area_fraction_rrtmg_sw_HEARTBEAT', & @@ -1057,6 +1087,10 @@ subroutine SetServices ( GC, RC ) DIMS = MAPL_DimsHorzVert, & VLOCATION = MAPL_VLocationCenter, __RC__) +! Note: The following four TAUxxPAR are REFRESH-frequency fields. As such, all +! the important provisos given in the comment on CLDxxSW above apply to these +! fields as well. Please read those provisos. + call MAPL_AddExportSpec(GC, & LONG_NAME = 'in_cloud_optical_thickness_of_low_clouds_RRTMG_PAR_REFRESH', & UNITS = '1' , &