!> @file biometeorology_pt_mod.f90 !--------------------------------------------------------------------------------! ! This file is part of PALM-4U. ! ! PALM-4U is free software: you can redistribute it and/or modify it under the ! terms of the GNU General Public License as published by the Free Software ! Foundation, either version 3 of the License, or (at your option) any later ! version. ! ! PALM-4U is distributed in the hope that it will be useful, but WITHOUT ANY ! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR ! A PARTICULAR PURPOSE. See the GNU General Public License for more details. ! ! You should have received a copy of the GNU General Public License along with ! PALM. If not, see . ! ! Copyright 1997-2019, Deutscher Wetterdienst (DWD) / ! German Meteorological Service (DWD) and ! Leibniz Universitaet Hannover !--------------------------------------------------------------------------------! ! ! Current revisions: ! ------------------ ! ! ! Former revisions: ! ----------------- ! $Id$ ! From branch resler@3462, pavelkrc: ! replace specific precision DLOG, DMOD intrinsics with generic precision for ! compatibility with 'wp' ! ! ! Initial revision ! ! ! ! Authors: ! -------- ! @author Dominik Froehlich ! ! ! Description: ! ------------ !> Module for the calculation of the thermal index perceived temperature (PT). !> The Perceived Temperature is the air temperature in a reference environ- !> ment that results in the same Predicted Mean Vote (Fanger 1972), PMV, as !> the real environment. In the reference environment the mean radiant !> temperature is equal to the air temperature and the wind speed is reduced !> to a value set at 0.1 m/s. The Perceived Temperature is linked to a !> reference individual: male aged 35 years, body height 1.75 m, weight !> 75 kg, walking with 4 km per hour on a horizontal plain (related to an !> internal heat production of 172.5 W). The person varies the thermal !> resistance of clothing (clo) in the range Icl=1.75 clo (winter) and !> Icl=0.50 clo (summer) to achieve thermal comfort (PMV = 0) if possible. !> If not possible the individual perceives cold stress (winter, 1.75 clo) !> or heat load (summer, 0.5 clo). !> !> @todo !> !> @note nothing now !> !> @bug no known bugs by now !------------------------------------------------------------------------------! MODULE biometeorology_pt_mod !-- Load required variables from existing modules USE kinds !< to set precision of INTEGER and REAL arrays according to PALM IMPLICIT NONE PRIVATE !-- no private interfaces --! PUBLIC calculate_pt_static, saturation_vapor_pressure, deltapmv, dpmv_adj, & dpmv_cold, fanger, calc_sultr, pt_regression, ireq_neutral, iso_ridder !-- PALM interfaces: INTERFACE calculate_pt_static MODULE PROCEDURE calculate_pt_static END INTERFACE calculate_pt_static INTERFACE saturation_vapor_pressure MODULE PROCEDURE saturation_vapor_pressure END INTERFACE saturation_vapor_pressure INTERFACE deltapmv MODULE PROCEDURE deltapmv END INTERFACE deltapmv INTERFACE dpmv_adj MODULE PROCEDURE dpmv_adj END INTERFACE dpmv_adj INTERFACE dpmv_cold MODULE PROCEDURE dpmv_cold END INTERFACE dpmv_cold INTERFACE fanger MODULE PROCEDURE fanger END INTERFACE fanger INTERFACE calc_sultr MODULE PROCEDURE calc_sultr END INTERFACE calc_sultr INTERFACE pt_regression MODULE PROCEDURE pt_regression END INTERFACE pt_regression INTERFACE ireq_neutral MODULE PROCEDURE ireq_neutral END INTERFACE ireq_neutral INTERFACE iso_ridder MODULE PROCEDURE iso_ridder END INTERFACE iso_ridder CONTAINS !------------------------------------------------------------------------------! ! Description: ! ------------ !> PT_BASIC.F90 Version of perceived temperature (PT, �C) for !> - standard measured/predicted meteorological values and TMRT !> as input; !> - regressions for determination of PT; !> - adjustment to Gagge's PMV* (2-node-model, 1986) as base of PT !> under warm/humid conditions (Icl= 0.50 clo) and under cold !> conditions (Icl= 1.75 clo) !> !> Function value is the Perceived Temperature, degree centigrade !------------------------------------------------------------------------------! SUBROUTINE calculate_pt_static( tt2m, el2m, vau1m, tmrt, pb, clo, pt_basic ) IMPLICIT NONE !-- Type of input of the argument list REAL(wp), INTENT ( IN ) :: tt2m !< Local air temperature (�C) REAL(wp), INTENT ( IN ) :: el2m !< Local vapour pressure (hPa) REAL(wp), INTENT ( IN ) :: tmrt !< Local mean radiant temperature (�C) REAL(wp), INTENT ( IN ) :: vau1m !< Local wind velocitry (m/s) REAL(wp), INTENT ( IN ) :: pb !< Local barometric air pressure (hPa) !-- Type of output of the argument list REAL(wp), INTENT ( OUT ) :: pt_basic !< Perceived temperature (�C) REAL(wp), INTENT ( OUT ) :: clo !< Clothing index (dimensionless) !-- Parameters for standard "Klima-Michel" REAL(wp), PARAMETER :: eta = 0._wp !< Mechanical work efficiency for walking on flat ground (compare to Fanger (1972) pp 24f) REAL(wp), PARAMETER :: actlev = 134.6862_wp !< Workload by activity per standardized surface (A_Du) !-- Type of program variables REAL(wp), PARAMETER :: eps = 0.0005 !< Accuracy in clothing insulation (clo) for evaluation the root of Fanger's PMV (pmva=0) REAL(wp) :: uclo REAL(wp) :: oclo REAL(wp) :: d_pmv REAL(wp) :: svp_tt REAL(wp) :: sult_lim REAL(wp) :: dgtcm REAL(wp) :: dgtcstd REAL(wp) :: clon REAL(wp) :: ireq_minimal REAL(wp) :: clo_fanger REAL(wp) :: pmv_o REAL(wp) :: pmv_u REAL(wp) :: pmva REAL(wp) :: ptc REAL(wp) :: d_std REAL(wp) :: pmvs REAL(wp) :: top INTEGER(iwp) :: nzaehl, nerr_kalt, nerr LOGICAL :: sultrieness !-- Initialise pt_basic = 9999.0_wp nerr = 0_iwp nzaehl = 0_iwp sultrieness = .FALSE. !-- Tresholds: clothing insulation (account for model inaccuracies) ! summer clothing uclo = 0.44453_wp ! winter clothing oclo = 1.76267_wp !-- decision: firstly calculate for winter or summer clothing IF ( tt2m <= 10._wp ) THEN !-- First guess: winter clothing insulation: cold stress clo = oclo CALL fanger ( tt2m, tmrt, el2m, vau1m, pb, clo, actlev, eta, pmva, top ) pmv_o = pmva IF ( pmva > 0._wp ) THEN !-- Case summer clothing insulation: heat load ? clo = uclo CALL fanger ( tt2m, tmrt, el2m, vau1m, pb, clo, actlev, eta, pmva, & top ) pmv_u = pmva IF ( pmva <= 0._wp ) THEN !-- Case: comfort achievable by varying clothing insulation !-- Between winter and summer set values CALL iso_ridder ( tt2m, tmrt, el2m, vau1m, pb, actlev, eta, uclo, & pmv_u, oclo, pmv_o, eps, pmva, top, nzaehl, clo ) IF ( nzaehl < 0_iwp ) THEN nerr = -1_iwp RETURN ENDIF ELSE IF ( pmva > 0.06_wp ) THEN clo = 0.5_wp CALL fanger ( tt2m, tmrt, el2m, vau1m, pb, clo, actlev, eta, & pmva, top ) ENDIF ELSE IF ( pmva < -0.11_wp ) THEN clo = 1.75_wp CALL fanger ( tt2m, tmrt, el2m, vau1m, pb, clo, actlev, eta, pmva, & top ) ENDIF ELSE !-- First guess: summer clothing insulation: heat load clo = uclo CALL fanger ( tt2m, tmrt, el2m, vau1m, pb, clo, actlev, eta, pmva, top ) pmv_u = pmva IF ( pmva < 0._wp ) THEN !-- Case winter clothing insulation: cold stress ? clo = oclo CALL fanger ( tt2m, tmrt, el2m, vau1m, pb, clo, actlev, eta, pmva, & top ) pmv_o = pmva IF ( pmva >= 0._wp ) THEN !-- Case: comfort achievable by varying clothing insulation ! between winter and summer set values CALL iso_ridder ( tt2m, tmrt, el2m, vau1m, pb, actlev, eta, uclo, & pmv_u, oclo, pmv_o, eps, pmva, top, nzaehl, clo ) IF ( nzaehl < 0_iwp ) THEN nerr = -1_iwp RETURN ENDIF ELSE IF ( pmva < -0.11_wp ) THEN clo = 1.75_wp CALL fanger ( tt2m, tmrt, el2m, vau1m, pb, clo, actlev, eta, & pmva, top ) ENDIF ELSE IF ( pmva > 0.06_wp ) THEN clo = 0.5_wp CALL fanger ( tt2m, tmrt, el2m, vau1m, pb, clo, actlev, eta, pmva, & top ) ENDIF ENDIF !-- Determine perceived temperature by regression equation + adjustments pmvs = pmva CALL pt_regression ( pmva, clo, pt_basic ) ptc = pt_basic IF ( clo >= 1.75_wp .AND. pmva <= -0.11_wp ) THEN !-- Adjust for cold conditions according to Gagge 1986 CALL dpmv_cold ( pmva, tt2m, vau1m, tmrt, nerr_kalt, d_pmv ) IF ( nerr_kalt > 0_iwp ) nerr = -5_iwp pmvs = pmva - d_pmv IF ( pmvs > -0.11_wp ) THEN d_pmv = 0._wp pmvs = -0.11_wp ENDIF CALL pt_regression ( pmvs, clo, pt_basic ) ENDIF clo_fanger = clo clon = clo IF ( clo > 0.5_wp .AND. pt_basic <= 8.73_wp ) THEN !-- Required clothing insulation (ireq) is exclusively defined for ! operative temperatures (top) less 10 (C) for a ! reference wind of 0.2 m/s according to 8.73 (C) for 0.1 m/s clon = ireq_neutral ( pt_basic, ireq_minimal, nerr ) clo = clon ENDIF CALL calc_sultr ( ptc, dgtcm, dgtcstd, sult_lim ) sultrieness = .FALSE. d_std = -99._wp IF ( pmva > 0.06_wp .AND. clo <= 0.5_wp ) THEN !-- Adjust for warm/humid conditions according to Gagge 1986 CALL saturation_vapor_pressure ( tt2m, svp_tt ) d_pmv = deltapmv ( pmva, tt2m, el2m, svp_tt, tmrt, vau1m, nerr ) pmvs = pmva + d_pmv CALL pt_regression ( pmvs, clo, pt_basic ) IF ( sult_lim < 99._wp ) THEN IF ( (pt_basic - ptc) > sult_lim ) sultrieness = .TRUE. !-- Set factor to threshold for sultriness IF ( dgtcstd /= 0_iwp ) THEN d_std = ( ( pt_basic - ptc ) - dgtcm ) / dgtcstd ENDIF ENDIF ENDIF END SUBROUTINE calculate_pt_static !------------------------------------------------------------------------------! ! Description: ! ------------ !> The SUBROUTINE calculates the saturation water vapour pressure !> (hPa = hecto Pascal) for a given temperature tt (�C). !> For example, tt can be the air temperature or the dew point temperature. !------------------------------------------------------------------------------! SUBROUTINE saturation_vapor_pressure( tt, p_st ) IMPLICIT NONE REAL(wp), INTENT ( IN ) :: tt !< ambient air temperature (�C) REAL(wp), INTENT ( OUT ) :: p_st !< saturation water vapour pressure (hPa) REAL(wp) :: b REAL(wp) :: c IF ( tt < 0._wp ) THEN !-- tt < 0 (�C): saturation water vapour pressure over ice b = 17.84362_wp c = 245.425_wp ELSE !-- tt >= 0 (�C): saturation water vapour pressure over water b = 17.08085_wp c = 234.175_wp ENDIF !-- Saturation water vapour pressure p_st = 6.1078_wp * EXP ( b * tt / ( c + tt ) ) END SUBROUTINE saturation_vapor_pressure !------------------------------------------------------------------------------! ! Description: ! ------------ ! Find the clothing insulation value clo_res (clo) to make Fanger's Predicted ! Mean Vote (PMV) equal comfort (pmva=0) for actual meteorological conditions ! (tt2m,tmrt, el2m, vau1m, pb) and values of individual's activity level ! (Klima-Michel: 134.682 W/m2 = 2.3 met, mechanincal work load eta= 0 W/m2), ! i.e. find the root of Fanger's comfort equation !------------------------------------------------------------------------------! SUBROUTINE iso_ridder( tt2m, tmrt, el2m, vau1m, pb, actlev, eta, uclo, & pmv_u, oclo, pmv_o, eps, pmva, top, nerr, & clo_res ) IMPLICIT NONE !-- Input variables of argument list: REAL(wp), INTENT ( IN ) :: tt2m !< Ambient temperature (�C) REAL(wp), INTENT ( IN ) :: tmrt !< Mean radiant temperature (�C) REAL(wp), INTENT ( IN ) :: el2m !< Water vapour pressure (hPa) REAL(wp), INTENT ( IN ) :: vau1m !< Wind speed (m/s) 1 m above ground REAL(wp), INTENT ( IN ) :: pb !< Barometric pressure (hPa) REAL(wp), INTENT ( IN ) :: actlev !< Individuals activity level per unit surface area (W/m2) REAL(wp), INTENT ( IN ) :: eta !< Individuals work efficiency (dimensionless) REAL(wp), INTENT ( IN ) :: uclo !< Lower threshold of bracketing clothing insulation (clo) REAL(wp), INTENT ( IN ) :: oclo !< Upper threshold of bracketing clothing insulation (clo) REAL(wp), INTENT ( IN ) :: eps !< (0.05) accuracy in clothing insulation (clo) for ! evaluation the root of Fanger's PMV (pmva=0) REAL(wp), INTENT ( IN ) :: pmv_o !< Fanger's PMV corresponding to oclo REAL(wp), INTENT ( IN ) :: pmv_u !< Fanger's PMV corresponding to uclo ! Output variables of argument list: REAL(wp), INTENT ( OUT ) :: pmva !< 0 (set to zero, because clo is evaluated for comfort) REAL(wp), INTENT ( OUT ) :: top !< Operative temperature (�C) at found root of Fanger's PMV REAL(wp), INTENT ( OUT ) :: clo_res !< Resulting clothing insulation value (clo) INTEGER(iwp), INTENT ( OUT ) :: nerr !< Error status / quality flag ! nerr >= 0, o.k., and nerr is the number of iterations for ! convergence ! nerr = -1: error = malfunction of Ridder's convergence method ! nerr = -2: error = maximum iterations (max_iteration) exceeded ! nerr = -3: error = root not bracketed between uclo and oclo !-- Type of program variables INTEGER(iwp), PARAMETER :: max_iteration = 15_iwp REAL(wp), PARAMETER :: guess_0 = -1.11e30_wp REAL(wp) :: x_ridder REAL(wp) :: clo_u REAL(wp) :: clo_o REAL(wp) :: x_unten REAL(wp) :: x_oben REAL(wp) :: x_mittel REAL(wp) :: x_neu REAL(wp) :: y_unten REAL(wp) :: y_oben REAL(wp) :: y_mittel REAL(wp) :: y_neu REAL(wp) :: s INTEGER(iwp) :: j !-- Initialise nerr = 0_iwp !-- Set pmva = 0 (comfort): Root of PMV depending on clothing insulation pmva = 0._wp clo_u = uclo y_unten = pmv_u clo_o = oclo y_oben = pmv_o IF ( ( y_unten > 0._wp .AND. y_oben < 0._wp ) .OR. & ( y_unten < 0._wp .AND. y_oben > 0._wp ) ) THEN x_unten = clo_u x_oben = clo_o x_ridder = guess_0 DO j = 1_iwp, max_iteration x_mittel = 0.5_wp * ( x_unten + x_oben ) CALL fanger ( tt2m, tmrt, el2m, vau1m, pb, x_mittel, actlev, eta, & y_mittel, top ) s = SQRT ( y_mittel**2 - y_unten * y_oben ) IF ( s == 0._wp ) THEN clo_res = x_mittel nerr = j RETURN ENDIF x_neu = x_mittel + ( x_mittel - x_unten ) * & ( SIGN ( 1._wp, y_unten - y_oben ) * y_mittel / s ) IF ( ABS ( x_neu - x_ridder ) <= eps ) THEN clo_res = x_ridder nerr = j RETURN ENDIF x_ridder = x_neu CALL fanger ( tt2m, tmrt, el2m, vau1m, pb, x_ridder, actlev, eta, & y_neu, top ) IF ( y_neu == 0._wp ) THEN clo_res = x_ridder nerr = j RETURN ENDIF IF ( SIGN ( y_mittel, y_neu ) /= y_mittel ) THEN x_unten = x_mittel y_unten = y_mittel x_oben = x_ridder y_oben = y_neu ELSE IF ( SIGN ( y_unten, y_neu ) /= y_unten ) THEN x_oben = x_ridder y_oben = y_neu ELSE IF ( SIGN ( y_oben, y_neu ) /= y_oben ) THEN x_unten = x_ridder y_unten = y_neu ELSE !-- Never get here in x_ridder: singularity in y nerr = -1_iwp clo_res = x_ridder RETURN ENDIF IF ( ABS ( x_oben - x_unten ) <= eps ) THEN clo_res = x_ridder nerr = j RETURN ENDIF ENDDO !-- x_ridder exceed maximum iterations nerr = -2_iwp clo_res = y_neu RETURN ELSE IF ( y_unten == 0. ) THEN x_ridder = clo_u ELSE IF ( y_oben == 0. ) THEN x_ridder = clo_o ELSE !-- x_ridder not bracketed by u_clo and o_clo nerr = -3_iwp clo_res = x_ridder RETURN ENDIF END SUBROUTINE iso_ridder !------------------------------------------------------------------------------! ! Description: ! ------------ !> Regression relations between perceived temperature (gt) and (adjusted) !> PMV (21.11.1996, files F3PMVPGT.TXT und F3CLOGT.TXT). Three cases: !> !> - PMV < 0: gt = 5.805 + 12.6784*pmv !> - PMV = 0 (gt depends on clothing insulation, FUNCTION iso_ridder): !> gt = 21.258 - 9.558*clo !> - PMV > 0: 16.826 + 6.163*pmv !> !> The regression presumes the Klima-Michel settings for reference !> individual and reference environment. !------------------------------------------------------------------------------! SUBROUTINE pt_regression( pmv, clo, pt ) IMPLICIT NONE REAL(wp), INTENT ( IN ) :: pmv !< Fangers predicted mean vote (dimensionless) REAL(wp), INTENT ( IN ) :: clo !< clothing insulation index (clo) REAL(wp), INTENT ( OUT ) :: pt !< pt (�C) corresponding to given PMV / clo IF ( pmv <= -0.11_wp ) THEN pt = 5.805_wp + 12.6784_wp * pmv ELSE IF ( pmv >= + 0.01_wp ) THEN pt = 16.826_wp + 6.163_wp * pmv ELSE pt = 21.258_wp - 9.558_wp * clo ENDIF ENDIF END SUBROUTINE pt_regression !------------------------------------------------------------------------------! ! Description: ! ------------ !> FANGER.F90 !> !> SI-VERSION: ACTLEV W m-2, DAMPFDRUCK hPa !> Berechnet das aktuelle Predicted Mean Vote nach Fanger !> !> The case of free convection (vau < 0.1 m/s) is dealt with vau = 0.1 m/s !------------------------------------------------------------------------------! SUBROUTINE fanger( tt, tmrt, pa, in_vau, pb, in_clo, actlev, eta, pmva, top ) IMPLICIT NONE !-- Input variables of argument list: REAL(wp), INTENT ( IN ) :: tt !< Ambient air temperature (�C) REAL(wp), INTENT ( IN ) :: tmrt !< Mean radiant temperature (�C) REAL(wp), INTENT ( IN ) :: pa !< Water vapour pressure (hPa) REAL(wp), INTENT ( IN ) :: pb !< Barometric pressure (hPa) at site REAL(wp), INTENT ( IN ) :: in_vau !< Wind speed (m/s) 1 m above ground REAL(wp), INTENT ( IN ) :: in_clo !< Clothing insulation (clo) REAL(wp), INTENT ( IN ) :: actlev !< Individuals activity level per unit surface area (W/m2) REAL(wp), INTENT ( IN ) :: eta !< Individuals mechanical work efficiency (dimensionless) !-- Output variables of argument list: REAL(wp), INTENT ( OUT ) :: pmva !< Actual Predicted Mean Vote (dimensionless) according ! to Fanger corresponding to meteorological (tt,tmrt,pa,vau,pb) ! and individual variables (clo, actlev, eta) REAL(wp), INTENT ( OUT ) :: top !< operative temperature (�C) !-- Internal variables REAL(wp), PARAMETER :: cels_offs = 273.15_wp !< Kelvin Celsius Offset (K) REAL(wp) :: f_cl REAL(wp) :: heat_convection REAL(wp) :: activity REAL(wp) :: t_skin_aver REAL(wp) :: bc REAL(wp) :: cc REAL(wp) :: dc REAL(wp) :: ec REAL(wp) :: gc REAL(wp) :: t_clothing REAL(wp) :: hr REAL(wp) :: clo REAL(wp) :: vau REAL(wp) :: z1 REAL(wp) :: z2 REAL(wp) :: z3 REAL(wp) :: z4 REAL(wp) :: z5 REAL(wp) :: z6 INTEGER(iwp) :: i !-- Clo must be > 0. to avoid div. by 0! clo = in_clo IF ( clo <= 0._wp ) clo = .001_wp !-- f_cl = Increase in surface due to clothing f_cl = 1._wp + .15_wp * clo !-- Case of free convection (vau < 0.1 m/s ) not considered vau = in_vau IF ( vau < .1_wp ) THEN vau = .1_wp ENDIF !-- Heat_convection = forced convection heat_convection = 12.1_wp * SQRT ( vau * pb / 1013.25_wp ) !-- Activity = inner heat produktion per standardized surface activity = actlev * ( 1._wp - eta ) !-- T_skin_aver = average skin temperature t_skin_aver = 35.7_wp - .0275_wp * activity !-- Calculation of constants for evaluation below bc = .155_wp * clo * 3.96_wp * 10._wp**( -8 ) * f_cl cc = f_cl * heat_convection ec = .155_wp * clo dc = ( 1._wp + ec * cc ) / bc gc = ( t_skin_aver + bc * ( tmrt + cels_offs )**4 + ec * cc * tt ) / bc !-- Calculation of clothing surface temperature (t_clothing) based on ! Newton-approximation with air temperature as initial guess t_clothing = tt DO I = 1, 3 t_clothing = t_clothing - ( ( t_clothing + cels_offs )**4 + t_clothing & * dc - gc ) / ( 4._wp * ( t_clothing + cels_offs )**3 + dc ) ENDDO !-- Empiric factor for the adaption of the heat ballance equation ! to the psycho-physical scale (Equ. 40 in FANGER) z1 = ( .303_wp * EXP ( -.036_wp * actlev ) + .0275_wp ) !-- Water vapour diffution through the skin z2 = .31_wp * ( 57.3_wp - .07_wp * activity-pa ) !-- Sweat evaporation from the skin surface z3 = .42_wp * ( activity - 58._wp ) !-- Loss of latent heat through respiration z4 = .0017_wp * actlev * ( 58.7_wp - pa ) + .0014_wp * actlev * & ( 34._wp - tt ) !-- Loss of radiational heat z5 = 3.96e-8_wp * f_cl * ( ( t_clothing + cels_offs )**4 - ( tmrt + & cels_offs )**4 ) IF ( ABS ( t_clothing - tmrt ) > 0._wp ) THEN hr = z5 / f_cl / ( t_clothing - tmrt ) ELSE hr = 0._wp ENDIF !-- Heat loss through forced convection cc*(t_clothing-TT) z6 = cc * ( t_clothing - tt ) !-- Predicted Mean Vote pmva = z1 * ( activity - z2 - z3 - z4 - z5 - z6 ) !-- Operative temperatur top = ( hr * tmrt + heat_convection * tt ) / ( hr + heat_convection ) END SUBROUTINE fanger !------------------------------------------------------------------------------! ! Description: ! ------------ !> For pmva > 0 and clo =0.5 the increment (deltapmv) is calculated !> that converts pmva into Gagge's et al. (1986) PMV*. !------------------------------------------------------------------------------! REAL(wp) FUNCTION deltapmv( pmva, tt2m, el2m, svp_tt, tmrt, vau1m, nerr ) IMPLICIT NONE !-- Input variables of argument list: REAL(wp), INTENT ( IN ) :: pmva !< Actual Predicted Mean Vote (PMV) according to Fanger REAL(wp), INTENT ( IN ) :: tt2m !< Ambient temperature (�C) at screen level REAL(wp), INTENT ( IN ) :: el2m !< Water vapour pressure (hPa) at screen level REAL(wp), INTENT ( IN ) :: svp_tt !< Saturation water vapour pressure (hPa) at tt2m REAL(wp), INTENT ( IN ) :: tmrt !< Mean radiant temperature (�C) at screen level REAL(wp), INTENT ( IN ) :: vau1m !< Wind speed (m/s) 1 m above ground !-- Output variables of argument list: INTEGER(iwp), INTENT ( OUT ) :: nerr !< Error status / quality flag ! 0 = o.k. ! -2 = pmva outside valid regression range ! -3 = rel. humidity set to 5 % or 95 %, respectively ! -4 = deltapmv set to avoid pmvs < 0 !-- Internal variable types: REAL(wp) :: pmv REAL(wp) :: pa_p50 REAL(wp) :: pa REAL(wp) :: apa REAL(wp) :: dapa REAL(wp) :: sqvel REAL(wp) :: ta REAL(wp) :: dtmrt REAL(wp) :: p10 REAL(wp) :: p95 REAL(wp) :: gew REAL(wp) :: gew2 REAL(wp) :: dpmv_1 REAL(wp) :: dpmv_2 REAL(wp) :: pmvs REAL(wp) :: bpmv(0:7) REAL(wp) :: bpa_p50(0:7) REAL(wp) :: bpa(0:7) REAL(wp) :: bapa(0:7) REAL(wp) :: bdapa(0:7) REAL(wp) :: bsqvel(0:7) REAL(wp) :: bta(0:7) REAL(wp) :: bdtmrt(0:7) REAL(wp) :: aconst(0:7) INTEGER(iwp) :: nreg DATA bpmv / & -0.0556602_wp, -0.1528680_wp, -0.2336104_wp, -0.2789387_wp, -0.3551048_wp,& -0.4304076_wp, -0.4884961_wp, -0.4897495_wp / DATA bpa_p50 / & -0.1607154_wp, -0.4177296_wp, -0.4120541_wp, -0.0886564_wp, +0.4285938_wp,& +0.6281256_wp, +0.5067361_wp, +0.3965169_wp / DATA bpa / & +0.0580284_wp, +0.0836264_wp, +0.1009919_wp, +0.1020777_wp, +0.0898681_wp,& +0.0839116_wp, +0.0853258_wp, +0.0866589_wp / DATA bapa / & -1.7838788_wp, -2.9306231_wp, -1.6350334_wp, +0.6211547_wp, +3.3918083_wp,& +5.5521025_wp, +8.4897418_wp, +16.6265851_wp / DATA bdapa / & +1.6752720_wp, +2.7379504_wp, +1.2940526_wp, -1.0985759_wp, -3.9054732_wp,& -6.0403012_wp, -8.9437119_wp, -17.0671201_wp / DATA bsqvel / & -0.0315598_wp, -0.0286272_wp, -0.0009228_wp, +0.0483344_wp, +0.0992366_wp,& +0.1491379_wp, +0.1951452_wp, +0.2133949_wp / DATA bta / & +0.0953986_wp, +0.1524760_wp, +0.0564241_wp, -0.0893253_wp, -0.2398868_wp,& -0.3515237_wp, -0.5095144_wp, -0.9469258_wp / DATA bdtmrt / & -0.0004672_wp, -0.0000514_wp, -0.0018037_wp, -0.0049440_wp, -0.0069036_wp,& -0.0075844_wp, -0.0079602_wp, -0.0089439_wp / DATA aconst / & +1.8686215_wp, +3.4260713_wp, +2.0116185_wp, -0.7777552_wp, -4.6715853_wp,& -7.7314281_wp, -11.7602578_wp, -23.5934198_wp / !-- Test for compliance with regression range IF ( pmva < -1.0_wp .OR. pmva > 7.0_wp ) THEN nerr = -2_iwp ELSE nerr = 0_iwp ENDIF !-- Initialise classic PMV pmv = pmva !-- Water vapour pressure of air p10 = 0.05_wp * svp_tt p95 = 1.00_wp * svp_tt IF ( el2m >= p10 .AND. el2m <= p95 ) THEN pa = el2m ELSE nerr = -3_iwp IF ( el2m < p10 ) THEN !-- Due to conditions of regress�on: r.H. >= 5 % pa = p10 ELSE !-- Due to conditions of regress�on: r.H. <= 95 % pa = p95 ENDIF ENDIF IF ( pa > 0._wp ) THEN !-- Natural logarithm of pa apa = LOG ( pa ) ELSE apa = -5._wp ENDIF !-- Ratio actual water vapour pressure to that of a r.H. of 50 % pa_p50 = 0.5_wp * svp_tt IF ( pa_p50 > 0._wp .AND. pa > 0._wp ) THEN dapa = apa - LOG ( pa_p50 ) pa_p50 = pa / pa_p50 ELSE dapa = -5._wp pa_p50 = 0._wp ENDIF !-- Square root of wind velocity IF ( vau1m >= 0._wp ) THEN sqvel = SQRT ( vau1m ) ELSE sqvel = 0._wp ENDIF !-- Air temperature ta = tt2m !-- Difference mean radiation to air temperature dtmrt = tmrt - ta !-- Select the valid regression coefficients nreg = INT ( pmv ) IF ( nreg < 0_iwp ) THEN !-- value of the FUNCTION in the case pmv <= -1 deltapmv = 0._wp RETURN ENDIF gew = MOD ( pmv, 1._wp ) IF ( gew < 0._wp ) gew = 0._wp IF ( nreg > 5_iwp ) THEN ! nreg=6 nreg = 5_iwp gew = pmv - 5._wp gew2 = pmv - 6._wp IF ( gew2 > 0_iwp ) THEN gew = ( gew - gew2 ) / gew ENDIF ENDIF !-- Regression valid for 0. <= pmv <= 6. dpmv_1 = & + bpa ( nreg ) * pa & + bpmv ( nreg ) * pmv & + bapa ( nreg ) * apa & + bta ( nreg ) * ta & + bdtmrt ( nreg ) * dtmrt & + bdapa ( nreg ) * dapa & + bsqvel ( nreg ) * sqvel & + bpa_p50 ( nreg ) * pa_p50 & + aconst ( nreg ) dpmv_2 = 0._wp IF ( nreg < 6_iwp ) THEN dpmv_2 = & + bpa ( nreg + 1_iwp ) * pa & + bpmv ( nreg + 1_iwp ) * pmv & + bapa ( nreg + 1_iwp ) * apa & + bta ( nreg + 1_iwp ) * ta & + bdtmrt ( nreg + 1_iwp ) * dtmrt & + bdapa ( nreg + 1_iwp ) * dapa & + bsqvel ( nreg + 1_iwp ) * sqvel & + bpa_p50 ( nreg + 1_iwp ) * pa_p50 & + aconst ( nreg + 1_iwp ) ENDIF !-- Calculate pmv modification deltapmv = ( 1._wp - gew ) * dpmv_1 + gew * dpmv_2 pmvs = pmva + deltapmv IF ( ( pmvs ) < 0._wp ) THEN !-- Prevent negative pmv* due to problems with clothing insulation nerr = -4_iwp IF ( pmvs > -0.11_wp ) THEN !-- Threshold from FUNCTION pt_regression for winter clothing insulation deltapmv = deltapmv + 0.11_wp ELSE !-- Set pmvs to "0" for compliance with summer clothing insulation deltapmv = -1._wp * pmva ENDIF ENDIF END FUNCTION deltapmv !------------------------------------------------------------------------------! ! Description: ! ------------ !> The subroutine "calc_sultr" returns a threshold value to perceived !> temperature allowing to decide whether the actual perceived temperature !> is linked to perecption of sultriness. The threshold values depends !> on the Fanger's classical PMV, expressed here as perceived temperature !> gtc. It is derived from all synoptic observations 1966-1995 of the sites !> 10170 Rostock-Warnemuende and 10803 Freiburg showing heat load. The !> threshold value is valid for subjects acclimatised to Central European !> climate conditions. It includes all effective influences as air !> temperature, relative humidity, mean radiant temperature, and wind !> velocity: !> - defined only under warm conditions: gtc based on 0.5 (clo) and !> gtc > 16.826 (�C) !> - undefined: sultr_res set at +99. !------------------------------------------------------------------------------! SUBROUTINE calc_sultr( gtc, dgtcm, dgtcstd, sultr_res ) IMPLICIT NONE !-- Input of the argument list: REAL(wp), INTENT ( IN ) :: gtc !< Classical perceived temperature: Base is Fanger's PMV !-- Additional output variables of argument list: REAL(wp), INTENT ( OUT ) :: dgtcm !< Mean deviation gtc (classical gt) to gt* (rational gt ! calculated based on Gagge's rational PMV*) REAL(wp), INTENT ( OUT ) :: dgtcstd !< dgtcm plus its standard deviation times a factor ! determining the significance to perceive sultriness REAL(wp), INTENT ( OUT ) :: sultr_res !-- Types of coefficients mean deviation: third order polynomial REAL(wp), PARAMETER :: dgtcka = +7.5776086_wp REAL(wp), PARAMETER :: dgtckb = -0.740603_wp REAL(wp), PARAMETER :: dgtckc = +0.0213324_wp REAL(wp), PARAMETER :: dgtckd = -0.00027797237_wp !-- Types of coefficients mean deviation plus standard deviation ! regression coefficients: third order polynomial REAL(wp), PARAMETER :: dgtcsa = +0.0268918_wp REAL(wp), PARAMETER :: dgtcsb = +0.0465957_wp REAL(wp), PARAMETER :: dgtcsc = -0.00054709752_wp REAL(wp), PARAMETER :: dgtcsd = +0.0000063714823_wp !-- Factor to mean standard deviation defining SIGNificance for ! sultriness REAL(wp), PARAMETER :: faktor = 1._wp !-- Initialise sultr_res = +99._wp dgtcm = 0._wp dgtcstd = 999999._wp IF ( gtc < 16.826_wp .OR. gtc > 56._wp ) THEN !-- Unallowed classical PMV/gtc RETURN ENDIF !-- Mean deviation dependent on gtc dgtcm = dgtcka + dgtckb * gtc + dgtckc * gtc**2._wp + dgtckd * gtc**3._wp !-- Mean deviation plus its standard deviation dgtcstd = dgtcsa + dgtcsb * gtc + dgtcsc * gtc**2._wp + dgtcsd * gtc**3._wp !-- Value of the FUNCTION sultr_res = dgtcm + faktor * dgtcstd IF ( ABS ( sultr_res ) > 99._wp ) sultr_res = +99._wp END SUBROUTINE calc_sultr !------------------------------------------------------------------------------! ! Description: ! ------------ !> Multiple linear regression to calculate an increment delta_kalt, !> to adjust Fanger's classical PMV (pmva) by Gagge's 2 node model, !> applying Fanger's convective heat transfer coefficient, hcf. !> Wind velocitiy of the reference environment is 0.10 m/s !------------------------------------------------------------------------------! SUBROUTINE dpmv_cold( pmva, tt2m, vau1m, tmrt, nerr, dpmv_cold_res ) IMPLICIT NONE !-- Type of input arguments REAL(wp), INTENT ( IN ) :: pmva !< Fanger's classical predicted mean vote REAL(wp), INTENT ( IN ) :: tt2m !< Air temperature (�C) 2 m above ground REAL(wp), INTENT ( IN ) :: vau1m !< Relative wind velocity 1 m above ground (m/s) REAL(wp), INTENT ( IN ) :: tmrt !< Mean radiant temperature (�C) !-- Type of output argument INTEGER(iwp), INTENT ( OUT ) :: nerr !< Error indicator: 0 = o.k., +1 = denominator for intersection = 0 REAL(wp), INTENT ( OUT ) :: dpmv_cold_res !< Increment to adjust pmva according to the results of Gagge's ! 2 node model depending on the input !-- Type of program variables REAL(wp) :: delta_kalt(3) REAL(wp) :: pmv_cross(2) REAL(wp) :: reg_a(3) REAL(wp) :: coeff(3,5) REAL(wp) :: r_nenner REAL(wp) :: pmvc REAL(wp) :: dtmrt REAL(wp) :: SQRT_v1m INTEGER(iwp) :: i INTEGER(iwp) :: j INTEGER(iwp) :: i_bin !-- Coefficient of the 3 regression lines ! 1:const 2:*pmvc 3:*tt2m 4:*SQRT_v1m 5:*dtmrt coeff(1,1:5) = & (/ +0.161_wp, +0.130_wp, -1.125E-03_wp, +1.106E-03_wp, -4.570E-04_wp /) coeff(2,1:5) = & (/ 0.795_wp, 0.713_wp, -8.880E-03_wp, -1.803E-03_wp, -2.816E-03_wp/) coeff(3,1:5) = & (/ +0.05761_wp, +0.458_wp, -1.829E-02_wp, -5.577E-03_wp, -1.970E-03_wp /) !-- Initialise nerr = 0_iwp dpmv_cold_res = 0._wp pmvc = pmva dtmrt = tmrt - tt2m SQRT_v1m = vau1m IF ( SQRT_v1m < 0.10_wp ) THEN SQRT_v1m = 0.10_wp ELSE SQRT_v1m = SQRT ( SQRT_v1m ) ENDIF DO i = 1, 3 delta_kalt (i) = 0._wp IF ( i < 3 ) THEN pmv_cross (i) = pmvc ENDIF ENDDO DO i = 1, 3 !-- Regression constant for given meteorological variables reg_a(i) = coeff(i, 1) + coeff(i, 3) * tt2m + coeff(i, 4) * SQRT_v1m + & coeff(i,5)*dtmrt delta_kalt(i) = reg_a(i) + coeff(i, 2) * pmvc ENDDO !-- Intersection points of regression lines in terms of Fanger's PMV DO i = 1, 2 r_nenner = coeff (i, 2) - coeff (i + 1, 2) IF ( ABS ( r_nenner ) > 0.00001_wp ) THEN pmv_cross(i) = ( reg_a (i + 1) - reg_a (i) ) / r_nenner ELSE nerr = 1_iwp RETURN ENDIF ENDDO i_bin = 3 DO i = 1, 2 IF ( pmva > pmv_cross (i) ) THEN i_bin = i EXIT ENDIF ENDDO !-- Adjust to operative temperature scaled according ! to classical PMV (Fanger) dpmv_cold_res = delta_kalt(i_bin) - dpmv_adj(pmva) END SUBROUTINE dpmv_cold !------------------------------------------------------------------------------! ! Description: ! ------------ !> Calculates the summand dpmv_adj adjusting to the operative temperature !> scaled according to classical PMV (Fanger) !> Reference environment: v_1m = 0.10 m/s, dTMRT = 0 K, r.h. = 50 % !------------------------------------------------------------------------------! REAL(wp) FUNCTION dpmv_adj( pmva ) IMPLICIT NONE REAL(wp), INTENT ( IN ) :: pmva INTEGER(iwp), PARAMETER :: n_bin = 3 INTEGER(iwp), PARAMETER :: n_ca = 0 INTEGER(iwp), PARAMETER :: n_ce = 3 REAL(wp), dimension (n_bin, n_ca:n_ce) :: coef REAL(wp) :: pmv INTEGER(iwp) :: i, i_bin ! range_1 range_2 range_3 DATA (coef(i, 0), i = 1, n_bin) /+0.0941540_wp, -0.1506620_wp, -0.0871439_wp/ DATA (coef(i, 1), i = 1, n_bin) /+0.0783162_wp, -1.0612651_wp, +0.1695040_wp/ DATA (coef(i, 2), i = 1, n_bin) /+0.1350144_wp, -1.0049144_wp, -0.0167627_wp/ DATA (coef(i, 3), i = 1, n_bin) /+0.1104037_wp, -0.2005277_wp, -0.0003230_wp/ IF ( pmva <= -2.1226_wp ) THEN i_bin = 3_iwp ELSE IF ( pmva <= -1.28_wp ) THEN i_bin = 2_iwp ELSE i_bin = 1_iwp ENDIF dpmv_adj = coef( i_bin, n_ca ) pmv = 1._wp DO i = n_ca + 1, n_ce pmv = pmv * pmva dpmv_adj = dpmv_adj + coef(i_bin, i) * pmv ENDDO RETURN END FUNCTION dpmv_adj !------------------------------------------------------------------------------! ! Description: ! ------------ !> Based on perceived temperature (gt) as input, ireq_neutral determines !> the required clothing insulation (clo) for thermally neutral conditions !> (neither body cooling nor body heating). It is related to the Klima- !> Michel activity level (134.682 W/m2). IREQ_neutral is only defined !> for gt < 10 (�C) !------------------------------------------------------------------------------! REAL(wp) FUNCTION ireq_neutral( gt, ireq_minimal, nerr ) IMPLICIT NONE !-- Type declaration of arguments REAL(wp), INTENT ( IN ) :: gt REAL(wp), INTENT ( OUT ) :: ireq_minimal INTEGER(iwp), INTENT ( OUT ) :: nerr !-- Type declaration for internal varables REAL(wp) :: top02 !-- Initialise nerr = 0_iwp !-- Convert perceived temperature from basis 0.1 m/s to basis 0.2 m/s top02 = 1.8788_wp + 0.9296_wp * gt !-- IREQ neutral conditions (thermal comfort) ireq_neutral = 1.62_wp - 0.0564_wp * top02 !-- Regression only defined for gt <= 10 (�C) IF ( ireq_neutral < 0.5_wp ) THEN IF ( ireq_neutral < 0.48_wp ) THEN nerr = 1_iwp ENDIF ireq_neutral = 0.5_wp ENDIF !-- Minimal required clothing insulation: maximal acceptable body cooling ireq_minimal = 1.26_wp - 0.0588_wp * top02 IF ( nerr > 0_iwp ) THEN ireq_minimal = ireq_neutral ENDIF RETURN END FUNCTION ireq_neutral END MODULE biometeorology_pt_mod