SURFEX v8.1
General documentation of Surfex
isba_lwnet_meb.F90
Go to the documentation of this file.
1 !SFX_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
2 !SFX_LIC This is part of the SURFEX software governed by the CeCILL-C licence
3 !SFX_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
4 !SFX_LIC for details. version 1.
5 ! ##########################################################################
6  SUBROUTINE isba_lwnet_meb(PLAI,PPSN,PPSNA,PEMIS_N,PEMIS_F,PFF, &
7  PTV,PTG,PTN,PLW_RAD,PLWNET_N,PLWNET_V,PLWNET_G, &
8  PLWNET_V_DTV,PLWNET_V_DTG,PLWNET_V_DTN, &
9  PLWNET_G_DTV,PLWNET_G_DTG,PLWNET_G_DTN, &
10  PLWNET_N_DTV,PLWNET_N_DTG,PLWNET_N_DTN, &
11  PSIGMA_F,PSIGMA_FN, &
12  PLWDOWN_GN )
13 
14 !
15 !!**** *ISBA_LWNET_MEB*
16 !!
17 !! PURPOSE
18 !! -------
19 !
20 ! Calculates the net longwave radiation budget terms for fully
21 ! coupled snow, soil-understory vegetation and canopy vegetation.
22 ! Flux derrivatives also herein.
23 !
24 !
25 !!** METHOD
26 !! ------
27 !
28 !! EXTERNAL
29 !! --------
30 !!
31 !! none
32 !!
33 !! IMPLICIT ARGUMENTS
34 !! ------------------
35 !!
36 !!
37 !! REFERENCE
38 !! ---------
39 !!
40 !! Noilhan and Planton (1989)
41 !! Belair (1995)
42 !! * to be done * (2011)
43 !!
44 !! AUTHOR
45 !! ------
46 !!
47 !! A. Boone * Meteo-France *
48 !! P. Samuelsson * SMHI *
49 !! S. Gollvik * SMHI *
50 !!
51 !! MODIFICATIONS
52 !! -------------
53 !! Original 22/01/11
54 !! A. Boone 04/03/17 add factor to dLWnet_n/dT_n term to cause it to vanish
55 !! under some extreme circumstances.
56 !!
57 !-------------------------------------------------------------------------------
58 !
59 !* 0. DECLARATIONS
60 ! ------------
61 !
62 USE modd_isba_par, ONLY : xemissoil, xemisveg
63 !
64 USE mode_meb, ONLY : meb_shield_factor
65 !
66 USE yomhook ,ONLY : lhook, dr_hook
67 USE parkind1 ,ONLY : jprb
68 !
69 IMPLICIT NONE
70 !
71 !* 0.1 declarations of arguments
72 !
73 REAL, DIMENSION(:), INTENT(IN) :: PLAI, PPSN, PPSNA, PLW_RAD, PSIGMA_F
74 !
75 REAL, DIMENSION(:), INTENT(IN) :: PTV, PTG, PTN
76 !
77 REAL, DIMENSION(:), INTENT(IN) :: PEMIS_N, PEMIS_F, PFF
78 !
79 REAL, DIMENSION(:), INTENT(OUT) :: PLWNET_N, PLWNET_V, PLWNET_G
80 !
81 REAL, DIMENSION(:), INTENT(OUT) :: PLWDOWN_GN
82 !
83 REAL, DIMENSION(:), INTENT(OUT) :: PLWNET_V_DTV, PLWNET_V_DTG, PLWNET_V_DTN
84 ! PLWNET_V_DTV, PLWNET_V_DTG, PLWNET_V_DTN = Vegetation canopy net radiation
85 ! derivatives w/r/t surface temperature(s) (W m-2 K-1)
86 !
87 REAL, DIMENSION(:), INTENT(OUT) :: PLWNET_G_DTV, PLWNET_G_DTG, PLWNET_G_DTN
88 ! PLWNET_G_DTV, PLWNET_G_DTG, PLWNET_G_DTN = Understory-ground net radiation
89 ! derivatives w/r/t surface temperature(s) (W m-2 K-1)
90 !
91 REAL, DIMENSION(:), INTENT(OUT) :: PLWNET_N_DTV, PLWNET_N_DTG, PLWNET_N_DTN
92 ! PLWNET_N_DTV, PLWNET_N_DTG, PLWNET_N_DTN = Ground-based snow net radiation
93 ! derivatives w/r/t surface temperature(s) (W m-2 K-1)
94 !
95 REAL, DIMENSION(:), INTENT(OUT) :: PSIGMA_FN
96 !
97 !* 0.2 declarations of local variables
98 !
99 !
100 REAL, DIMENSION(SIZE(PLAI)) :: ZLWUP
101 !
102 REAL, DIMENSION(SIZE(PLAI)) :: ZSIGMA_FA, ZPN, ZFRAC, ZEMIS
103 !
104 REAL, DIMENSION(SIZE(PLAI)) :: ZLW_G_A, ZLW_G_B, ZLW_G_C, &
105  ZLW_G_D, ZLW_G_E, ZLW_G_F, ZLW_G_G, &
106  ZLW_G_H, ZLW_G_I, ZLW_G_J, ZLW_G_K, &
107  ZLW_G_L
108 !
109 REAL, DIMENSION(SIZE(PLAI)) :: ZLW_N_A, ZLW_N_B, ZLW_N_C, &
110  ZLW_N_D, ZLW_N_E, ZLW_N_F, ZLW_N_G, &
111  ZLW_N_H, ZLW_N_I, ZLW_N_J, ZLW_N_K, &
112  ZLW_N_L
113 !
114 REAL, DIMENSION(SIZE(PEMIS_N)) :: ZEMIS_G
115 !
116 !* 0.3 declarations of local parameters
117 !
118 REAL, PARAMETER :: ZTGRAD_MAX = 60. ! K When Tn-Tv approaches this difference, the derrivative term
119  ! dLWnet_n/dTn vanishes. See below for details.
120 REAL, PARAMETER :: ZTGRAD_DIF = 10. ! K This is the range in Tn-Tv over which dLWnet_n/dTn
121  ! linearly vanishes. See below for details.
122 !
123 REAL(KIND=JPRB) :: ZHOOK_HANDLE
124 !-------------------------------------------------------------------------------
125 !
126 !* 0. Initialization:
127 ! ---------------
128 !
129 IF (lhook) CALL dr_hook('ISBA_LWNET_MEB',0,zhook_handle)
130 !-------------------------------------------------------------------------------
131 !
132 ! Soil with flooded part:
133 !
134 zemis_g(:) = xemissoil*(1.-pff(:)) + pff(:)*pemis_f(:)
135 !
136 !* 1. View factors: transmission
137 ! --------------------------
138 !
139 psigma_fn(:) = 1.0 - meb_shield_factor(plai,ppsna) ! NOTE: Effects of intercepted snow on the
140 ! ! canopy is neglected.
141 !
142 !* 2. Longwave radiation terms
143 ! ------------------------
144 
145 ! - over snow-free fraction:
146 
147 zfrac(:) = 1.-ppsn(:)
148 zpn(:) = ppsn(:)*(1.-ppsna(:))
149 zsigma_fa(: ) = (1.-zpn(:))*psigma_f(:) + zpn(:)*psigma_fn(:)
150 
151  CALL lw_flux_comp(zpn,plw_rad,zfrac,psigma_f,zsigma_fa, &
152  zemis_g,ptv,ptg, &
153  zlw_g_a,zlw_g_b,zlw_g_c,zlw_g_d,zlw_g_e,zlw_g_f, &
154  zlw_g_g,zlw_g_h,zlw_g_i,zlw_g_j,zlw_g_k,zlw_g_l )
155 
156 ! - over snow-covered fraction:
157 
158 zfrac(:) = ppsn(:)
159 zpn(:) = ppsn(:) + ppsna(:)*(1.-ppsn(:))
160 zsigma_fa(: ) = (1.-zpn(:))*psigma_f(:) + zpn(:)*psigma_fn(:)
161 
162  CALL lw_flux_comp(zpn,plw_rad,zfrac,psigma_fn,zsigma_fa, &
163  pemis_n,ptv,ptn, &
164  zlw_n_a,zlw_n_b,zlw_n_c,zlw_n_d,zlw_n_e,zlw_n_f, &
165  zlw_n_g,zlw_n_h,zlw_n_i,zlw_n_j,zlw_n_k,zlw_n_l )
166 
167 !------------------------------------------------------------------
168 ! Diagnostics
169 !------------------------------------------------------------------
170 
171 ! Total LW energy flux reaching ground/snow surface: (W m-2)
172 ! (explicit part: the implicit flux needs to have the derrivative terms added)
173 
174 plwdown_gn(:) = zlw_g_c(:) + zlw_g_f(:) + zlw_n_c(:) + zlw_n_f(:) + &
175  zlw_g_j(:) + zlw_g_k(:) + zlw_n_j(:) + zlw_n_k(:)
176 
177 !------------------------------------------------------------------
178 ! - compute derivatives: W m-2 K-1
179 !------------------------------------------------------------------
180 
181 plwnet_v_dtv(:) = ( zlw_g_g(:) - zlw_g_h(:) - 2*zlw_g_f(:) &
182  + zlw_n_g(:) - zlw_n_h(:) - 2*zlw_n_f(:) )*4/ptv(:)
183 plwnet_v_dtg(:) = ( zlw_g_i(:) - zlw_g_j(:) - zlw_g_k(:) &
184  - zlw_g_l(:) )*4/ptg(:)
185 plwnet_v_dtn(:) = ( zlw_n_i(:) - zlw_n_j(:) - zlw_n_k(:) &
186  - zlw_n_l(:) )*4/ptn(:)
187 
188 plwnet_g_dtv(:) = ( zlw_g_f(:) - zlw_g_g(:) )*4/ptv(:)
189 plwnet_g_dtg(:) = ( zlw_g_j(:) - zlw_g_i(:) )*4/ptg(:)
190 plwnet_g_dtn(:) = zlw_n_j(:) *4/ptn(:)
191 
192 plwnet_n_dtv(:) = ( zlw_n_f(:) - zlw_n_g(:) )*4/ptv(:)
193 plwnet_n_dtg(:) = zlw_g_k(:) *4/ptg(:)
194 plwnet_n_dtn(:) = ( zlw_n_k(:) - zlw_n_i(:) )*4/ptn(:)
195 !
196 ! Note that for very thin snow combined with extremely cold
197 ! conditions, very weak LWdown forcing and strong Tn-Tc differences,
198 ! the derrivative of Lwnet_n/dT_n can pose problems (i.e. the assumption of
199 ! a linear change over the given time step is a poor approximation),
200 ! thus in these rare instances, this derrivative is forced
201 ! to vanish (over some small temperature difference range). This has no impact
202 ! on numerical stability (since the LW linerization has no effect on this) and has
203 ! no impact on results in the general sense (since it is temporary and the assumed gradieints
204 ! are so large that they are rarely, if ever during a run, encountered).
205 ! Use PTV as a proxy for air T (under these conditions, it should be close to PTC for example)
206 
207 plwnet_n_dtn(:) = plwnet_n_dtn(:)*min(1., max(0., ztgrad_max-max(0.,ptn(:)-ptv(:)))/ztgrad_dif)
208 
209 !------------------------------------------------------------------
210 ! - Compute *explicit* net budgets (at time t): W m-2
211 ! NOTE: fully implicit budgets (at time t+dt) are computed
212 ! *after* energy budget
213 !------------------------------------------------------------------
214 
215 zlwup(:) = zlw_g_b(:) + zlw_g_e(:) + zlw_g_f(:) + zlw_g_h(:) &
216  + zlw_g_l(:) &
217  + zlw_n_b(:) + zlw_n_e(:) + zlw_n_f(:) + zlw_n_h(:) &
218  + zlw_n_l(:)
219 
220 plwnet_g(:) = zlw_g_c(:) + zlw_g_f(:) + zlw_g_j(:) &
221  - zlw_g_i(:) - zlw_g_d(:) - zlw_g_g(:) &
222  + zlw_n_j(:)
223 
224 plwnet_n(:) = zlw_n_c(:) + zlw_n_f(:) + zlw_n_k(:) &
225  - zlw_n_i(:) - zlw_n_d(:) - zlw_n_g(:) &
226  + zlw_g_k(:)
227 
228 plwnet_v(:) = plw_rad(:) - zlwup(:) - plwnet_g(:) - plwnet_n(:)
229 
230 IF (lhook) CALL dr_hook('ISBA_LWNET_MEB',1,zhook_handle)
231 
232 CONTAINS
233 !=========================================================
234 SUBROUTINE lw_flux_comp(PPN,PLW_RAD,PFRAC,PSIGMA_F,PSIGMA_FA, &
235  PEMIS_S,PTV,PTEMP_S, &
236  PLW_A,PLW_B,PLW_C,PLW_D,PLW_E,PLW_F,PLW_G,PLW_H, &
237  PLW_I,PLW_J,PLW_K,PLW_L )
238 
239 USE modd_csts, ONLY : xstefan
240 
241 IMPLICIT NONE
242 
243 REAL, DIMENSION(:), INTENT(IN) :: PPN, PLW_RAD, PSIGMA_F, PSIGMA_FA, PFRAC
244 REAL, DIMENSION(:), INTENT(IN) :: PTEMP_S, PTV
245 REAL, DIMENSION(:), INTENT(IN) :: PEMIS_S
246 REAL, DIMENSION(:), INTENT(OUT) :: PLW_A, PLW_B, PLW_C, PLW_D, PLW_E, PLW_F, &
247  PLW_G, PLW_H, PLW_I, PLW_J, PLW_K, PLW_L
248 
249 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
250 
251 REAL, DIMENSION(SIZE(PLW_RAD)) :: ZWORK
252 REAL(KIND=JPRB) :: ZHOOK_HANDLE
253 !--------------------------------------------------------------
254 IF (lhook) CALL dr_hook('LW_FLUX_COMP',0,zhook_handle)
255 
256 plw_a(:) = plw_rad(:)*pfrac(:)
257 plw_b(:) = plw_a(:)* psigma_f(:) *(1.-xemisveg)
258 plw_c(:) = plw_a(:)*(1.-psigma_f(:))
259 plw_d(:) = plw_c(:) *(1.-pemis_s(:))
260 plw_e(:) = plw_d(:)*(1.-psigma_fa(:))
261 
262 plw_f(:) = psigma_fa(:) * xemisveg * pfrac(:) *xstefan*(ptv(:)**4)
263 plw_g(:) = plw_f(:) *(1.-pemis_s(:))
264 plw_h(:) = plw_g(:)*(1.-psigma_fa(:))
265 
266 zwork(:) = (1.-xemisveg)*psigma_fa(:)
267 plw_i(:) = pemis_s(:) * pfrac(:) *xstefan*(ptemp_s(:)**4)
268 plw_j(:) = plw_i(:)*zwork(:) *(1.-ppn(:))
269 plw_k(:) = plw_i(:)*zwork(:)* ppn(:)
270 plw_l(:) = plw_i(:)*(1.-psigma_fa(:))
271 
272 IF (lhook) CALL dr_hook('LW_FLUX_COMP',1,zhook_handle)
273 
274 END SUBROUTINE lw_flux_comp
275 !=========================================================
276 
277 END SUBROUTINE isba_lwnet_meb
yomhook::dr_hook
Definition: yomhook.F90:20
mode_meb
Definition: mode_meb.F90:6
lw_flux_comp
subroutine lw_flux_comp(PPN, PLW_RAD, PFRAC, PSIGMA_F, PSIGMA_FA, EMIS_S, PTV, PTEMP_S, LW_A, PLW_B, PLW_C, PLW_D, PLW_E, PLW_F, PLW_G, PLW_H, LW_I, PLW_J, PLW_K, PLW_L)
Definition: isba_lwnet_meb.F90:238
modd_csts::xstefan
real, save xstefan
Definition: modd_csts.F90:59
mode_meb::meb_shield_factor
Definition: mode_meb.F90:65
parkind1::jprb
integer, parameter jprb
Definition: parkind1.F90:32
yomhook::lhook
logical lhook
Definition: yomhook.F90:15
modd_csts
Definition: modd_csts.F90:6
isba_lwnet_meb
subroutine isba_lwnet_meb(PLAI, PPSN, PPSNA, PEMIS_N, PEMIS_F, PFF,
Definition: isba_lwnet_meb.F90:7
Generated on Tue Jan 16 2018 16:23:29 for SURFEX v8.1 by   doxygen 1.8.13