!WRF:MODEL_LAYER:PHYSICS
!
SUBROUTINE radiation_driver (itimestep,dt, & 2,28
RTHRATENLW,RTHRATENSW,RTHRATEN, &
GLW,GSW, &
XLAT,XLONG,ALBEDO,CLDFRA,EMISS, &
rho_phy,moist,n_moist, &
p8w,p_phy,Pb,pi_phy,dz8w,t_phy,t8w,GMT, &
!#$ JULDAY,config_flags,RADT,STEPRA,ICLOUD, &
JULDAY,RADT, &
NRADS,NRADL,RA_SW_PHYSICS,RA_LW_PHYSICS, &
ICLOUD, &
taucldi,taucldc,warm_rain, &
! for eta radiation
XLAND,TSK,HTOP,HBOT,CUPPT,VEGFRA,SNOW,LTOP, &
julyr,NPHS, &
RSWTOA,RLWTOA,CZMEAN, & !Added
CFRACL,CFRACM,CFRACH, & !Added
ACFRST,NCFRST,ACFRCV,NCFRCV, & !Added
ids,ide, jds,jde, kds,kde, &
ims,ime, jms,jme, kms,kme, &
its,ite, jts,jte, kts,kte )
!-------------------------------------------------------------------------
!#$USE module_bc
!#$USE module_state_description
USE module_model_constants
!#$USE module_wrf_error
! *** add new modules of schemes here
!#$USE module_ra_sw
!#$USE module_ra_gsfcsw
!#$USE module_ra_rrtm
USE module_ra_gfdleta
! This driver calls subroutines for the radiation parameterizations.
!
! short wave radiation choices:
! 1. swrad (19??)
!
! long wave radiation choices:
! 1. rrtmlwrad
!
!----------------------------------------------------------------------
IMPLICIT NONE
!======================================================================
! Grid structure in physics part of WRF
!----------------------------------------------------------------------
! The horizontal velocities used in the physics are unstaggered
! relative to temperature/moisture variables. All predicted
! variables are carried at half levels except w, which is at full
! levels. Some arrays with names (*8w) are at w (full) levels.
!
!----------------------------------------------------------------------
! In WRF, kms (smallest number) is the bottom level and kme (largest
! number) is the top level. In your scheme, if 1 is at the top level,
! then you have to reverse the order in the k direction.
!
! kme - half level (no data at this level)
! kme ----- full level
! kme-1 - half level
! kme-1 ----- full level
! .
! .
! .
! kms+2 - half level
! kms+2 ----- full level
! kms+1 - half level
! kms+1 ----- full level
! kms - half level
! kms ----- full level
!
!======================================================================
! Grid structure in physics part of WRF
!-------------------------------------
! The horizontal velocities used in the physics are unstaggered
! relative to temperature/moisture variables. All predicted
! variables are carried at half levels except w, which is at full
! levels. Some arrays with names (*8w) are at w (full) levels.
!
!==================================================================
! Definitions
!-----------
! Rho_d dry density (kg/m^3)
! Theta_m moist potential temperature (K)
! Qv water vapor mixing ratio (kg/kg)
! Qc cloud water mixing ratio (kg/kg)
! Qr rain water mixing ratio (kg/kg)
! Qi cloud ice mixing ratio (kg/kg)
! Qs snow mixing ratio (kg/kg)
!-----------------------------------------------------------------
!-- RTHRATEN Rho_dTheta_m tendency
! due to radiation (kg/m^3 . K)
!-- RTHRATENLW Rho_dTheta_m tendency
! due to long wave radiation (kg/m^3 . K)
!-- RTHRATENSW Rho_dTheta_m temperature tendency
! due to short wave radiation (kg/m^3 . K)
!-- dt time step (s)
!-- itimestep number of time steps
!-- GLW downward long wave flux at ground surface (W/m^2)
!-- GSW downward short wave flux at ground surface (W/m^2)
!-- XLAT latitude, south is negative (degree)
!-- XLONG longitude, west is negative (degree)
!-- ALBEDO albedo (between 0 and 1)
!-- CLDFRA cloud fraction (between 0 and 1)
!-- EMISS surface emissivity (between 0 and 1)
!-- rho_phy density (kg/m^3)
!-- rr dry air density (kg/m^3)
!-- moist moisture array (4D - last index is species) (kg/kg)
!-- n_moist number of moisture species
!-- p8w pressure at full levels (Pa)
!-- p_phy pressure (Pa)
!-- Pb base-state pressure (Pa)
!-- pi_phy exner function (dimensionless)
!-- dz8w dz between full levels (m)
!-- t_phy temperature (K)
!-- t8w temperature at full levels (K)
!-- GMT Greenwich Mean Time Hour of model start (hour)
!-- JULDAY the initial day (Julian day)
!-- config_flags boundary condition flag
!-- RADT time for calling radiation (min)
!-- DEGRAD conversion factor for
! degrees to radians (pi/180.) (rad/deg)
!-- DPD degrees per day for earth's
! orbital position (deg/day)
!-- R_d gas constant for dry air (J/kg/K)
!-- CP heat capacity at constant pressure for dry air (J/kg/K)
!-- G acceleration due to gravity (m/s^2)
!-- rvovrd R_v divided by R_d (dimensionless)
!-- XTIME time since simulation start (min)
!-- DECLIN solar declination angle (deg)
!-- SOLCON solar constant (W/m^2)
!-- P_QV species index for water vapor
!-- P_QC species index for cloud water
!-- P_QR species index for rain water
!-- P_QI species index for cloud ice
!-- P_QS species index for snow
!-- P_QG species index for graupel
!-- ids start index for i in domain
!-- ide end index for i in domain
!-- jds start index for j in domain
!-- jde end index for j in domain
!-- kds start index for k in domain
!-- kde end index for k in domain
!-- ims start index for i in memory
!-- ime end index for i in memory
!-- jms start index for j in memory
!-- jme end index for j in memory
!-- kms start index for k in memory
!-- kme end index for k in memory
!-- its start index for i in tile
!-- ite end index for i in tile
!-- jts start index for j in tile
!-- jte end index for j in tile
!-- kts start index for k in tile
!-- kte end index for k in tile
!
!==================================================================
!#$TYPE(grid_config_rec_type), INTENT(IN ) :: config_flags
!
INTEGER, INTENT(IN ) :: ids,ide, jds,jde, kds,kde, &
ims,ime, jms,jme, kms,kme, &
its,ite, jts,jte, kts,kte, &
n_moist
INTEGER, INTENT(IN ) :: NRADS,NRADL,ICLOUD
CHARACTER(20),INTENT(IN) :: RA_LW_PHYSICS,RA_SW_PHYSICS
LOGICAL, INTENT(IN ) :: warm_rain
REAL, INTENT(IN ) :: RADT
REAL, DIMENSION( ims:ime, jms:jme ), &
INTENT(IN ) :: XLAND, &
TSK, &
VEGFRA, &
SNOW
REAL, DIMENSION( ims:ime, jms:jme ), &
INTENT(INOUT) :: HTOP, &
HBOT, &
CUPPT
INTEGER, INTENT(IN ) :: julyr,NPHS
!
REAL, DIMENSION( ims:ime, kms:kme, jms:jme ), &
INTENT(IN ) :: dz8w, &
p8w, &
p_phy, &
Pb, &
pi_phy, &
t_phy, &
t8w, &
rho_phy, &
CLDFRA
!
REAL, DIMENSION( ims:ime, kms:kme, jms:jme, n_moist ), &
INTENT(IN ) :: moist
!
REAL, DIMENSION( ims:ime, kms:kme, jms:jme ), &
INTENT(INOUT) :: RTHRATEN, &
RTHRATENLW, &
RTHRATENSW
REAL, DIMENSION( ims:ime, kms:kme, jms:jme ), &
INTENT(INOUT) :: taucldi,taucldc
!
REAL, DIMENSION( ims:ime, jms:jme ), &
INTENT(IN ) :: XLAT, &
XLONG, &
ALBEDO, &
EMISS
!
REAL, DIMENSION( ims:ime, jms:jme ),INTENT(OUT) :: CFRACH, & !Added
CFRACL, & !Added
CFRACM, & !Added
CZMEAN !Added
!
REAL, DIMENSION( ims:ime, jms:jme ), &
INTENT(INOUT) :: GSW, &
GLW, &
RLWTOA, & !Added
RSWTOA, & !Added
ACFRST, & !Added
ACFRCV !Added
!
REAL, INTENT(IN ) :: GMT,dt
!
INTEGER, INTENT(IN ) :: JULDAY, itimestep
!
INTEGER,DIMENSION(3),INTENT(IN) :: LTOP
!
INTEGER,DIMENSION( ims:ime, jms:jme ),INTENT(INOUT) :: NCFRST, & !Added
NCFRCV !Added
! LOCAL VAR
REAL, DIMENSION( its:ite, jts:jte ) :: GLAT,GLON
REAL :: XTIME,DECLIN,SOLCON
INTEGER :: i,j,k
!#$LOGICAL :: gfdl_lw,gfdl_sw
LOGICAL :: longwave,shortwave
REAL :: OBECL,SINOB,SXLONG,ARG,DECDEG, &
DJUL,RJUL,ECCFAC
INTEGER :: P_QV=1
INTEGER :: P_QC=2
CHARACTER(20),PARAMETER :: RRTMSCHEME='RRTMSCHEME'
CHARACTER(20),PARAMETER :: GFDLLWSCHEME='GFDLLWSCHEME'
CHARACTER(20),PARAMETER :: SWRADSCHEME='SWRADSCHEME'
CHARACTER(20),PARAMETER :: GSFCSWSCHEME='GSFCSWSCHEME'
CHARACTER(20),PARAMETER :: GFDLSWSCHEME='GFDLSWSCHEME'
!------------------------------------------------------------------
!#$if (config_flags%ra_lw_physics .eq. 0 .and. &
!#$ config_flags%ra_sw_physics .eq. 0) return
!#$gfdl_lw = .false.
!#$gfdl_sw = .false.
longwave=.false.
shortwave=.false.
!#@IF (itimestep .eq. 1 .or. mod(itimestep-1,STEPRA) .eq. 0) THEN
!---------------
! initialize data
DO j=jts,jte
DO i=its,ite
GSW(I,J)=0.
GLW(I,J)=0.
GLAT(I,J)=XLAT(I,J)*DEGRAD
GLON(I,J)=XLONG(I,J)*DEGRAD
ENDDO
ENDDO
DO j=jts,jte
DO k=kts,kte
DO i=its,ite
RTHRATEN(I,K,J)=0.
ENDDO
ENDDO
ENDDO
!---------------
! Calculate constant for short wave radiation
CALL radconst
(XTIME,DECLIN,SOLCON,GMT,JULDAY, &
DEGRAD,DPD,itimestep,dt )
! CALL cal_cldfra(CLDFRA,moist(ims,kms,jms,P_QC), &
! moist(ims,kms,jms,P_QI),P_QI, &
! ids,ide, jds,jde, kds,kde, &
! ims,ime, jms,jme, kms,kme, &
! its,ite, jts,jte, kts,kte )
!#$WRITE(wrf_err_message,*)'SOLCON=',SOLCON,DECLIN,XTIME
!#$CALL wrf_debug(50,wrf_err_message)
!
IF(mod(itimestep-1,NRADL).eq.0)THEN
longwave=.true.
ENDIF
!
IF(LONGWAVE)THEN
!#$lwrad_select: SELECT CASE(config_flags%ra_lw_physics)
lwrad_select: SELECT CASE(ra_lw_physics)
CASE (RRTMSCHEME)
!#$ CALL wrf_debug (100, 'CALL rrtm')
! CALL RRTMLWRAD(RTHRATEN,GLW,EMISS, &
! moist(ims,kms,jms,P_QV), &
! moist(ims,kms,jms,P_QC), &
! moist(ims,kms,jms,P_QR), &
! moist(ims,kms,jms,P_QI), &
! moist(ims,kms,jms,P_QS), &
! moist(ims,kms,jms,P_QG), &
! p8w,p_phy,pi_phy,dz8w,t_phy, &
! t8w, rho_phy, CLDFRA, &
! R_d,G,P_QI,P_QS,P_QG,ICLOUD, &
! PARAM_FIRST_SCALAR,warm_rain, &
! ids,ide, jds,jde, kds,kde, &
! ims,ime, jms,jme, kms,kme, &
! its,ite, jts,jte, kts,kte )
CASE (GFDLLWSCHEME)
!#$ CALL wrf_debug (100, 'CALL gfdllw')
!#$ gfdl_lw = .true.
CALL ETARA
(DT,RTHRATEN,RTHRATENLW,RTHRATENSW,pi_phy, &
XLAND,p8w,dz8w,rho_phy,p_phy,t_phy, &
moist(ims,kms,jms,P_QV), &
moist(ims,kms,jms,P_QC),TSK,GLW,GSW, &
RLWTOA,RSWTOA,CZMEAN, & !Added
GLAT,GLON,HTOP,HBOT,ALBEDO,CUPPT, &
VEGFRA,SNOW,G,GMT, & !Modified
SHORTWAVE,LONGWAVE,NRADS,NRADL,NPHS,itimestep, & !Modified
!#$ julyr,julday,gfdl_lw,gfdl_sw, &
julyr,julday,LTOP, & !Modified
CFRACL,CFRACM,CFRACH, &
ACFRST,NCFRST,ACFRCV,NCFRCV, &
ids,ide, jds,jde, kds,kde, &
ims,ime, jms,jme, kms,kme, &
its,ite, jts,jte, kts,kte )
CASE DEFAULT
!#$ CALL wrf_message ('The long wave scheme does not exist')
END SELECT lwrad_select
!#$IF (config_flags%ra_lw_physics .gt. 0 ) THEN
DO j=jts,jte
DO k=kts,kte
DO i=its,ite
RTHRATENLW(I,K,J)=RTHRATEN(I,K,J)
ENDDO
ENDDO
ENDDO
!#$ENDIF
longwave=.false.
ENDIF !End LONGWAVE block
!
IF(mod(itimestep-1,NRADS).eq.0)THEN
shortwave=.true.
ENDIF
!
IF(SHORTWAVE)THEN
!#$swrad_select: SELECT CASE(config_flags%ra_sw_physics)
swrad_select: SELECT CASE(ra_sw_physics)
CASE (SWRADSCHEME)
!#$ CALL wrf_debug(100, 'CALL swrad')
! CALL SWRAD(dt,RTHRATEN,GSW,XLAT,XLONG,ALBEDO, &
! rho_phy,t_phy, &
! moist(ims,kms,jms,P_QV), &
! moist(ims,kms,jms,P_QC), &
! moist(ims,kms,jms,P_QR), &
! moist(ims,kms,jms,P_QI), &
! moist(ims,kms,jms,P_QS), &
! moist(ims,kms,jms,P_QG), &
! p_phy,pi_phy,dz8w,GMT, &
! R_d,CP,G,JULDAY, &
! XTIME,DECLIN,SOLCON, &
! P_QI,P_QS,P_QG,RADT,ICLOUD,DEGRAD, &
! PARAM_FIRST_SCALAR,warm_rain, &
! ids,ide, jds,jde, kds,kde, &
! ims,ime, jms,jme, kms,kme, &
! its,ite, jts,jte, kts,kte )
CASE (GSFCSWSCHEME)
! CALL GSFCSWRAD(RTHRATEN,GSW,XLAT,XLONG, &
! ALBEDO,t_phy, &
! moist(ims,kms,jms,P_QV), &
! moist(ims,kms,jms,P_QC), &
! moist(ims,kms,jms,P_QR), &
! moist(ims,kms,jms,P_QI), &
! moist(ims,kms,jms,P_QS), &
! moist(ims,kms,jms,P_QG), &
! p_phy,p8w,pi_phy,CLDFRA, &
! GMT,CP,G,JULDAY,XTIME,DECLIN,SOLCON, &
! P_QI,P_QS,P_QG,RADT,DEGRAD, &
! taucldi,taucldc,PARAM_FIRST_SCALAR,warm_rain, &
! ids,ide, jds,jde, kds,kde, &
! ims,ime, jms,jme, kms,kme, &
! its,ite, jts,jte, kts,kte )
CASE (GFDLSWSCHEME)
!#$ CALL wrf_debug (100, 'CALL gfdlsw')
!#$ gfdl_sw = .true.
CALL ETARA
(DT,RTHRATEN,RTHRATENLW,RTHRATENSW,pi_phy, &
XLAND,p8w,dz8w,rho_phy,p_phy,t_phy, &
moist(ims,kms,jms,P_QV), &
moist(ims,kms,jms,P_QC),TSK,GLW,GSW, &
RLWTOA,RSWTOA,CZMEAN, & !Added
GLAT,GLON,HTOP,HBOT,ALBEDO,CUPPT, &
VEGFRA,SNOW,G,GMT, & !Modified
SHORTWAVE,LONGWAVE,NRADS,NRADL,NPHS,itimestep, & !Modified
!#$ julyr,julday,gfdl_lw,gfdl_sw, &
julyr,julday,LTOP, & !Modified
CFRACL,CFRACM,CFRACH, &
ACFRST,NCFRST,ACFRCV,NCFRCV, &
ids,ide, jds,jde, kds,kde, &
ims,ime, jms,jme, kms,kme, &
its,ite, jts,jte, kts,kte )
CASE DEFAULT
!#$ CALL wrf_message('The short wave scheme does not exist')
END SELECT swrad_select
!#$IF (config_flags%ra_sw_physics .gt. 0) THEN
DO j=jts,jte
DO k=kts,kte
DO i=its,ite
RTHRATENSW(I,K,J)=RTHRATEN(I,K,J)-RTHRATENLW(I,K,J)
ENDDO
ENDDO
ENDDO
!#$ENDIF
shortwave=.false.
ENDIF !End SHORTWAVE block
END SUBROUTINE radiation_driver
!---------------------------------------------------------------------
SUBROUTINE radconst(XTIME,DECLIN,SOLCON,GMT,JULDAY, & 2,2
DEGRAD,DPD,step,dt )
!---------------------------------------------------------------------
!#$USE module_wrf_error
IMPLICIT NONE
!---------------------------------------------------------------------
INTEGER, INTENT(IN ) :: JULDAY, step
REAL, INTENT(IN ) :: GMT,dt,DEGRAD,DPD
REAL, INTENT(OUT ) :: XTIME,DECLIN,SOLCON
REAL :: OBECL,SINOB,SXLONG,ARG,JULIAN, &
DECDEG,DJUL,RJUL,ECCFAC
! for short wave radiation
DECLIN=0.
SOLCON=0.
!-----OBECL : OBLIQUITY = 23.5 DEGREE.
OBECL=23.5*DEGRAD
SINOB=SIN(OBECL)
XTIME=float(step)*dt/60.
!-----CALCULATE LONGITUDE OF THE SUN FROM VERNAL EQUINOX:
JULIAN=FLOAT(JULDAY-1)+(XTIME/60.+GMT)/24.
IF(JULIAN.GE.80.)SXLONG=DPD*(JULIAN-80.)
IF(JULIAN.LT.80.)SXLONG=DPD*(JULIAN+285.)
SXLONG=SXLONG*DEGRAD
ARG=SINOB*SIN(SXLONG)
DECLIN=ASIN(ARG)
DECDEG=DECLIN/DEGRAD
!----SOLAR CONSTANT ECCENTRICITY FACTOR (PALTRIDGE AND PLATT 1976)
DJUL=JULIAN*360./365.
RJUL=DJUL*DEGRAD
ECCFAC=1.000110+0.034221*COS(RJUL)+0.001280*SIN(RJUL)+0.000719* &
COS(2*RJUL)+0.000077*SIN(2*RJUL)
SOLCON=1370.*ECCFAC
!#$write(wrf_err_message,10)DECDEG,SOLCON
10 FORMAT(1X,'*** SOLAR DECLINATION ANGLE = ',F6.2,' DEGREES.', &
' SOLAR CONSTANT = ',F8.2,' W/M**2 ***')
!#$CALL wrf_debug (50, wrf_err_message)
END SUBROUTINE radconst
!---------------------------------------------------------------------
SUBROUTINE cal_cldfra(CLDFRA,QC,QI,P_QI, & 1
ids,ide, jds,jde, kds,kde, &
ims,ime, jms,jme, kms,kme, &
its,ite, jts,jte, kts,kte )
!---------------------------------------------------------------------
IMPLICIT NONE
!---------------------------------------------------------------------
INTEGER, INTENT(IN ) :: ids,ide, jds,jde, kds,kde, &
ims,ime, jms,jme, kms,kme, &
its,ite, jts,jte, kts,kte
INTEGER, INTENT(IN ) :: P_QI
!
REAL, DIMENSION( ims:ime, kms:kme, jms:jme ), INTENT(OUT ) :: &
CLDFRA
REAL, DIMENSION( ims:ime, kms:kme, jms:jme ), INTENT(IN ) :: &
QI, &
QC
REAL thresh
INTEGER:: i,j,k
!---------------------------------------------------------------------
thresh=1.0e-6
IF (P_QI .gt. 0) THEN
DO j = jts,jte
DO k = kts,kte
DO i = its,ite
IF ( QC(i,k,j)+QI(I,k,j) .gt. thresh) THEN
CLDFRA(i,k,j)=1.
ELSE
CLDFRA(i,k,j)=0.
ENDIF
ENDDO
ENDDO
ENDDO
ELSE
DO j = jts,jte
DO k = kts,kte
DO i = its,ite
IF ( QC(i,k,j) .gt. thresh) THEN
CLDFRA(i,k,j)=1.
ELSE
CLDFRA(i,k,j)=0.
ENDIF
ENDDO
ENDDO
ENDDO
ENDIF
END SUBROUTINE cal_cldfra