SURFEX v8.1
General documentation of Surfex
bem.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 bem(BOP, T, B, DMT, PTSTEP, PSUNTIME, KDAY, PPS, PRHOA, PT_CAN, &
7  PQ_CAN, PU_CAN, PHU_BLD, PT_RAD_IND, PFLX_BLD_FL, PFLX_BLD_MA,&
8  PRADHT_IN, PRAD_RF_MA, PRAD_RF_FL, PRAD_WL_MA, PRAD_WL_FL,&
9  PRAD_WIN_MA, PRAD_WIN_FL, PCONV_RF_BLD, PCONV_WL_BLD, &
10  PCONV_WIN_BLD, PLOAD_IN_FL, PLOAD_IN_MA )
11 ! ##########################################################################
12 !
13 !!**** *BEM*
14 !!
15 !! PURPOSE
16 !! -------
17 !
18 ! Computes the temperature and humidity evolution of indoor air,
19 ! building energy demand, HVAC energy consumption,
20 ! waste heat from HVAC systems, and heat fluxes from indoor to building surfaces.
21 !
22 !
23 !!** METHOD
24 ! ------
25 !
26 ! NOMENCLATURE: bld - refers to building plant area;
27 ! floor- refers to building plant area multiplied
28 ! by the number of floors;
29 ! wall - refers to wall area (excluding windows).
30 ! win - refers to window area.
31 ! mass - refers to internal mass area.
32 !
33 !
34 ! solar radiation transmitted through windows
35 ! *******************************************
36 !
37 ! Qsol_tr_win = Qsol_facade * tr_win * GR
38 !
39 !
40 ! indoor wall conv/rad heat transfer coefficients
41 ! ***********************************************
42 !
43 ! The calculation of CHTC accounts for favorable or unfavorable convection
44 ! depending on the relative position between the hot layer and cold layer
45 !
46 !
47 ! building energy demand
48 ! **********************
49 !
50 ! Calculation of the cooling and heating, sensible and latent building energy demand.
51 ! The sensible demand includes the convective heat transfer from indoor surfaces, the
52 ! convective fraction of internal heat gains, and sensible infiltration/ventilation heat
53 ! gains. The latent demand includes the latent fraction of internal heat gains and latent
54 ! infiltration/ventilation heat gains.
55 !
56 ! surface areas and volummes (referred to m2_bld)
57 ! ***********************************************
58 !
59 ! Awall = WALL_O_HOR * (1 - GR) / BLD [m2_wall/m2_bld]
60 ! Awin = WALL_O_HOR * GR / BLD [m2_win/m2_bld]
61 ! Amass = 2 * N_FLOOR [m2_mass/m2_bld]
62 ! N_FLOOR = BLD_HEIGHT / FLOOR_HEIGHT [#]
63 ! Aroof = 1 [m2_roof/m2_bld]
64 ! Afloor = 1 [m2_floor/m2_bld]
65 ! Vol_air = BLD_HEIGHT [m3_bld/m2_bld]
66 !
67 !
68 ! evolution of the internal temperature
69 ! *************************************
70 !
71 
72 ! dTin
73 ! Vol_air * ro_air * cp_air * ---- = h_wall * Awall * (Twall - Tin)
74 ! dt + h_roof * Aroof * (Troof -Tin)
75 ! + h_floor * Afloor *(Tfloor - Tin)
76 ! + h_mass * Amass * (Tmass - Tin)
77 ! + h_win * Awin * (Twin - Tin)
78 ! + Qig * (1 - fig_rad) * (1-fig_lat)
79 ! + Vinf * ro_air * cp_air * (Tout - Tin)
80 ! + Vsys * ro_air * cp_air * (Tsys - Tin)
81 !
82 !
83 ! evolution of the internal specific humidity
84 ! *******************************************
85 !
86 ! dQin
87 ! Vol_air * ro_air * lv_air * ---- = Qig * fig_lat
88 ! dt + Vinf * ro_air * lv_air * (Qout - Qin)
89 ! + Vsys * ro_air * lv_air * (Qsys - Qin)
90 !
91 !
92 ! heat fluxes from indoor to surfaces
93 ! ***********************************
94 !
95 ! Qin_wall = h_wall * (Tin - Twall) [W/m2_wall]
96 ! Qin_roof = h_roof * (Tin - Troof) [W/m2_roof]
97 ! Qin_floor = h_floor * (Tin - Tfloor) [W/m2_floor]
98 ! Qin_mass = h_wall * (Tin - Tmass)
99 ! + Qig * fig_rad * (1-fig_lat)/ 2
100 ! + Qsol_tr_win [W/m2_mass]
101 !
102 !
103 ! energy consumption and waste heat from cooling system
104 ! *****************************************************
105 !
106 ! Qhvac = Qbld / COP
107 ! Qwaste = Qbld + Qhvac
108 !
109 !
110 ! energy consumption and waste heat from heating system
111 ! *****************************************************
112 !
113 ! Qhvac = Qbld / Eff
114 ! Qwaste = Qhvac - Qbld
115 !
116 !
117 !! EXTERNAL
118 !! --------
119 !!
120 !!
121 !! IMPLICIT ARGUMENTS
122 !! ------------------
123 !!
124 !!
125 !! REFERENCE
126 !! ---------
127 !!
128 !!
129 !! AUTHOR
130 !! ------
131 !!
132 !! B. Bueno * Meteo-France *
133 !!
134 !!! MODIFICATIONS
135 !! -------------
136 !! Original 2010
137 !! G. Pigeon nov. 2011: inclusion floor/mass budget inside
138 !! add automatic/manual ventilation
139 !! conserve exchanges with the different surfaces inside 1 time step
140 !! G. Pigeon sept. 2012: use of TARP/DOE coef for indoor convection
141 !! use of both T%XT_WALL_A and T%XT_WALL_B for calculations
142 !! the internal mass depth is 1/2 of the floor depth
143 !! add the option of no atmospheric heat releases by HVAC system (B%XF_WATER_COND < 0)
144 !! G. Pigeon oct. 2012: use indoor air density + new solar heat gain distribution
145 !! V. Masson May 2013 implicitation of internal building temperature evolution
146 !-------------------------------------------------------------------------------
147 !
148 !* 0. DECLARATIONS
149 ! ------------
150 !
151 USE modd_bem_option_n, ONLY : bem_options_t
152 USE modd_teb_n, ONLY : teb_t
153 USE modd_bem_n, ONLY : bem_t
154 USE modd_diag_misc_teb_n, ONLY : diag_misc_teb_t
155 !
156 USE modd_csts,ONLY : xcpd,xstefan,xlvtt,xg, xrv, xrd
157 
158 USE mode_thermos
159 USE mode_psychro
160 USE modi_dx_air_cooling_coil_cv
161 USE modi_floor_layer_e_budget
162 USE modi_mass_layer_e_budget
163 USE mode_conv_doe
164 !
165 USE yomhook ,ONLY : lhook, dr_hook
166 USE parkind1 ,ONLY : jprb
167 !
168 IMPLICIT NONE
169 !
170 !* 0.1 Declarations of arguments
171 !
172 TYPE(bem_options_t), INTENT(INOUT) :: BOP
173 TYPE(teb_t), INTENT(INOUT) :: T
174 TYPE(bem_t), INTENT(INOUT) :: B
175 TYPE(diag_misc_teb_t), INTENT(INOUT) :: DMT
176 !
177 REAL, INTENT(IN) :: PTSTEP ! Time step
178 REAL, DIMENSION(:), INTENT(IN) :: PSUNTIME ! current solar time since midnight (solar time, s)
179 INTEGER, INTENT(IN) :: KDAY ! Simulation day
180 !
181 REAL, DIMENSION(:), INTENT(IN) :: PPS ! Canyon air pressure [Pa]
182 REAL, DIMENSION(:), INTENT(IN) :: PRHOA ! Air density at the lowest level [kg m-3]
183 REAL, DIMENSION(:), INTENT(IN) :: PT_CAN ! Canyon air temperature [K]
184 REAL, DIMENSION(:), INTENT(IN) :: PQ_CAN ! Canyon air specific humidity [kg kg-1]
185 REAL, DIMENSION(:), INTENT(IN) :: PU_CAN ! Canyon wind speed (m s-1)
186 !
187 REAL, DIMENSION(:), INTENT(OUT) :: PHU_BLD ! Indoor relative humidity 0 < (-) < 1
188 REAL, DIMENSION(:), INTENT(IN) :: PT_RAD_IND ! Indoor mean radiant temperature [K]
189 !
190 REAL, DIMENSION(:), INTENT(OUT) :: PFLX_BLD_FL! Heat flux from indoor air to floor
191  ! [W m-2(bld)]
192 REAL, DIMENSION(:), INTENT(OUT) :: PFLX_BLD_MA ! Heat flux from indoor air to mass
193  ! [W m-2(bld)]
194 REAL, DIMENSION(:), INTENT(IN) :: PRADHT_IN ! Indoor radiant heat transfer coefficient
195  ! [W K-1 m-2]
196 REAL, DIMENSION(:) , INTENT(IN) :: PRAD_RF_MA ! Rad. fluxes between roof and mass
197 REAL, DIMENSION(:) , INTENT(IN) :: PRAD_RF_FL ! Rad. fluxes between roof and floor
198 REAL, DIMENSION(:) , INTENT(IN) :: PRAD_WL_MA ! Rad. fluxes between wall and mass
199 REAL, DIMENSION(:) , INTENT(IN) :: PRAD_WL_FL ! Rad. fluxes between wall and floor
200 REAL, DIMENSION(:) , INTENT(IN) :: PRAD_WIN_MA ! Rad. fluxes between wind. and mass
201 REAL, DIMENSION(:) , INTENT(IN) :: PRAD_WIN_FL ! Rad. fluxes between wind. and floor
202 REAL, DIMENSION(:) , INTENT(IN) :: PCONV_RF_BLD ! Conv. fluxes between roof and indoor air
203 REAL, DIMENSION(:) , INTENT(IN) :: PCONV_WL_BLD ! Conv. fluxes between wall and indoor air
204 REAL, DIMENSION(:) , INTENT(IN) :: PCONV_WIN_BLD ! Conv. fluxes between wind. and indoor air
205 REAL, DIMENSION(:) , INTENT(IN) :: PLOAD_IN_FL ! solar + int heat gain on floor W/m2 [floor]
206 REAL, DIMENSION(:) , INTENT(IN) :: PLOAD_IN_MA ! solar + int heat gain on floor W/m2 [mass]
207 !
208 !* 0.2 Declarations of local variables
209 !
210 INTEGER :: IRF ! Number of roof layers
211 INTEGER :: IWL ! Number of wall layers
212 !REAL :: ZTCOMF_MAX ! Maximum comfort temperature for nat.vent [K]
213 !
214 REAL, DIMENSION(SIZE(B%XTI_BLD)) :: ZFAN_AP ! Fan design pressure increase [Pa]
215 REAL, DIMENSION(SIZE(B%XTI_BLD)) :: ZFAN_EFF ! Fan total efficiency
216 !
217 LOGICAL, DIMENSION(SIZE(B%XTI_BLD)):: GSCHED ! Day-night schedule flag
218  ! *to be transported to inputs*
219 !
220 REAL, DIMENSION(SIZE(B%XTI_BLD)) :: ZF_NIGHT ! Reduction factor of int.gains at night
221 REAL, DIMENSION(SIZE(B%XTI_BLD)) :: ZF_DAY ! Amplification factor of int.gains at daytime
222 !
223 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZAC_IN_MA_COOL, ZAC_IN_FL_COOL, &
224  ZAC_IN_RF_COOL, ZAC_IN_WL_A_COOL, &
225  ZAC_IN_WL_B_COOL, ZAC_IN_WIN_COOL
226 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZAC_IN_MA_HEAT, ZAC_IN_FL_HEAT, &
227  ZAC_IN_RF_HEAT, ZAC_IN_WL_A_HEAT, &
228  ZAC_IN_WL_B_HEAT, ZAC_IN_WIN_HEAT
229 !
230 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZQIN ! Internal heat gains [W m-2(bld)]
231 !
232 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZV_VENT ! Ventilation flow rate [m3 s-1 m-2(bld)]
233 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZINF ! Infiltration flow rate [m3 s-1 m-2(bld)]
234 !
235 LOGICAL, DIMENSION(SIZE(B%XTI_BLD)):: GNAT_VENT ! Is Natural ventilation active ?
236 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZNAT_VENT ! Nat.vent airflow rate [m3 s-1 m-2(bld)]
237 REAL,DIMENSION(SIZE(B%XTI_BLD)) :: ZTI_BLD ! Indoor air temperature at time step t + dt [K]
238 REAL,DIMENSION(SIZE(B%XTI_BLD)) :: ZTI_BLD_OPEN ! Indoor air temperature if windows opened
239 REAL,DIMENSION(SIZE(B%XTI_BLD)) :: ZTI_BLD_CLOSED! Indoor air temperature if windows closed
240 !
241 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZQCOOL_TRGT ! Specific humidity cooling setpoing [kg kg-1]
242 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZQHEAT_TRGT ! Specific humidity heating setpoing [kg kg-1]
243 !
244 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZSHR ! Rated sensible heat rate
245 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZM_SYS_RAT ! Auxiliar mass flow rate [kg s-1 m-2(bld)]
246 !
247 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZXMIX ! Outdoor mixing fraction
248 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZT_MIX ! Mixing air temperature [K]
249 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZQ_MIX ! Mixing air specific humidity [kg kg-1]
250 !
251 REAL,DIMENSION(SIZE(B%XTI_BLD)) :: ZQI_BLD ! Indoor air humidity at time step t + dt [K}
252 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZWASTE
253 !
254 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZDQS_FL
255 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZIMB_FL
256 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZDQS_MA
257 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZIMB_MA
258 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZLOAD_FL ! sum of solar and internal loads on floor
259 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZLOAD_MA ! sum of solar and internal loads on mass
260 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZRAD_FL_MA ! Rad. fluxes from floor to mass
261 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZCONV_FL_BLD ! Conv. fluxes from floor to indoor air
262 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZCONV_MA_BLD ! Conv. fluxes from mass to indoor air
263 REAL, DIMENSION(SIZE(B%XTI_BLD)):: ZRHOI ! indoor air density
264 !
265 INTEGER :: JJ ! Loop counter
266 REAL(KIND=JPRB) :: ZHOOK_HANDLE
267 !
268 !!REAL :: ZEXPL = 0.5 !explicit coefficient for internal temperature evol.
269 !!REAL :: ZIMPL = 0.5 !implicit coef..
270 !
271 !-------------------------------------------------------------------------------
272 IF (lhook) CALL dr_hook('BEM',0,zhook_handle)
273 !
274 !* 1. Initializations
275 ! ---------------
276 !
277 zrhoi(:) = pps(:) / (xrd * b%XTI_BLD(:) * ( 1.+((xrv/xrd)-1.)*b%XQI_BLD(:) ) )
278 ! *Temperal definitions for nat.vent*
279 !ZTCOMF_MAX = 26. + 273.16
280 !
281 ! *Definitions
282 zfan_ap(:) = 600.0
283 zfan_eff(:) = 0.7
284 !
285 ! *Other calcs
286 irf = SIZE(t%XT_ROOF,2)
287 iwl = SIZE(t%XT_WALL_A,2)
288 !
289 !
290 ! initial condition of QI_BLD equivalent to 50% RH
291 IF (any(b%XQI_BLD(:) <= 1e-6)) b%XQI_BLD = 0.5 * qsat(b%XTI_BLD, pps)
292 !
293 ! *Temperal definitions for shedule*
294 gsched(:) = .false.
295 WHERE (gsched(:))
296  zf_night(:) = 0.8
297  zf_day(:) = 1.2
298 ELSE WHERE
299  zf_night(:) = 1.
300  zf_day(:) = 1.
301 END WHERE
302 !
303 ! *Int.gains schedule
304 !
305 zqin = dmt%XQIN * b%XN_FLOOR
306 WHERE (psuntime(:) > 0. .AND. psuntime(:) < 25200.) ! night between 0000 and 0700
307  zqin(:) = zqin(:) * zf_night(:)
308 ELSEWHERE
309  zqin(:) = zqin(:) * zf_day(:)
310 END WHERE
311 
312 ! *Change of units AC/H -> [m3 s-1 m-2(bld)]
313 zv_vent(:) = b%XV_VENT(:) * t%XBLD_HEIGHT(:) / 3600.
314 zinf(:) = b%XINF (:) * t%XBLD_HEIGHT(:) / 3600.
315 !
316 !* 2. heat balance for building floor and mass
317 ! ----------------------------------------
318 !
319 !* 2.1 total load on the internal mass or floor
320 zload_fl(:) = (zqin(:) * b%XQIN_FRAD(:) * (1.-b%XQIN_FLAT(:)) + dmt%XTR_SW_WIN(:)) / (b%XMASS_O_BLD(:)+1.)
321 WHERE (b%XN_FLOOR(:) > 1.)
322  zload_ma(:) = zload_fl(:)
323 ELSEWHERE
324  zload_ma(:) = 0.
325 ENDWHERE
326 !
327 !* 2.2 FLOOR HEAT BALANCE
328 !
329  CALL floor_layer_e_budget(b, ptstep, pflx_bld_fl, zdqs_fl, zimb_fl, pradht_in, &
330  prad_wl_fl, prad_rf_fl, prad_win_fl, pload_in_fl, &
331  zrad_fl_ma, zconv_fl_bld)
332 !
333 !* 2.3 MASS HEAT BALANCE
334 !
335  CALL mass_layer_e_budget(b, ptstep, pflx_bld_ma, zdqs_ma, zimb_ma, pradht_in, &
336  prad_wl_ma, prad_rf_ma, prad_win_ma, pload_in_ma, &
337  zrad_fl_ma, zconv_ma_bld )
338 !
339 !
340 zac_in_wl_a_cool = chtc_vert_doe(t%XT_WALL_A(:,iwl), dmt%XTCOOL_TARGET)
341 zac_in_wl_b_cool = chtc_vert_doe(t%XT_WALL_B(:,iwl), dmt%XTCOOL_TARGET)
342 zac_in_win_cool = chtc_vert_doe(b%XT_WIN2 , dmt%XTCOOL_TARGET)
343 zac_in_ma_cool = chtc_vert_doe(b%XT_MASS(:,1) , dmt%XTCOOL_TARGET)
344 zac_in_rf_cool = chtc_down_doe(t%XT_ROOF(:,irf) , dmt%XTCOOL_TARGET)
345 zac_in_fl_cool = chtc_up_doe(b%XT_FLOOR(:,1) , dmt%XTCOOL_TARGET)
346 
347 zac_in_wl_a_heat = chtc_vert_doe(t%XT_WALL_A(:,iwl), dmt%XTHEAT_TARGET)
348 zac_in_wl_b_heat = chtc_vert_doe(t%XT_WALL_B(:,iwl), dmt%XTHEAT_TARGET)
349 zac_in_win_heat = chtc_vert_doe(b%XT_WIN2 , dmt%XTHEAT_TARGET)
350 zac_in_ma_heat = chtc_vert_doe(b%XT_MASS(:,1) , dmt%XTHEAT_TARGET)
351 zac_in_rf_heat = chtc_down_doe(t%XT_ROOF(:,irf) , dmt%XTHEAT_TARGET)
352 zac_in_fl_heat = chtc_up_doe(b%XT_FLOOR(:,1) , dmt%XTHEAT_TARGET)
353 
354 DO jj=1,SIZE(zac_in_win_cool)
355  zac_in_wl_a_cool(jj) = max(1., zac_in_wl_a_cool(jj))
356  zac_in_wl_b_cool(jj) = max(1., zac_in_wl_b_cool(jj))
357  zac_in_win_cool(jj) = max(1., zac_in_win_cool(jj))
358  zac_in_ma_cool(jj) = max(1., zac_in_ma_cool(jj))
359  zac_in_rf_cool(jj) = max(1., zac_in_rf_cool(jj))
360  zac_in_fl_cool(jj) = max(1., zac_in_fl_cool(jj))
361 
362  zac_in_wl_a_heat(jj) = max(1., zac_in_wl_a_heat(jj))
363  zac_in_wl_b_heat(jj) = max(1., zac_in_wl_b_heat(jj))
364  zac_in_win_heat(jj) = max(1., zac_in_win_heat(jj))
365  zac_in_ma_heat(jj) = max(1., zac_in_ma_heat(jj))
366  zac_in_rf_heat(jj) = max(1., zac_in_rf_heat(jj))
367  zac_in_fl_heat(jj) = max(1., zac_in_fl_heat(jj))
368 ENDDO
369 
370 !* 4. Indoor energy balance calculation
371 ! ---------------------------------
372 !
373 DO jj=1,SIZE(pt_can)
374  ! *first guess of indoor temperature
375 
376  zti_bld(jj) = b%XTI_BLD(jj) + ptstep/(zrhoi(jj) * xcpd * t%XBLD_HEIGHT(jj)) &
377  * ( t%XWALL_O_BLD(jj) * pconv_wl_bld(jj) + b%XGLAZ_O_BLD (jj) * pconv_win_bld(jj) &
378  + b%XMASS_O_BLD(jj) * zconv_ma_bld(jj) + pconv_rf_bld(jj) + zconv_fl_bld(jj) &
379  + zqin(jj) * (1 - b%XQIN_FRAD(jj)) * (1 - b%XQIN_FLAT(jj)) )
380  !
381  !################################################################################
382  ! *is natural surventilation active at the current time step ?
383  !---------------------------------------------------------------------------------
384  !
385  ! *no surventilation possible
386 
387  IF (b%CNATVENT(jj)=='NONE') THEN
388  !
389  gnat_vent(jj) = .false.
390  !
391  ! *automatic management of surventilation
392  ELSEIF (b%CNATVENT(jj)=='AUTO' .OR. b%CNATVENT(jj)=='MECH') THEN
393  !
394  IF (mod(psuntime(jj), 3600.) .LT. ptstep) THEN
395  !
396  IF ( b%XTI_BLD(jj).GT. pt_can(jj) + 1 ) THEN ! condition to enable the
397  IF (b%CNATVENT(jj)=='AUTO') THEN
398  ! natural surventilation rate calculation (window opening)
399  CALL get_nat_vent(b%XTI_BLD(jj), pt_can(jj), pu_can(jj), b%XGR(jj), &
400  b%XFLOOR_HW_RATIO(jj), t%XBLD_HEIGHT(jj), znat_vent(jj))
401  ELSE IF (b%CNATVENT(jj)=='MECH') THEN
402  ! mechanical surventilation rate calculation : 5 volumes/hour
403  znat_vent(jj) = 5.0*t%XBLD_HEIGHT(jj)/3600.
404  END IF
405  !
406  zti_bld_open(jj) = zti_bld(jj) &
407  + znat_vent(jj) * ptstep/t%XBLD_HEIGHT(jj) * (pt_can(jj) - b%XTI_BLD(jj))
408  zti_bld_closed(jj) = zti_bld(jj) &
409  + (zinf(jj) + zv_vent(jj)) * ptstep/t%XBLD_HEIGHT(jj) * (pt_can(jj) - b%XTI_BLD(jj))
410  !
411  gnat_vent(jj) = (zti_bld_open(jj) <= dmt%XTCOOL_TARGET (jj) .AND. &
412  zti_bld_open(jj) < zti_bld_closed(jj) .AND. &
413  zti_bld_open(jj) > dmt%XTHEAT_TARGET (jj) + 4.)
414  !
415  ELSE
416  gnat_vent(jj) = .false.
417  ENDIF
418  b%LNATVENT_NIGHT(jj) = gnat_vent(jj)
419  ELSE
420  gnat_vent(jj) = b%LNATVENT_NIGHT(jj)
421  ENDIF
422  !
423  ! *manual management of surventilation
424  ELSEIF (b%CNATVENT(jj)=='MANU') THEN
425  !
426  b%LNATVENT_NIGHT(jj) = b%LNATVENT_NIGHT(jj) .AND. &
427  .NOT. ( psuntime(jj) > 5.*3600 .AND. psuntime(jj) < 18.*3600 )
428  !
429  gnat_vent(jj) = ( psuntime(jj) > 18.*3600. .AND. psuntime(jj) < 21.*3600. &
430  .AND. pt_can(jj) < b%XTI_BLD(jj)+2. &
431  .AND. pt_can(jj) > dmt%XTHEAT_TARGET(jj) &
432  .AND. ( b%XTI_BLD(jj) > dmt%XTHEAT_TARGET(jj)+5. &
433  .OR. b%XTI_BLD(jj) == dmt%XTCOOL_TARGET(jj) ) )
434  gnat_vent(jj) = gnat_vent(jj) .OR. b%LNATVENT_NIGHT(jj)
435  !
436  ENDIF
437  !
438  ! Decicion about natural surventilation OK
439  !################################################################################
440  !
441  !
442  !################################################################################
443  ! COMPUTE ENERGY DEMAND
444  !---------------------------------------------------------------------------------
445 
446  ! *If natural surventilation ACTIVE
447  IF (gnat_vent(jj)) THEN
448  !
449  CALL get_nat_vent(b%XTI_BLD(jj), pt_can(jj), pu_can(jj), b%XGR(jj), &
450  b%XFLOOR_HW_RATIO(jj), t%XBLD_HEIGHT(jj), znat_vent(jj) )
451  !
452  zv_vent(jj) = 0.
453  zinf(jj) = 0.
454  !
455  dmt%XH_BLD_COOL (jj) = 0.0 ! No HVAC consumption
456  dmt%XH_BLD_HEAT (jj) = 0.0
457  dmt%XLE_BLD_COOL(jj) = 0.0 ! No HVAC consumption
458  dmt%XLE_BLD_HEAT(jj) = 0.0
459  !
460  dmt%XT_BLD_COOL (jj) = 0.0 ! No HVAC consumption
461  dmt%XHVAC_COOL (jj) = 0.0
462  dmt%XT_SYS (jj) = b%XTI_BLD(jj) ! No mechanical ventilation
463  dmt%XQ_SYS (jj) = b%XQI_BLD(jj) !
464  dmt%XH_WASTE (jj) = 0.0
465  dmt%XLE_WASTE (jj) = 0.0
466  dmt%XFAN_POWER (jj) = 0.0
467  dmt%XHVAC_HEAT (jj) = 0.0
468  !
469  dmt%XM_SYS (jj) = 0.0
470  dmt%XCOP (jj) = 0.0
471  dmt%XCAP_SYS(jj) = 0.0
472  !
473  ! *If natural surventilation INACTIVE
474  ELSE
475  !
476  znat_vent(jj) = 0.
477  !
478  ! ------------------------------------------------
479  ! * Building energy demand for heating and cooling
480  ! ------------------------------------------------
481  !
482  dmt%XH_BLD_COOL(jj) = t%XWALL_O_BLD(jj)/2. * (zac_in_wl_a_cool(jj) * (t%XT_WALL_A(jj,iwl) - dmt%XTCOOL_TARGET(jj)) &
483  + zac_in_wl_b_cool(jj) * (t%XT_WALL_B(jj,iwl) - dmt%XTCOOL_TARGET(jj))) &
484  + b%XGLAZ_O_BLD(jj) * zac_in_win_cool(jj) * (b%XT_WIN2(jj) - dmt%XTCOOL_TARGET(jj)) &
485  + zac_in_ma_cool(jj)* b%XMASS_O_BLD(jj) * (b%XT_MASS(jj,1) - dmt%XTCOOL_TARGET(jj)) &
486  + zac_in_rf_cool(jj) * (t%XT_ROOF(jj,irf) - dmt%XTCOOL_TARGET(jj)) &
487  + zac_in_fl_cool(jj) * (b%XT_FLOOR(jj,1) - dmt%XTCOOL_TARGET(jj)) &
488  + zqin(jj) * (1 - b%XQIN_FRAD(jj)) * (1 - b%XQIN_FLAT(jj)) &
489  + (zinf(jj) + zv_vent(jj)) * zrhoi(jj) * xcpd * (pt_can(jj) - dmt%XTCOOL_TARGET(jj))
490  !
491  dmt%XH_BLD_HEAT(jj) = - ( t%XWALL_O_BLD(jj)/2. * (zac_in_wl_a_heat(jj) * (t%XT_WALL_A(jj,iwl) - dmt%XTHEAT_TARGET(jj)) &
492  + zac_in_wl_b_heat(jj) * (t%XT_WALL_B(jj,iwl) - dmt%XTHEAT_TARGET(jj))) &
493  + b%XGLAZ_O_BLD(jj) * zac_in_win_heat(jj) * (b%XT_WIN2(jj) - dmt%XTHEAT_TARGET(jj)) &
494  + zac_in_ma_heat(jj)* b%XMASS_O_BLD(jj) * (b%XT_MASS(jj,1) - dmt%XTHEAT_TARGET(jj)) &
495  + zac_in_rf_heat(jj) * (t%XT_ROOF(jj,irf) - dmt%XTHEAT_TARGET(jj)) &
496  + zac_in_fl_heat(jj) * (b%XT_FLOOR(jj,1) - dmt%XTHEAT_TARGET(jj)) &
497  + zqin(jj) * (1 - b%XQIN_FRAD(jj))* (1 - b%XQIN_FLAT(jj)) &
498  + (zinf(jj) + zv_vent(jj)) * zrhoi(jj) * xcpd * (pt_can(jj) - dmt%XTHEAT_TARGET(jj)))
499  !
500  zqcool_trgt(jj) = 0.62198 * b%XHR_TARGET(jj) * psat(dmt%XTCOOL_TARGET(jj)) / &
501  (pps(jj)- b%XHR_TARGET(jj) * psat(dmt%XTCOOL_TARGET(jj)))
502  !
503  dmt%XLE_BLD_COOL(jj) = zqin(jj) * b%XQIN_FLAT(jj) &
504  + (zinf(jj) + zv_vent(jj)) * zrhoi(jj) * xlvtt * (pq_can(jj) - zqcool_trgt(jj))
505  !
506  zqheat_trgt(jj) = 0.62198 * b%XHR_TARGET(jj) * psat(dmt%XTHEAT_TARGET(jj)) / &
507  (pps(jj)- b%XHR_TARGET(jj) * psat(dmt%XTHEAT_TARGET(jj)))
508  !
509  dmt%XLE_BLD_HEAT(jj) = zqin(jj) * b%XQIN_FLAT(jj) &
510  + (zinf(jj) + zv_vent(jj)) * zrhoi(jj) * xlvtt * (pq_can(jj) - zqheat_trgt(jj))
511  !
512  ! * Autosize calculations
513  !
514  IF (bop%LAUTOSIZE .AND. kday==15) THEN
515  !
516  IF (dmt%XH_BLD_COOL(jj) > b%XAUX_MAX(jj)) THEN
517  !
518  b%XAUX_MAX (jj) = dmt%XH_BLD_COOL(jj)
519  !
520  ! Cooling coil sensible heat rate
521  zshr(jj) = min(xcpd * (dmt%XTCOOL_TARGET(jj) - b%XT_ADP(jj)) / &
522  (enth_fn_t_q(dmt%XTCOOL_TARGET(jj),zqcool_trgt(jj)) - &
523  enth_fn_t_q(b%XT_ADP(jj),qsat(b%XT_ADP(jj),pps(jj)))), 1.)
524  ! Cooling Coil Capacity [W m-2(bld)]
525  b%XCAP_SYS_RAT(jj) = dmt%XH_BLD_COOL(jj) / zshr(jj)
526  !
527  ! Cooling rated air flow rate [kg s-1 m-2(bld)]
528  zm_sys_rat(jj) = dmt%XH_BLD_COOL(jj) / xcpd / (dmt%XTCOOL_TARGET(jj)-(14.0+273.16))
529  IF (zm_sys_rat(jj) > b%XM_SYS_RAT(jj)) b%XM_SYS_RAT(jj) = zm_sys_rat(jj)
530  !
531  ! Impose condition
532  IF (b%XM_SYS_RAT(jj)/zrhoi(jj)/b%XCAP_SYS_RAT(jj) < 0.00004027) THEN
533  b%XCAP_SYS_RAT(jj) = b%XM_SYS_RAT(jj)/zrhoi(jj)/0.00004027
534  ELSE IF (b%XM_SYS_RAT(jj)/zrhoi(jj)/b%XCAP_SYS_RAT(jj) > 0.00006041) THEN
535  b%XCAP_SYS_RAT(jj) = b%XM_SYS_RAT(jj)/zrhoi(jj)/0.00006041
536  END IF
537  !
538  END IF
539  !
540  END IF
541  !
542  ! * END Autosize calculations
543  !
544  ! * system efficiency
545  ! ...................
546  !
547  dmt%XM_SYS (jj) = b%XM_SYS_RAT (jj)
548  dmt%XCOP (jj) = b%XCOP_RAT (jj)
549  dmt%XCAP_SYS(jj) = b%XCAP_SYS_RAT(jj)
550  !
551  ! * Mixing conditions
552  ! .................
553  !
554  zxmix(jj) = zv_vent(jj) * zrhoi(jj) / dmt%XM_SYS(jj)
555  zt_mix(jj) = zxmix(jj) * pt_can(jj) + (1.-zxmix(jj)) * b%XTI_BLD(jj)
556  zq_mix(jj) = zxmix(jj) * pq_can(jj) + (1.-zxmix(jj)) * b%XQI_BLD(jj)
557  !
558  ! ---------------------------------------------
559  ! * COOLING system : Performance and Waste heat
560  ! ---------------------------------------------
561  !
562  IF (dmt%XH_BLD_COOL(jj) >= 0.0) THEN
563  !
564  ! *ideal system
565  IF (bop%CCOOL_COIL=='IDEAL') THEN
566  !
567  dmt%XT_BLD_COOL(jj) = dmt%XH_BLD_COOL(jj) + dmt%XLE_BLD_COOL(jj)
568  !desactivation of LE_BLD_COOL impact on HVAC_COOL calculation
569  !following too much impact in VURCA simulation (23/01/2012)
570  !this would be the case for a vaporization system !
571  !DMT%XHVAC_COOL (JJ) = DMT%XT_BLD_COOL(JJ) / B%XCOP_RAT(JJ)
572  dmt%XHVAC_COOL (jj) = dmt%XH_BLD_COOL(jj) / b%XCOP_RAT(jj)
573  IF (dmt%XHVAC_COOL(jj) < 0.0) dmt%XHVAC_COOL(jj) = 0.0
574  !
575  dmt%XT_SYS(jj) = zt_mix(jj) - dmt%XH_BLD_COOL (jj) /dmt%XM_SYS(jj) / xcpd
576  !DMT%XQ_SYS(JJ) = ZQ_MIX(JJ) - DMT%XLE_BLD_COOL(JJ) / DMT%XM_SYS(JJ)/ XLVTT
577  !desactivation following too much impact in VURCA simulation
578  !(23/01/2012)
579  dmt%XQ_SYS(jj) = zq_mix(jj)
580  !
581  dmt%XH_WASTE(jj) = dmt%XHVAC_COOL(jj) * (1.+b%XCOP_RAT(jj)) * (1. - b%XF_WATER_COND(jj))
582  dmt%XLE_WASTE(jj) = dmt%XHVAC_COOL(jj) * (1.+b%XCOP_RAT(jj)) * b%XF_WATER_COND(jj)
583  !
584  ! *real system
585  ELSEIF (bop%CCOOL_COIL=='DXCOIL') THEN
586  !
587  CALL dx_air_cooling_coil_cv(pt_can(jj), pq_can(jj), pps(jj), zrhoi(jj), zt_mix(jj), &
588  zq_mix(jj), b%XCOP_RAT(jj), b%XCAP_SYS_RAT(jj), &
589  b%XT_ADP(jj), b%XF_WATER_COND(jj), dmt%XM_SYS(jj), &
590  dmt%XH_BLD_COOL(jj), dmt%XH_WASTE(jj), dmt%XLE_WASTE(jj), &
591  dmt%XCOP(jj), dmt%XCAP_SYS(jj), dmt%XT_SYS(jj), &
592  dmt%XQ_SYS(jj), dmt%XHVAC_COOL(jj), dmt%XT_BLD_COOL(jj) )
593  !
594  ENDIF !end type of cooling system
595 
596  !!! case of system without atmospheric releases. I-e releases in soil/water F_WATER_COND < 0
597  IF (b%XF_WATER_COND(jj) < 0) THEN
598  dmt%XH_WASTE(jj) = 0.
599  dmt%XLE_WASTE(jj) = 0.
600  ENDIF
601  !!!!
602  !
603  ! From EP Engineering Reference (p. 647)
604  dmt%XFAN_POWER(jj) = dmt%XM_SYS(jj) * zfan_ap(jj) * zfan_eff(jj) * zrhoi(jj)
605  !
606  dmt%XH_BLD_HEAT (jj) = 0.0
607  dmt%XLE_BLD_HEAT(jj) = 0.0
608  dmt%XHVAC_HEAT (jj) = 0.0
609  !
610  ! ---------------------------------------------
611  ! * HEATING system : Performance and Waste heat
612  ! ---------------------------------------------
613  !
614  ELSE IF (dmt%XH_BLD_HEAT(jj) > 0.0) THEN
615  !
616  ! *specific computation for real heating system
617  IF (bop%CHEAT_COIL .EQ. 'FINCAP') THEN
618  IF (dmt%XH_BLD_HEAT(jj) > b%XCAP_SYS_HEAT(jj)) dmt%XH_BLD_HEAT(jj) = b%XCAP_SYS_HEAT(jj)
619  END IF
620  !
621  dmt%XT_SYS(jj) = zt_mix(jj) + dmt%XH_BLD_HEAT(jj) / dmt%XM_SYS(jj) / xcpd
622  dmt%XQ_SYS(jj) = zq_mix(jj)
623  !
624  dmt%XHVAC_HEAT (jj) = dmt%XH_BLD_HEAT(jj) / b%XEFF_HEAT(jj)
625  dmt%XH_WASTE (jj) = dmt%XHVAC_HEAT(jj) - dmt%XH_BLD_HEAT(jj)
626  dmt%XLE_WASTE (jj) = 0.0
627  dmt%XH_BLD_COOL (jj) = 0.0
628  dmt%XLE_BLD_COOL(jj) = 0.0
629  dmt%XT_BLD_COOL (jj) = 0.0
630  dmt%XHVAC_COOL (jj) = 0.0
631 ! From EP Engineering Reference (p. 647)
632  dmt%XFAN_POWER(jj) = dmt%XM_SYS(jj)*zfan_ap(jj)*(zfan_eff(jj)*zrhoi(jj))
633  !
634  ! ------------------------------
635  ! * NEITHEIR COOLING NOR HEATING
636  ! ------------------------------
637  !
638  ELSE
639  !
640  dmt%XH_BLD_COOL (jj) = 0.0
641  dmt%XH_BLD_HEAT (jj) = 0.0
642  dmt%XLE_BLD_COOL(jj) = 0.0
643  dmt%XLE_BLD_HEAT(jj) = 0.0
644  !
645  dmt%XT_BLD_COOL (jj) = 0.0
646  dmt%XHVAC_COOL (jj) = 0.0
647  dmt%XT_SYS (jj) = zt_mix(jj)
648  dmt%XQ_SYS (jj) = zq_mix(jj)
649  dmt%XH_WASTE (jj) = 0.0
650  dmt%XLE_WASTE (jj) = 0.0
651  dmt%XFAN_POWER (jj) = 0.0
652  dmt%XHVAC_HEAT (jj) = 0.0
653  !
654  END IF !end for heating/cooling sytem
655  !
656  END IF
657  !
658  !---------------------------------------------------------------------------------
659  ! ENERGY DEMAND COMPUTED
660  !################################################################################
661 ENDDO
662 !
663 !---------------------------------------------------
664 ! EVOLUTION OF THE INTERNAL TEMPERATURE AND HUMIDITY
665 !###################################################
666 !
667 zti_bld(:) = ( zti_bld(:) + ptstep/t%XBLD_HEIGHT(:) * &
668  ((zinf(:) + znat_vent(:)) * pt_can(:) + dmt%XM_SYS(:) / zrhoi(:) * (dmt%XT_SYS(:) ) )) &
669  / (1. + ptstep/t%XBLD_HEIGHT(:)*(zinf(:) + znat_vent(:) + dmt%XM_SYS(:) / zrhoi(:)) )
670 zqi_bld(:) = ( b%XQI_BLD(:) + ptstep/t%XBLD_HEIGHT(:) * &
671  (zqin(:) * b%XQIN_FLAT(:) / (zrhoi(:) * xlvtt) + (zinf(:) + znat_vent(:)) * (pq_can(:)) &
672  + dmt%XM_SYS(:) / zrhoi(:) * (dmt%XQ_SYS(:) ) ))&
673  / (1. + ptstep/t%XBLD_HEIGHT(:)* (zinf(:) + znat_vent(:) + dmt%XM_SYS(:) / zrhoi(:)) )
674 !
675 ! Update variables
676 b%XTI_BLD(:) = zti_bld(:)
677 b%XQI_BLD(:) = zqi_bld(:)
678 !
679 ! Waste heat due to infiltration/ventilation
680 zwaste(:) = (zinf(:)+zv_vent(:)+znat_vent(:)) * zrhoi(:)
681 dmt%XH_WASTE (:) = dmt%XH_WASTE (:) + zwaste(:) * xcpd * (b%XTI_BLD(:) - pt_can(:))
682 dmt%XLE_WASTE(:) = dmt%XLE_WASTE(:) + zwaste(:) * xlvtt * (b%XQI_BLD(:) - pq_can(:))
683 !
684 !
685 IF (lhook) CALL dr_hook('BEM',1,zhook_handle)
686 !
687 CONTAINS
688 !
689 SUBROUTINE get_nat_vent(PTI_BLD, PPT_CAN, PPU_CAN, PGR, PF_AUX, PPBLD_HEIGHT, PNAT_VENT)
690 !
691 IMPLICIT NONE
692 !
693 REAL, INTENT(IN) :: PTI_BLD
694 REAL, INTENT(IN) :: PPT_CAN
695 REAL, INTENT(IN) :: PPU_CAN
696 REAL, INTENT(IN) :: PGR
697 REAL, INTENT(IN) :: PF_AUX
698 REAL, INTENT(IN) :: PPBLD_HEIGHT
699 REAL, INTENT(OUT) :: PNAT_VENT
700 REAL(KIND=JPRB) :: ZHOOK_HANDLE
701 !
702 IF (lhook) CALL dr_hook('BEM:GET_NAT_VENT',0,zhook_handle)
703 !
704 pnat_vent = xg * (pti_bld - ppt_can)
705 IF (pnat_vent .LT. 0.) THEN ! exceptional case with MANU ventilation system
706  pnat_vent= ppbld_height/3600. !minimum value
707 ELSE
708  pnat_vent = 1./3. * (pnat_vent/ppt_can)**(1./2.) &
709  * (1.5 + pti_bld/pnat_vent * 1./2. * ppu_can**2*0.1)**(3./2.) &
710  * pgr * pf_aux / 1.5 / 2.
711  pnat_vent = min(pnat_vent, 5.0*ppbld_height/3600.)
712 ENDIF
713 !
714 IF (lhook) CALL dr_hook('BEM:GET_NAT_VENT',1,zhook_handle)
715 !
716 END SUBROUTINE get_nat_vent
717 !
718 END SUBROUTINE bem
modd_bem_option_n
Definition: modd_bem_optionn.F90:6
mass_layer_e_budget
subroutine mass_layer_e_budget(B, PTSTEP, PFLX_BLD_MA, PDQS_MA, PIMB_MA, PRADHT_IN, PRAD_WL_MA, PRAD_RF_MA, PRAD_WIN_MA, PLOAD_MA, PRAD_FL_MA, PCONV_MA_BLD)
Definition: mass_layer_e_budget.F90:9
mode_conv_doe::chtc_down_doe
real function, dimension(size(pts)) chtc_down_doe(PTS, PTA)
Definition: mode_conv_DOE.F90:319
modd_diag_misc_teb_n
Definition: modd_diag_misc_tebn.F90:6
modd_csts::xcpd
real, save xcpd
Definition: modd_csts.F90:63
yomhook::dr_hook
Definition: yomhook.F90:20
mode_conv_doe::chtc_vert_doe
Definition: mode_conv_DOE.F90:40
modd_teb_n::teb_t
Definition: modd_tebn.F90:46
bem
subroutine bem(BOP, T, B, DMT, PTSTEP, PSUNTIME, KDAY, PPS, PRHOA, PT_CAN, PQ_CAN, PU_CAN, PHU_BLD, PT_RAD_IND, PFLX_BLD_FL, PFLX_BLD_MA, PRADHT_IN, PRAD_RF_MA, PRAD_RF_FL, PRAD_WL_MA, PRAD_WL_FL, PRAD_WIN_MA, PRAD_WIN_FL, PCONV_RF_BLD, PCONV_WL_BLD, PCONV_WIN_BLD, PLOAD_IN_FL, PLOAD_IN_MA)
Definition: bem.F90:11
modd_csts::xstefan
real, save xstefan
Definition: modd_csts.F90:59
modd_bem_n::bem_t
Definition: modd_bemn.F90:48
modd_csts::xlvtt
real, save xlvtt
Definition: modd_csts.F90:70
mode_conv_doe
Definition: mode_conv_DOE.F90:6
modd_csts::xrd
real, save xrd
Definition: modd_csts.F90:62
modd_bem_option_n::bem_options_t
Definition: modd_bem_optionn.F90:46
modd_csts::xg
real, save xg
Definition: modd_csts.F90:55
parkind1::jprb
integer, parameter jprb
Definition: parkind1.F90:32
modd_diag_misc_teb_n::diag_misc_teb_t
Definition: modd_diag_misc_tebn.F90:43
get_nat_vent
subroutine get_nat_vent(PTI_BLD, PPT_CAN, PPU_CAN, PGR, PF_AUX, PPBLD_HEIGHT, PNAT_VENT)
Definition: bem.F90:690
mode_conv_doe::chtc_up_doe
Definition: mode_conv_DOE.F90:44
mode_thermos
Definition: mode_thermos.F90:6
modd_csts::xrv
real, save xrv
Definition: modd_csts.F90:62
modd_teb_n
Definition: modd_tebn.F90:6
mode_thermos::qsat
Definition: mode_thermos.F90:57
mode_psychro
Definition: mode_psychro.F90:6
modd_bem_n
Definition: modd_bemn.F90:6
yomhook::lhook
logical lhook
Definition: yomhook.F90:15
dx_air_cooling_coil_cv
subroutine dx_air_cooling_coil_cv(PT_CANYON, PQ_CANYON, PPS, PRHOA, PT_IN, PQ_IN, PCOP_RAT, PCAP_SYS_RAT, PT_ADP, PF_WATER_COND, PM_SYS, PH_BLD_COOL, PH_WASTE, PLE_WASTE, PCOP, PCAP_SYS, PT_OUT, PQ_OUT, PDX_POWER, PT_BLD_COOL)
Definition: dx_air_cooling_coil_cv.F90:11
modd_csts
Definition: modd_csts.F90:6
floor_layer_e_budget
subroutine floor_layer_e_budget(B, PTSTEP, PFLX_BLD_FL, PDQS_FL, PIMB_FL, PRADHT_IN, PRAD_WL_FL, PRAD_RF_FL, PRAD_WIN_FL, PLOAD_FL, PRAD_FL_MA, PCONV_FL_BLD)
Definition: floor_layer_e_budget.F90:9
mode_psychro::enth_fn_t_q
Definition: mode_psychro.F90:52
mode_thermos::psat
Definition: mode_thermos.F90:48
Generated on Tue Jan 16 2018 16:23:22 for SURFEX v8.1 by   doxygen 1.8.13