[1817] | 1 | !> @file land_surface_model_mod.f90 |
---|
[2000] | 2 | !------------------------------------------------------------------------------! |
---|
[1496] | 3 | ! This file is part of PALM. |
---|
| 4 | ! |
---|
[2000] | 5 | ! PALM is free software: you can redistribute it and/or modify it under the |
---|
| 6 | ! terms of the GNU General Public License as published by the Free Software |
---|
| 7 | ! Foundation, either version 3 of the License, or (at your option) any later |
---|
| 8 | ! version. |
---|
[1496] | 9 | ! |
---|
| 10 | ! PALM is distributed in the hope that it will be useful, but WITHOUT ANY |
---|
| 11 | ! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR |
---|
| 12 | ! A PARTICULAR PURPOSE. See the GNU General Public License for more details. |
---|
| 13 | ! |
---|
| 14 | ! You should have received a copy of the GNU General Public License along with |
---|
| 15 | ! PALM. If not, see <http://www.gnu.org/licenses/>. |
---|
| 16 | ! |
---|
[1818] | 17 | ! Copyright 1997-2016 Leibniz Universitaet Hannover |
---|
[2000] | 18 | !------------------------------------------------------------------------------! |
---|
[1496] | 19 | ! |
---|
| 20 | ! Current revisions: |
---|
| 21 | ! ----------------- |
---|
[1789] | 22 | ! |
---|
[2032] | 23 | ! |
---|
[1789] | 24 | ! Former revisions: |
---|
| 25 | ! ----------------- |
---|
| 26 | ! $Id: land_surface_model_mod.f90 2032 2016-10-21 15:13:51Z gronemeier $ |
---|
| 27 | ! |
---|
[2032] | 28 | ! 2031 2016-10-21 15:11:58Z knoop |
---|
| 29 | ! renamed variable rho to rho_ocean |
---|
| 30 | ! |
---|
[2001] | 31 | ! 2000 2016-08-20 18:09:15Z knoop |
---|
| 32 | ! Forced header and separation lines into 80 columns |
---|
| 33 | ! |
---|
[1979] | 34 | ! 1978 2016-07-29 12:08:31Z maronga |
---|
| 35 | ! Bugfix: initial values of pave_surface and water_surface were not set. |
---|
| 36 | ! |
---|
[1977] | 37 | ! 1976 2016-07-27 13:28:04Z maronga |
---|
| 38 | ! Parts of the code have been reformatted. Use of radiation model output is |
---|
| 39 | ! generalized and simplified. Added more output quantities due to modularization |
---|
| 40 | ! |
---|
[1973] | 41 | ! 1972 2016-07-26 07:52:02Z maronga |
---|
| 42 | ! Further modularization: output of cross sections and 3D data is now done in this |
---|
| 43 | ! module. Moreover, restart data is written and read directly within this module. |
---|
| 44 | ! |
---|
| 45 | ! |
---|
[1967] | 46 | ! 1966 2016-07-18 11:54:18Z maronga |
---|
| 47 | ! Bugfix: calculation of m_total in soil model was not set to zero at model start |
---|
| 48 | ! |
---|
[1950] | 49 | ! 1949 2016-06-17 07:19:16Z maronga |
---|
| 50 | ! Bugfix: calculation of qsws_soil_eb with precipitation = .TRUE. gave |
---|
| 51 | ! qsws_soil_eb = 0 due to a typo |
---|
| 52 | ! |
---|
[1857] | 53 | ! 1856 2016-04-13 12:56:17Z maronga |
---|
| 54 | ! Bugfix: for water surfaces, the initial water surface temperature is set equal |
---|
| 55 | ! to the intital skin temperature. Moreover, the minimum value of r_a is now |
---|
| 56 | ! 1.0 to avoid too large fluxes at the first model time step |
---|
| 57 | ! |
---|
[1851] | 58 | ! 1849 2016-04-08 11:33:18Z hoffmann |
---|
| 59 | ! prr moved to arrays_3d |
---|
[1852] | 60 | ! |
---|
[1827] | 61 | ! 1826 2016-04-07 12:01:39Z maronga |
---|
| 62 | ! Cleanup after modularization |
---|
| 63 | ! |
---|
[1818] | 64 | ! 1817 2016-04-06 15:44:20Z maronga |
---|
| 65 | ! Added interface for lsm_init_arrays. Added subroutines for check_parameters, |
---|
| 66 | ! header, and parin. Renamed some subroutines. |
---|
| 67 | ! |
---|
[1789] | 68 | ! 1788 2016-03-10 11:01:04Z maronga |
---|
[1788] | 69 | ! Bugfix: calculate lambda_surface based on temperature gradient between skin |
---|
| 70 | ! layer and soil layer instead of Obukhov length |
---|
| 71 | ! Changed: moved calculation of surface specific humidity to energy balance solver |
---|
| 72 | ! New: water surfaces are available by using a fixed sea surface temperature. |
---|
| 73 | ! The roughness lengths are calculated dynamically using the Charnock |
---|
| 74 | ! parameterization. This involves the new roughness length for moisture z0q. |
---|
| 75 | ! New: modified solution of the energy balance solver and soil model for |
---|
| 76 | ! paved surfaces (i.e. asphalt concrete). |
---|
| 77 | ! Syntax layout improved. |
---|
| 78 | ! Changed: parameter dewfall removed. |
---|
[1758] | 79 | ! |
---|
[1784] | 80 | ! 1783 2016-03-06 18:36:17Z raasch |
---|
| 81 | ! netcdf variables moved to netcdf module |
---|
| 82 | ! |
---|
[1758] | 83 | ! 1757 2016-02-22 15:49:32Z maronga |
---|
[1757] | 84 | ! Bugfix: set tm_soil_m to zero after allocation. Added parameter |
---|
| 85 | ! unscheduled_radiation_calls to control calls of the radiation model based on |
---|
| 86 | ! the skin temperature change during one time step (preliminary version). Set |
---|
| 87 | ! qsws_soil_eb to zero at model start (previously set to qsws_eb). Removed MAX |
---|
| 88 | ! function as it cannot be vectorized. |
---|
[1710] | 89 | ! |
---|
| 90 | ! 1709 2015-11-04 14:47:01Z maronga |
---|
[1709] | 91 | ! Renamed pt_1 and qv_1 to pt1 and qv1. |
---|
| 92 | ! Bugfix: set initial values for t_surface_p in case of restart runs |
---|
| 93 | ! Bugfix: zero resistance caused crash when using radiation_scheme = 'clear-sky' |
---|
| 94 | ! Bugfix: calculation of rad_net when using radiation_scheme = 'clear-sky' |
---|
| 95 | ! Added todo action |
---|
[1698] | 96 | ! |
---|
| 97 | ! 1697 2015-10-28 17:14:10Z raasch |
---|
| 98 | ! bugfix: misplaced cpp-directive |
---|
| 99 | ! |
---|
| 100 | ! 1695 2015-10-27 10:03:11Z maronga |
---|
[1692] | 101 | ! Bugfix: REAL constants provided with KIND-attribute in call of |
---|
[1696] | 102 | ! Replaced rif with ol |
---|
| 103 | ! |
---|
[1692] | 104 | ! 1691 2015-10-26 16:17:44Z maronga |
---|
[1691] | 105 | ! Added skip_time_do_lsm to allow for spin-ups without LSM. Various bugfixes: |
---|
| 106 | ! Soil temperatures are now defined at the edges of the layers, calculation of |
---|
| 107 | ! shb_eb corrected, prognostic equation for skin temperature corrected. Surface |
---|
| 108 | ! fluxes are now directly transfered to atmosphere |
---|
[1552] | 109 | ! |
---|
[1683] | 110 | ! 1682 2015-10-07 23:56:08Z knoop |
---|
| 111 | ! Code annotations made doxygen readable |
---|
| 112 | ! |
---|
[1591] | 113 | ! 1590 2015-05-08 13:56:27Z maronga |
---|
| 114 | ! Bugfix: definition of character strings requires same length for all elements |
---|
| 115 | ! |
---|
[1586] | 116 | ! 1585 2015-04-30 07:05:52Z maronga |
---|
| 117 | ! Modifications for RRTMG. Changed tables to PARAMETER type. |
---|
| 118 | ! |
---|
[1572] | 119 | ! 1571 2015-03-12 16:12:49Z maronga |
---|
| 120 | ! Removed upper-case variable names. Corrected distribution of precipitation to |
---|
| 121 | ! the liquid water reservoir and the bare soil fractions. |
---|
| 122 | ! |
---|
[1556] | 123 | ! 1555 2015-03-04 17:44:27Z maronga |
---|
| 124 | ! Added output of r_a and r_s |
---|
| 125 | ! |
---|
[1554] | 126 | ! 1553 2015-03-03 17:33:54Z maronga |
---|
| 127 | ! Improved better treatment of roughness lengths. Added default soil temperature |
---|
| 128 | ! profile |
---|
| 129 | ! |
---|
[1552] | 130 | ! 1551 2015-03-03 14:18:16Z maronga |
---|
[1551] | 131 | ! Flux calculation is now done in prandtl_fluxes. Added support for data output. |
---|
| 132 | ! Vertical indices have been replaced. Restart runs are now possible. Some |
---|
| 133 | ! variables have beem renamed. Bugfix in the prognostic equation for the surface |
---|
| 134 | ! temperature. Introduced z0_eb and z0h_eb, which overwrite the setting of |
---|
| 135 | ! roughness_length and z0_factor. Added Clapp & Hornberger parametrization for |
---|
| 136 | ! the hydraulic conductivity. Bugfix for root fraction and extraction |
---|
| 137 | ! calculation |
---|
[1496] | 138 | ! |
---|
[1514] | 139 | ! intrinsic function MAX and MIN |
---|
[1496] | 140 | ! |
---|
[1501] | 141 | ! 1500 2014-12-03 17:42:41Z maronga |
---|
| 142 | ! Corrected calculation of aerodynamic resistance (r_a). |
---|
| 143 | ! Precipitation is now added to liquid water reservoir using LE_liq. |
---|
| 144 | ! Added support for dry runs. |
---|
| 145 | ! |
---|
[1497] | 146 | ! 1496 2014-12-02 17:25:50Z maronga |
---|
| 147 | ! Initial revision |
---|
| 148 | ! |
---|
[1496] | 149 | ! |
---|
| 150 | ! Description: |
---|
| 151 | ! ------------ |
---|
[1682] | 152 | !> Land surface model, consisting of a solver for the energy balance at the |
---|
| 153 | !> surface and a four layer soil scheme. The scheme is similar to the TESSEL |
---|
| 154 | !> scheme implemented in the ECMWF IFS model, with modifications according to |
---|
| 155 | !> H-TESSEL. The implementation is based on the formulation implemented in the |
---|
| 156 | !> DALES and UCLA-LES models. |
---|
| 157 | !> |
---|
[1691] | 158 | !> @todo Consider partial absorption of the net shortwave radiation by the |
---|
[1709] | 159 | !> skin layer. |
---|
[1788] | 160 | !> @todo Improve surface water parameterization |
---|
[1691] | 161 | !> @todo Invert indices (running from -3 to 0. Currently: nzb_soil=0, |
---|
| 162 | !> nzt_soil=3)). |
---|
| 163 | !> @todo Implement surface runoff model (required when performing long-term LES |
---|
| 164 | !> with considerable precipitation. |
---|
[1709] | 165 | !> @todo Fix crashes with radiation_scheme == 'constant' |
---|
[1682] | 166 | !> |
---|
[1709] | 167 | !> @note No time step criterion is required as long as the soil layers do not |
---|
| 168 | !> become too thin. |
---|
[1496] | 169 | !------------------------------------------------------------------------------! |
---|
[1682] | 170 | MODULE land_surface_model_mod |
---|
| 171 | |
---|
[1691] | 172 | USE arrays_3d, & |
---|
[1849] | 173 | ONLY: hyp, ol, pt, pt_p, prr, q, q_p, ql, qsws, shf, ts, us, vpt, z0, & |
---|
| 174 | z0h, z0q |
---|
[1496] | 175 | |
---|
[1691] | 176 | USE cloud_parameters, & |
---|
[1849] | 177 | ONLY: cp, hyrho, l_d_cp, l_d_r, l_v, pt_d_t, rho_l, r_d, r_v |
---|
[1496] | 178 | |
---|
[1691] | 179 | USE control_parameters, & |
---|
| 180 | ONLY: cloud_physics, dt_3d, humidity, intermediate_timestep_count, & |
---|
| 181 | initializing_actions, intermediate_timestep_count_max, & |
---|
[1817] | 182 | max_masks, precipitation, pt_surface, & |
---|
| 183 | rho_surface, roughness_length, surface_pressure, & |
---|
| 184 | timestep_scheme, tsc, z0h_factor, time_since_reference_point |
---|
[1496] | 185 | |
---|
[1691] | 186 | USE indices, & |
---|
| 187 | ONLY: nbgp, nxlg, nxrg, nyng, nysg, nzb, nzb_s_inner |
---|
[1496] | 188 | |
---|
[1691] | 189 | USE kinds |
---|
[1496] | 190 | |
---|
[1691] | 191 | USE pegrid |
---|
[1496] | 192 | |
---|
[1691] | 193 | USE radiation_model_mod, & |
---|
[1976] | 194 | ONLY: force_radiation_call, rad_net, rad_sw_in, rad_lw_out, & |
---|
| 195 | rad_lw_out_change_0, unscheduled_radiation_calls |
---|
[1691] | 196 | |
---|
| 197 | USE statistics, & |
---|
| 198 | ONLY: hom, statistic_regions |
---|
| 199 | |
---|
[1496] | 200 | IMPLICIT NONE |
---|
| 201 | |
---|
| 202 | ! |
---|
| 203 | !-- LSM model constants |
---|
[1682] | 204 | INTEGER(iwp), PARAMETER :: nzb_soil = 0, & !< bottom of the soil model (to be switched) |
---|
| 205 | nzt_soil = 3, & !< top of the soil model (to be switched) |
---|
| 206 | nzs = 4 !< number of soil layers (fixed for now) |
---|
[1496] | 207 | |
---|
| 208 | REAL(wp), PARAMETER :: & |
---|
[1551] | 209 | b_ch = 6.04_wp, & ! Clapp & Hornberger exponent |
---|
| 210 | lambda_h_dry = 0.19_wp, & ! heat conductivity for dry soil |
---|
[1496] | 211 | lambda_h_sm = 3.44_wp, & ! heat conductivity of the soil matrix |
---|
| 212 | lambda_h_water = 0.57_wp, & ! heat conductivity of water |
---|
| 213 | psi_sat = -0.388_wp, & ! soil matrix potential at saturation |
---|
[1551] | 214 | rho_c_soil = 2.19E6_wp, & ! volumetric heat capacity of soil |
---|
| 215 | rho_c_water = 4.20E6_wp, & ! volumetric heat capacity of water |
---|
[1496] | 216 | m_max_depth = 0.0002_wp ! Maximum capacity of the water reservoir (m) |
---|
| 217 | |
---|
| 218 | |
---|
| 219 | ! |
---|
| 220 | !-- LSM variables |
---|
[1788] | 221 | INTEGER(iwp) :: veg_type = 2, & !< NAMELIST veg_type_2d |
---|
| 222 | soil_type = 3 !< NAMELIST soil_type_2d |
---|
[1496] | 223 | |
---|
[1788] | 224 | INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: soil_type_2d, & !< soil type, 0: user-defined, 1-7: generic (see list) |
---|
| 225 | veg_type_2d !< vegetation type, 0: user-defined, 1-19: generic (see list) |
---|
| 226 | |
---|
| 227 | LOGICAL, DIMENSION(:,:), ALLOCATABLE :: water_surface, & !< flag parameter for water surfaces (classes 14+15) |
---|
| 228 | pave_surface, & !< flag parameter for pavements (asphalt etc.) (class 20) |
---|
| 229 | building_surface !< flag parameter indicating that the surface element is covered by buildings (no LSM actions, not implemented yet) |
---|
| 230 | |
---|
[1691] | 231 | LOGICAL :: conserve_water_content = .TRUE., & !< open or closed bottom surface for the soil model |
---|
| 232 | force_radiation_call_l = .FALSE., & !< flag parameter for unscheduled radiation model calls |
---|
| 233 | land_surface = .FALSE. !< flag parameter indicating wheather the lsm is used |
---|
[1496] | 234 | |
---|
| 235 | ! value 9999999.9_wp -> generic available or user-defined value must be set |
---|
| 236 | ! otherwise -> no generic variable and user setting is optional |
---|
[1682] | 237 | REAL(wp) :: alpha_vangenuchten = 9999999.9_wp, & !< NAMELIST alpha_vg |
---|
| 238 | canopy_resistance_coefficient = 9999999.9_wp, & !< NAMELIST g_d |
---|
| 239 | c_surface = 20000.0_wp, & !< Surface (skin) heat capacity |
---|
| 240 | drho_l_lv, & !< (rho_l * l_v)**-1 |
---|
| 241 | exn, & !< value of the Exner function |
---|
| 242 | e_s = 0.0_wp, & !< saturation water vapour pressure |
---|
| 243 | field_capacity = 9999999.9_wp, & !< NAMELIST m_fc |
---|
| 244 | f_shortwave_incoming = 9999999.9_wp, & !< NAMELIST f_sw_in |
---|
| 245 | hydraulic_conductivity = 9999999.9_wp, & !< NAMELIST gamma_w_sat |
---|
| 246 | ke = 0.0_wp, & !< Kersten number |
---|
[1691] | 247 | lambda_h_sat = 0.0_wp, & !< heat conductivity for saturated soil |
---|
[1682] | 248 | lambda_surface_stable = 9999999.9_wp, & !< NAMELIST lambda_surface_s |
---|
| 249 | lambda_surface_unstable = 9999999.9_wp, & !< NAMELIST lambda_surface_u |
---|
| 250 | leaf_area_index = 9999999.9_wp, & !< NAMELIST lai |
---|
| 251 | l_vangenuchten = 9999999.9_wp, & !< NAMELIST l_vg |
---|
| 252 | min_canopy_resistance = 9999999.9_wp, & !< NAMELIST r_canopy_min |
---|
| 253 | min_soil_resistance = 50.0_wp, & !< NAMELIST r_soil_min |
---|
| 254 | m_total = 0.0_wp, & !< weighted total water content of the soil (m3/m3) |
---|
| 255 | n_vangenuchten = 9999999.9_wp, & !< NAMELIST n_vg |
---|
[1788] | 256 | pave_depth = 9999999.9_wp, & !< depth of the pavement |
---|
| 257 | pave_heat_capacity = 1.94E6_wp, & !< volumetric heat capacity of pavement (e.g. roads) |
---|
| 258 | pave_heat_conductivity = 1.00_wp, & !< heat conductivity for pavements (e.g. roads) |
---|
[1682] | 259 | q_s = 0.0_wp, & !< saturation specific humidity |
---|
| 260 | residual_moisture = 9999999.9_wp, & !< NAMELIST m_res |
---|
| 261 | rho_cp, & !< rho_surface * cp |
---|
[2031] | 262 | rho_lv, & !< rho_ocean * l_v |
---|
[1682] | 263 | rd_d_rv, & !< r_d / r_v |
---|
| 264 | saturation_moisture = 9999999.9_wp, & !< NAMELIST m_sat |
---|
[1691] | 265 | skip_time_do_lsm = 0.0_wp, & !< LSM is not called before this time |
---|
[1682] | 266 | vegetation_coverage = 9999999.9_wp, & !< NAMELIST c_veg |
---|
| 267 | wilting_point = 9999999.9_wp, & !< NAMELIST m_wilt |
---|
| 268 | z0_eb = 9999999.9_wp, & !< NAMELIST z0 (lsm_par) |
---|
[1788] | 269 | z0h_eb = 9999999.9_wp, & !< NAMELIST z0h (lsm_par) |
---|
| 270 | z0q_eb = 9999999.9_wp !< NAMELIST z0q (lsm_par) |
---|
[1496] | 271 | |
---|
[1551] | 272 | REAL(wp), DIMENSION(nzb_soil:nzt_soil) :: & |
---|
[1691] | 273 | ddz_soil_stag, & !< 1/dz_soil_stag |
---|
| 274 | dz_soil_stag, & !< soil grid spacing (center-center) |
---|
| 275 | root_extr = 0.0_wp, & !< root extraction |
---|
[1551] | 276 | root_fraction = (/9999999.9_wp, 9999999.9_wp, & |
---|
[1682] | 277 | 9999999.9_wp, 9999999.9_wp /), & !< distribution of root surface area to the individual soil layers |
---|
| 278 | zs = (/0.07_wp, 0.28_wp, 1.00_wp, 2.89_wp/), & !< soil layer depths (m) |
---|
| 279 | soil_moisture = 0.0_wp !< soil moisture content (m3/m3) |
---|
[1496] | 280 | |
---|
[1551] | 281 | REAL(wp), DIMENSION(nzb_soil:nzt_soil+1) :: & |
---|
[1691] | 282 | soil_temperature = (/290.0_wp, 287.0_wp, 285.0_wp, 283.0_wp, & !< soil temperature (K) |
---|
| 283 | 283.0_wp /), & |
---|
| 284 | ddz_soil, & !< 1/dz_soil |
---|
| 285 | dz_soil !< soil grid spacing (edge-edge) |
---|
[1496] | 286 | |
---|
| 287 | #if defined( __nopointer ) |
---|
[1682] | 288 | REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET :: t_surface, & !< surface temperature (K) |
---|
| 289 | t_surface_p, & !< progn. surface temperature (K) |
---|
| 290 | m_liq_eb, & !< liquid water reservoir (m) |
---|
| 291 | m_liq_eb_av, & !< liquid water reservoir (m) |
---|
| 292 | m_liq_eb_p !< progn. liquid water reservoir (m) |
---|
[1496] | 293 | #else |
---|
[1551] | 294 | REAL(wp), DIMENSION(:,:), POINTER :: t_surface, & |
---|
| 295 | t_surface_p, & |
---|
| 296 | m_liq_eb, & |
---|
| 297 | m_liq_eb_p |
---|
[1496] | 298 | |
---|
[1551] | 299 | REAL(wp), DIMENSION(:,:), ALLOCATABLE, TARGET :: t_surface_1, t_surface_2, & |
---|
| 300 | m_liq_eb_av, & |
---|
| 301 | m_liq_eb_1, m_liq_eb_2 |
---|
[1496] | 302 | #endif |
---|
| 303 | |
---|
| 304 | ! |
---|
| 305 | !-- Temporal tendencies for time stepping |
---|
[1682] | 306 | REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tt_surface_m, & !< surface temperature tendency (K) |
---|
| 307 | tm_liq_eb_m !< liquid water reservoir tendency (m) |
---|
[1496] | 308 | |
---|
| 309 | ! |
---|
| 310 | !-- Energy balance variables |
---|
[1691] | 311 | REAL(wp), DIMENSION(:,:), ALLOCATABLE :: & |
---|
[1682] | 312 | alpha_vg, & !< coef. of Van Genuchten |
---|
| 313 | c_liq, & !< liquid water coverage (of vegetated area) |
---|
| 314 | c_liq_av, & !< average of c_liq |
---|
| 315 | c_soil_av, & !< average of c_soil |
---|
| 316 | c_veg, & !< vegetation coverage |
---|
| 317 | c_veg_av, & !< average of c_veg |
---|
| 318 | f_sw_in, & !< fraction of absorbed shortwave radiation by the surface layer (not implemented yet) |
---|
| 319 | ghf_eb, & !< ground heat flux |
---|
| 320 | ghf_eb_av, & !< average of ghf_eb |
---|
| 321 | gamma_w_sat, & !< hydraulic conductivity at saturation |
---|
| 322 | g_d, & !< coefficient for dependence of r_canopy on water vapour pressure deficit |
---|
| 323 | lai, & !< leaf area index |
---|
| 324 | lai_av, & !< average of lai |
---|
[1691] | 325 | lambda_surface_s, & !< coupling between surface and soil (depends on vegetation type) |
---|
| 326 | lambda_surface_u, & !< coupling between surface and soil (depends on vegetation type) |
---|
[1682] | 327 | l_vg, & !< coef. of Van Genuchten |
---|
| 328 | m_fc, & !< soil moisture at field capacity (m3/m3) |
---|
| 329 | m_res, & !< residual soil moisture |
---|
| 330 | m_sat, & !< saturation soil moisture (m3/m3) |
---|
| 331 | m_wilt, & !< soil moisture at permanent wilting point (m3/m3) |
---|
| 332 | n_vg, & !< coef. Van Genuchten |
---|
| 333 | qsws_eb, & !< surface flux of latent heat (total) |
---|
| 334 | qsws_eb_av, & !< average of qsws_eb |
---|
| 335 | qsws_liq_eb, & !< surface flux of latent heat (liquid water portion) |
---|
| 336 | qsws_liq_eb_av, & !< average of qsws_liq_eb |
---|
| 337 | qsws_soil_eb, & !< surface flux of latent heat (soil portion) |
---|
| 338 | qsws_soil_eb_av, & !< average of qsws_soil_eb |
---|
| 339 | qsws_veg_eb, & !< surface flux of latent heat (vegetation portion) |
---|
| 340 | qsws_veg_eb_av, & !< average of qsws_veg_eb |
---|
| 341 | rad_net_l, & !< local copy of rad_net (net radiation at surface) |
---|
| 342 | r_a, & !< aerodynamic resistance |
---|
[1709] | 343 | r_a_av, & !< average of r_a |
---|
[1682] | 344 | r_canopy, & !< canopy resistance |
---|
[1691] | 345 | r_soil, & !< soil resistance |
---|
[1682] | 346 | r_soil_min, & !< minimum soil resistance |
---|
| 347 | r_s, & !< total surface resistance (combination of r_soil and r_canopy) |
---|
[1691] | 348 | r_s_av, & !< average of r_s |
---|
[1682] | 349 | r_canopy_min, & !< minimum canopy (stomatal) resistance |
---|
| 350 | shf_eb, & !< surface flux of sensible heat |
---|
| 351 | shf_eb_av !< average of shf_eb |
---|
[1496] | 352 | |
---|
[1691] | 353 | |
---|
[1496] | 354 | REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: & |
---|
[1788] | 355 | lambda_h, & !< heat conductivity of soil (W/m/K) |
---|
[1682] | 356 | lambda_w, & !< hydraulic diffusivity of soil (?) |
---|
[1788] | 357 | gamma_w, & !< hydraulic conductivity of soil (W/m/K) |
---|
[1682] | 358 | rho_c_total !< volumetric heat capacity of the actual soil matrix (?) |
---|
[1496] | 359 | |
---|
| 360 | #if defined( __nopointer ) |
---|
| 361 | REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: & |
---|
[1682] | 362 | t_soil, & !< Soil temperature (K) |
---|
| 363 | t_soil_av, & !< Average of t_soil |
---|
| 364 | t_soil_p, & !< Prog. soil temperature (K) |
---|
| 365 | m_soil, & !< Soil moisture (m3/m3) |
---|
| 366 | m_soil_av, & !< Average of m_soil |
---|
| 367 | m_soil_p !< Prog. soil moisture (m3/m3) |
---|
[1496] | 368 | #else |
---|
| 369 | REAL(wp), DIMENSION(:,:,:), POINTER :: & |
---|
[1551] | 370 | t_soil, t_soil_p, & |
---|
[1496] | 371 | m_soil, m_soil_p |
---|
| 372 | |
---|
| 373 | REAL(wp), DIMENSION(:,:,:), ALLOCATABLE, TARGET :: & |
---|
[1551] | 374 | t_soil_av, t_soil_1, t_soil_2, & |
---|
| 375 | m_soil_av, m_soil_1, m_soil_2 |
---|
[1496] | 376 | #endif |
---|
| 377 | |
---|
| 378 | |
---|
| 379 | REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: & |
---|
[1757] | 380 | tt_soil_m, & !< t_soil storage array |
---|
| 381 | tm_soil_m, & !< m_soil storage array |
---|
[1682] | 382 | root_fr !< root fraction (sum=1) |
---|
[1496] | 383 | |
---|
[1551] | 384 | |
---|
[1496] | 385 | ! |
---|
[1551] | 386 | !-- Predefined Land surface classes (veg_type) |
---|
[1788] | 387 | CHARACTER(26), DIMENSION(0:20), PARAMETER :: veg_type_name = (/ & |
---|
[1590] | 388 | 'user defined ', & ! 0 |
---|
| 389 | 'crops, mixed farming ', & ! 1 |
---|
| 390 | 'short grass ', & ! 2 |
---|
[1585] | 391 | 'evergreen needleleaf trees', & ! 3 |
---|
| 392 | 'deciduous needleleaf trees', & ! 4 |
---|
[1590] | 393 | 'evergreen broadleaf trees ', & ! 5 |
---|
| 394 | 'deciduous broadleaf trees ', & ! 6 |
---|
| 395 | 'tall grass ', & ! 7 |
---|
| 396 | 'desert ', & ! 8 |
---|
| 397 | 'tundra ', & ! 9 |
---|
| 398 | 'irrigated crops ', & ! 10 |
---|
| 399 | 'semidesert ', & ! 11 |
---|
| 400 | 'ice caps and glaciers ', & ! 12 |
---|
| 401 | 'bogs and marshes ', & ! 13 |
---|
| 402 | 'inland water ', & ! 14 |
---|
| 403 | 'ocean ', & ! 15 |
---|
| 404 | 'evergreen shrubs ', & ! 16 |
---|
| 405 | 'deciduous shrubs ', & ! 17 |
---|
| 406 | 'mixed forest/woodland ', & ! 18 |
---|
[1788] | 407 | 'interrupted forest ', & ! 19 |
---|
| 408 | 'pavements/roads ' & ! 20 |
---|
[1585] | 409 | /) |
---|
[1496] | 410 | |
---|
| 411 | ! |
---|
[1551] | 412 | !-- Soil model classes (soil_type) |
---|
[1585] | 413 | CHARACTER(12), DIMENSION(0:7), PARAMETER :: soil_type_name = (/ & |
---|
| 414 | 'user defined', & ! 0 |
---|
[1590] | 415 | 'coarse ', & ! 1 |
---|
| 416 | 'medium ', & ! 2 |
---|
| 417 | 'medium-fine ', & ! 3 |
---|
| 418 | 'fine ', & ! 4 |
---|
| 419 | 'very fine ', & ! 5 |
---|
| 420 | 'organic ', & ! 6 |
---|
| 421 | 'loamy (CH) ' & ! 7 |
---|
[1585] | 422 | /) |
---|
[1551] | 423 | ! |
---|
| 424 | !-- Land surface parameters according to the respective classes (veg_type) |
---|
| 425 | |
---|
| 426 | ! |
---|
| 427 | !-- Land surface parameters I |
---|
| 428 | !-- r_canopy_min, lai, c_veg, g_d |
---|
[1788] | 429 | REAL(wp), DIMENSION(0:3,1:20), PARAMETER :: veg_pars = RESHAPE( (/ & |
---|
[1585] | 430 | 180.0_wp, 3.00_wp, 0.90_wp, 0.00_wp, & ! 1 |
---|
| 431 | 110.0_wp, 2.00_wp, 0.85_wp, 0.00_wp, & ! 2 |
---|
| 432 | 500.0_wp, 5.00_wp, 0.90_wp, 0.03_wp, & ! 3 |
---|
| 433 | 500.0_wp, 5.00_wp, 0.90_wp, 0.03_wp, & ! 4 |
---|
| 434 | 175.0_wp, 5.00_wp, 0.90_wp, 0.03_wp, & ! 5 |
---|
| 435 | 240.0_wp, 6.00_wp, 0.99_wp, 0.13_wp, & ! 6 |
---|
| 436 | 100.0_wp, 2.00_wp, 0.70_wp, 0.00_wp, & ! 7 |
---|
| 437 | 250.0_wp, 0.05_wp, 0.00_wp, 0.00_wp, & ! 8 |
---|
| 438 | 80.0_wp, 1.00_wp, 0.50_wp, 0.00_wp, & ! 9 |
---|
| 439 | 180.0_wp, 3.00_wp, 0.90_wp, 0.00_wp, & ! 10 |
---|
| 440 | 150.0_wp, 0.50_wp, 0.10_wp, 0.00_wp, & ! 11 |
---|
| 441 | 0.0_wp, 0.00_wp, 0.00_wp, 0.00_wp, & ! 12 |
---|
| 442 | 240.0_wp, 4.00_wp, 0.60_wp, 0.00_wp, & ! 13 |
---|
| 443 | 0.0_wp, 0.00_wp, 0.00_wp, 0.00_wp, & ! 14 |
---|
| 444 | 0.0_wp, 0.00_wp, 0.00_wp, 0.00_wp, & ! 15 |
---|
| 445 | 225.0_wp, 3.00_wp, 0.50_wp, 0.00_wp, & ! 16 |
---|
| 446 | 225.0_wp, 1.50_wp, 0.50_wp, 0.00_wp, & ! 17 |
---|
| 447 | 250.0_wp, 5.00_wp, 0.90_wp, 0.03_wp, & ! 18 |
---|
[1788] | 448 | 175.0_wp, 2.50_wp, 0.90_wp, 0.03_wp, & ! 19 |
---|
| 449 | 0.0_wp, 0.00_wp, 0.00_wp, 0.00_wp & ! 20 |
---|
| 450 | /), (/ 4, 20 /) ) |
---|
[1496] | 451 | |
---|
| 452 | ! |
---|
[1788] | 453 | !-- Land surface parameters II z0, z0h, z0q |
---|
| 454 | REAL(wp), DIMENSION(0:2,1:20), PARAMETER :: roughness_par = RESHAPE( (/ & |
---|
| 455 | 0.25_wp, 0.25E-2_wp, 0.25E-2_wp, & ! 1 |
---|
| 456 | 0.20_wp, 0.20E-2_wp, 0.20E-2_wp, & ! 2 |
---|
| 457 | 2.00_wp, 2.00_wp, 2.00_wp, & ! 3 |
---|
| 458 | 2.00_wp, 2.00_wp, 2.00_wp, & ! 4 |
---|
| 459 | 2.00_wp, 2.00_wp, 2.00_wp, & ! 5 |
---|
| 460 | 2.00_wp, 2.00_wp, 2.00_wp, & ! 6 |
---|
| 461 | 0.47_wp, 0.47E-2_wp, 0.47E-2_wp, & ! 7 |
---|
| 462 | 0.013_wp, 0.013E-2_wp, 0.013E-2_wp, & ! 8 |
---|
| 463 | 0.034_wp, 0.034E-2_wp, 0.034E-2_wp, & ! 9 |
---|
| 464 | 0.5_wp, 0.50E-2_wp, 0.50E-2_wp, & ! 10 |
---|
| 465 | 0.17_wp, 0.17E-2_wp, 0.17E-2_wp, & ! 11 |
---|
| 466 | 1.3E-3_wp, 1.3E-4_wp, 1.3E-4_wp, & ! 12 |
---|
| 467 | 0.83_wp, 0.83E-2_wp, 0.83E-2_wp, & ! 13 |
---|
| 468 | 0.00_wp, 0.00_wp, 0.00_wp, & ! 14 |
---|
| 469 | 0.00_wp, 0.00_wp, 0.00_wp, & ! 15 |
---|
| 470 | 0.10_wp, 0.10E-2_wp, 0.10E-2_wp, & ! 16 |
---|
| 471 | 0.25_wp, 0.25E-2_wp, 0.25E-2_wp, & ! 17 |
---|
| 472 | 2.00_wp, 2.00E-2_wp, 2.00E-2_wp, & ! 18 |
---|
| 473 | 1.10_wp, 1.10E-2_wp, 1.10E-2_wp, & ! 19 |
---|
| 474 | 1.0E-4_wp, 1.0E-5_wp, 1.0E-5_wp & ! 20 |
---|
| 475 | /), (/ 3, 20 /) ) |
---|
[1496] | 476 | |
---|
| 477 | ! |
---|
[1551] | 478 | !-- Land surface parameters III lambda_surface_s, lambda_surface_u, f_sw_in |
---|
[1788] | 479 | REAL(wp), DIMENSION(0:2,1:20), PARAMETER :: surface_pars = RESHAPE( (/ & |
---|
[1585] | 480 | 10.0_wp, 10.0_wp, 0.05_wp, & ! 1 |
---|
| 481 | 10.0_wp, 10.0_wp, 0.05_wp, & ! 2 |
---|
| 482 | 20.0_wp, 15.0_wp, 0.03_wp, & ! 3 |
---|
| 483 | 20.0_wp, 15.0_wp, 0.03_wp, & ! 4 |
---|
| 484 | 20.0_wp, 15.0_wp, 0.03_wp, & ! 5 |
---|
| 485 | 20.0_wp, 15.0_wp, 0.03_wp, & ! 6 |
---|
| 486 | 10.0_wp, 10.0_wp, 0.05_wp, & ! 7 |
---|
| 487 | 15.0_wp, 15.0_wp, 0.00_wp, & ! 8 |
---|
| 488 | 10.0_wp, 10.0_wp, 0.05_wp, & ! 9 |
---|
| 489 | 10.0_wp, 10.0_wp, 0.05_wp, & ! 10 |
---|
| 490 | 10.0_wp, 10.0_wp, 0.05_wp, & ! 11 |
---|
| 491 | 58.0_wp, 58.0_wp, 0.00_wp, & ! 12 |
---|
| 492 | 10.0_wp, 10.0_wp, 0.05_wp, & ! 13 |
---|
[1788] | 493 | 1.0E10_wp, 1.0E10_wp, 0.00_wp, & ! 14 |
---|
| 494 | 1.0E10_wp, 1.0E10_wp, 0.00_wp, & ! 15 |
---|
[1585] | 495 | 10.0_wp, 10.0_wp, 0.05_wp, & ! 16 |
---|
| 496 | 10.0_wp, 10.0_wp, 0.05_wp, & ! 17 |
---|
| 497 | 20.0_wp, 15.0_wp, 0.03_wp, & ! 18 |
---|
[1788] | 498 | 20.0_wp, 15.0_wp, 0.03_wp, & ! 19 |
---|
| 499 | 0.0_wp, 0.0_wp, 0.00_wp & ! 20 |
---|
| 500 | /), (/ 3, 20 /) ) |
---|
[1496] | 501 | |
---|
| 502 | ! |
---|
| 503 | !-- Root distribution (sum = 1) level 1, level 2, level 3, level 4, |
---|
[1788] | 504 | REAL(wp), DIMENSION(0:3,1:20), PARAMETER :: root_distribution = RESHAPE( (/ & |
---|
[1585] | 505 | 0.24_wp, 0.41_wp, 0.31_wp, 0.04_wp, & ! 1 |
---|
| 506 | 0.35_wp, 0.38_wp, 0.23_wp, 0.04_wp, & ! 2 |
---|
| 507 | 0.26_wp, 0.39_wp, 0.29_wp, 0.06_wp, & ! 3 |
---|
| 508 | 0.26_wp, 0.38_wp, 0.29_wp, 0.07_wp, & ! 4 |
---|
| 509 | 0.24_wp, 0.38_wp, 0.31_wp, 0.07_wp, & ! 5 |
---|
| 510 | 0.25_wp, 0.34_wp, 0.27_wp, 0.14_wp, & ! 6 |
---|
| 511 | 0.27_wp, 0.27_wp, 0.27_wp, 0.09_wp, & ! 7 |
---|
| 512 | 1.00_wp, 0.00_wp, 0.00_wp, 0.00_wp, & ! 8 |
---|
| 513 | 0.47_wp, 0.45_wp, 0.08_wp, 0.00_wp, & ! 9 |
---|
| 514 | 0.24_wp, 0.41_wp, 0.31_wp, 0.04_wp, & ! 10 |
---|
| 515 | 0.17_wp, 0.31_wp, 0.33_wp, 0.19_wp, & ! 11 |
---|
| 516 | 0.00_wp, 0.00_wp, 0.00_wp, 0.00_wp, & ! 12 |
---|
| 517 | 0.25_wp, 0.34_wp, 0.27_wp, 0.11_wp, & ! 13 |
---|
| 518 | 0.00_wp, 0.00_wp, 0.00_wp, 0.00_wp, & ! 14 |
---|
| 519 | 0.00_wp, 0.00_wp, 0.00_wp, 0.00_wp, & ! 15 |
---|
| 520 | 0.23_wp, 0.36_wp, 0.30_wp, 0.11_wp, & ! 16 |
---|
| 521 | 0.23_wp, 0.36_wp, 0.30_wp, 0.11_wp, & ! 17 |
---|
| 522 | 0.19_wp, 0.35_wp, 0.36_wp, 0.10_wp, & ! 18 |
---|
[1788] | 523 | 0.19_wp, 0.35_wp, 0.36_wp, 0.10_wp, & ! 19 |
---|
| 524 | 0.00_wp, 0.00_wp, 0.00_wp, 0.00_wp & ! 20 |
---|
| 525 | /), (/ 4, 20 /) ) |
---|
[1496] | 526 | |
---|
| 527 | ! |
---|
| 528 | !-- Soil parameters according to the following porosity classes (soil_type) |
---|
[1551] | 529 | |
---|
[1496] | 530 | ! |
---|
[1551] | 531 | !-- Soil parameters I alpha_vg, l_vg, n_vg, gamma_w_sat |
---|
[1585] | 532 | REAL(wp), DIMENSION(0:3,1:7), PARAMETER :: soil_pars = RESHAPE( (/ & |
---|
[1496] | 533 | 3.83_wp, 1.250_wp, 1.38_wp, 6.94E-6_wp, & ! 1 |
---|
| 534 | 3.14_wp, -2.342_wp, 1.28_wp, 1.16E-6_wp, & ! 2 |
---|
| 535 | 0.83_wp, -0.588_wp, 1.25_wp, 0.26E-6_wp, & ! 3 |
---|
| 536 | 3.67_wp, -1.977_wp, 1.10_wp, 2.87E-6_wp, & ! 4 |
---|
| 537 | 2.65_wp, 2.500_wp, 1.10_wp, 1.74E-6_wp, & ! 5 |
---|
[1551] | 538 | 1.30_wp, 0.400_wp, 1.20_wp, 0.93E-6_wp, & ! 6 |
---|
| 539 | 0.00_wp, 0.00_wp, 0.00_wp, 0.57E-6_wp & ! 7 |
---|
| 540 | /), (/ 4, 7 /) ) |
---|
[1496] | 541 | |
---|
| 542 | ! |
---|
| 543 | !-- Soil parameters II m_sat, m_fc, m_wilt, m_res |
---|
[1585] | 544 | REAL(wp), DIMENSION(0:3,1:7), PARAMETER :: m_soil_pars = RESHAPE( (/ & |
---|
[1496] | 545 | 0.403_wp, 0.244_wp, 0.059_wp, 0.025_wp, & ! 1 |
---|
| 546 | 0.439_wp, 0.347_wp, 0.151_wp, 0.010_wp, & ! 2 |
---|
| 547 | 0.430_wp, 0.383_wp, 0.133_wp, 0.010_wp, & ! 3 |
---|
| 548 | 0.520_wp, 0.448_wp, 0.279_wp, 0.010_wp, & ! 4 |
---|
| 549 | 0.614_wp, 0.541_wp, 0.335_wp, 0.010_wp, & ! 5 |
---|
[1551] | 550 | 0.766_wp, 0.663_wp, 0.267_wp, 0.010_wp, & ! 6 |
---|
| 551 | 0.472_wp, 0.323_wp, 0.171_wp, 0.000_wp & ! 7 |
---|
| 552 | /), (/ 4, 7 /) ) |
---|
[1496] | 553 | |
---|
| 554 | |
---|
| 555 | SAVE |
---|
| 556 | |
---|
| 557 | |
---|
| 558 | PRIVATE |
---|
| 559 | |
---|
[1817] | 560 | |
---|
[1551] | 561 | ! |
---|
[1817] | 562 | !-- Public functions |
---|
| 563 | PUBLIC lsm_check_data_output, lsm_check_data_output_pr, & |
---|
[1972] | 564 | lsm_check_parameters, lsm_define_netcdf_grid, lsm_3d_data_averaging,& |
---|
| 565 | lsm_data_output_2d, lsm_data_output_3d, lsm_energy_balance, & |
---|
| 566 | lsm_header, lsm_init, lsm_init_arrays, lsm_parin, lsm_soil_model, & |
---|
| 567 | lsm_swap_timelevel, lsm_read_restart_data, lsm_last_actions |
---|
[1817] | 568 | ! |
---|
[1551] | 569 | !-- Public parameters, constants and initial values |
---|
[1826] | 570 | PUBLIC land_surface, skip_time_do_lsm |
---|
[1496] | 571 | |
---|
[1551] | 572 | ! |
---|
[1783] | 573 | !-- Public grid variables |
---|
| 574 | PUBLIC nzb_soil, nzs, nzt_soil, zs |
---|
[1496] | 575 | |
---|
[1551] | 576 | ! |
---|
| 577 | !-- Public 2D output variables |
---|
[1972] | 578 | PUBLIC ghf_eb, qsws_eb, qsws_liq_eb, qsws_soil_eb,qsws_veg_eb, r_a, r_s, & |
---|
| 579 | shf_eb |
---|
[1551] | 580 | |
---|
| 581 | ! |
---|
| 582 | !-- Public prognostic variables |
---|
[1972] | 583 | PUBLIC m_soil, t_soil |
---|
[1551] | 584 | |
---|
[1496] | 585 | |
---|
[1817] | 586 | INTERFACE lsm_check_data_output |
---|
| 587 | MODULE PROCEDURE lsm_check_data_output |
---|
| 588 | END INTERFACE lsm_check_data_output |
---|
| 589 | |
---|
| 590 | INTERFACE lsm_check_data_output_pr |
---|
| 591 | MODULE PROCEDURE lsm_check_data_output_pr |
---|
| 592 | END INTERFACE lsm_check_data_output_pr |
---|
| 593 | |
---|
| 594 | INTERFACE lsm_check_parameters |
---|
| 595 | MODULE PROCEDURE lsm_check_parameters |
---|
| 596 | END INTERFACE lsm_check_parameters |
---|
| 597 | |
---|
[1972] | 598 | INTERFACE lsm_3d_data_averaging |
---|
| 599 | MODULE PROCEDURE lsm_3d_data_averaging |
---|
| 600 | END INTERFACE lsm_3d_data_averaging |
---|
| 601 | |
---|
| 602 | INTERFACE lsm_data_output_2d |
---|
| 603 | MODULE PROCEDURE lsm_data_output_2d |
---|
| 604 | END INTERFACE lsm_data_output_2d |
---|
| 605 | |
---|
| 606 | INTERFACE lsm_data_output_3d |
---|
| 607 | MODULE PROCEDURE lsm_data_output_3d |
---|
| 608 | END INTERFACE lsm_data_output_3d |
---|
| 609 | |
---|
| 610 | INTERFACE lsm_define_netcdf_grid |
---|
| 611 | MODULE PROCEDURE lsm_define_netcdf_grid |
---|
| 612 | END INTERFACE lsm_define_netcdf_grid |
---|
| 613 | |
---|
[1496] | 614 | INTERFACE lsm_energy_balance |
---|
| 615 | MODULE PROCEDURE lsm_energy_balance |
---|
| 616 | END INTERFACE lsm_energy_balance |
---|
| 617 | |
---|
[1817] | 618 | INTERFACE lsm_header |
---|
| 619 | MODULE PROCEDURE lsm_header |
---|
| 620 | END INTERFACE lsm_header |
---|
| 621 | |
---|
| 622 | INTERFACE lsm_init |
---|
| 623 | MODULE PROCEDURE lsm_init |
---|
| 624 | END INTERFACE lsm_init |
---|
| 625 | |
---|
| 626 | INTERFACE lsm_init_arrays |
---|
| 627 | MODULE PROCEDURE lsm_init_arrays |
---|
| 628 | END INTERFACE lsm_init_arrays |
---|
| 629 | |
---|
| 630 | INTERFACE lsm_parin |
---|
| 631 | MODULE PROCEDURE lsm_parin |
---|
| 632 | END INTERFACE lsm_parin |
---|
| 633 | |
---|
[1496] | 634 | INTERFACE lsm_soil_model |
---|
| 635 | MODULE PROCEDURE lsm_soil_model |
---|
| 636 | END INTERFACE lsm_soil_model |
---|
| 637 | |
---|
[1551] | 638 | INTERFACE lsm_swap_timelevel |
---|
| 639 | MODULE PROCEDURE lsm_swap_timelevel |
---|
| 640 | END INTERFACE lsm_swap_timelevel |
---|
[1496] | 641 | |
---|
[1972] | 642 | INTERFACE lsm_read_restart_data |
---|
| 643 | MODULE PROCEDURE lsm_read_restart_data |
---|
| 644 | END INTERFACE lsm_read_restart_data |
---|
| 645 | |
---|
| 646 | INTERFACE lsm_last_actions |
---|
| 647 | MODULE PROCEDURE lsm_last_actions |
---|
| 648 | END INTERFACE lsm_last_actions |
---|
| 649 | |
---|
[1496] | 650 | CONTAINS |
---|
| 651 | |
---|
| 652 | !------------------------------------------------------------------------------! |
---|
| 653 | ! Description: |
---|
| 654 | ! ------------ |
---|
[1817] | 655 | !> Check data output for land surface model |
---|
[1496] | 656 | !------------------------------------------------------------------------------! |
---|
[1976] | 657 | SUBROUTINE lsm_check_data_output( var, unit, i, ilen, k ) |
---|
[1817] | 658 | |
---|
| 659 | |
---|
[1976] | 660 | USE control_parameters, & |
---|
| 661 | ONLY: data_output, message_string |
---|
[1496] | 662 | |
---|
[1976] | 663 | IMPLICIT NONE |
---|
[1496] | 664 | |
---|
[1976] | 665 | CHARACTER (LEN=*) :: unit !< |
---|
| 666 | CHARACTER (LEN=*) :: var !< |
---|
[1496] | 667 | |
---|
[1976] | 668 | INTEGER(iwp) :: i |
---|
| 669 | INTEGER(iwp) :: ilen |
---|
| 670 | INTEGER(iwp) :: k |
---|
[1496] | 671 | |
---|
[1976] | 672 | SELECT CASE ( TRIM( var ) ) |
---|
[1496] | 673 | |
---|
[1976] | 674 | CASE ( 'm_soil' ) |
---|
| 675 | IF ( .NOT. land_surface ) THEN |
---|
| 676 | message_string = 'output of "' // TRIM( var ) // '" requi' // & |
---|
| 677 | 'res land_surface = .TRUE.' |
---|
| 678 | CALL message( 'check_parameters', 'PA0404', 1, 2, 0, 6, 0 ) |
---|
| 679 | ENDIF |
---|
| 680 | unit = 'm3/m3' |
---|
| 681 | |
---|
| 682 | CASE ( 't_soil' ) |
---|
| 683 | IF ( .NOT. land_surface ) THEN |
---|
| 684 | message_string = 'output of "' // TRIM( var ) // '" requi' // & |
---|
| 685 | 'res land_surface = .TRUE.' |
---|
| 686 | CALL message( 'check_parameters', 'PA0404', 1, 2, 0, 6, 0 ) |
---|
| 687 | ENDIF |
---|
| 688 | unit = 'K' |
---|
[1817] | 689 | |
---|
[1976] | 690 | CASE ( 'lai*', 'c_liq*', 'c_soil*', 'c_veg*', 'ghf_eb*', 'm_liq_eb*',& |
---|
| 691 | 'qsws_eb*', 'qsws_liq_eb*', 'qsws_soil_eb*', 'qsws_veg_eb*', & |
---|
| 692 | 'r_a*', 'r_s*', 'shf_eb*' ) |
---|
| 693 | IF ( k == 0 .OR. data_output(i)(ilen-2:ilen) /= '_xy' ) THEN |
---|
| 694 | message_string = 'illegal value for data_output: "' // & |
---|
| 695 | TRIM( var ) // '" & only 2d-horizontal ' // & |
---|
| 696 | 'cross sections are allowed for this value' |
---|
| 697 | CALL message( 'check_parameters', 'PA0111', 1, 2, 0, 6, 0 ) |
---|
| 698 | ENDIF |
---|
| 699 | IF ( TRIM( var ) == 'lai*' .AND. .NOT. land_surface ) THEN |
---|
| 700 | message_string = 'output of "' // TRIM( var ) // '" requi' // & |
---|
| 701 | 'res land_surface = .TRUE.' |
---|
| 702 | CALL message( 'check_parameters', 'PA0404', 1, 2, 0, 6, 0 ) |
---|
| 703 | ENDIF |
---|
| 704 | IF ( TRIM( var ) == 'c_liq*' .AND. .NOT. land_surface ) THEN |
---|
| 705 | message_string = 'output of "' // TRIM( var ) // '" requi' // & |
---|
| 706 | 'res land_surface = .TRUE.' |
---|
| 707 | CALL message( 'check_parameters', 'PA0404', 1, 2, 0, 6, 0 ) |
---|
| 708 | ENDIF |
---|
| 709 | IF ( TRIM( var ) == 'c_soil*' .AND. .NOT. land_surface ) THEN |
---|
| 710 | message_string = 'output of "' // TRIM( var ) // '" requi' // & |
---|
| 711 | 'res land_surface = .TRUE.' |
---|
| 712 | CALL message( 'check_parameters', 'PA0404', 1, 2, 0, 6, 0 ) |
---|
| 713 | ENDIF |
---|
| 714 | IF ( TRIM( var ) == 'c_veg*' .AND. .NOT. land_surface ) THEN |
---|
| 715 | message_string = 'output of "' // TRIM( var ) // '" requi' // & |
---|
| 716 | 'res land_surface = .TRUE.' |
---|
| 717 | CALL message( 'check_parameters', 'PA0401', 1, 2, 0, 6, 0 ) |
---|
| 718 | ENDIF |
---|
| 719 | IF ( TRIM( var ) == 'ghf_eb*' .AND. .NOT. land_surface ) THEN |
---|
| 720 | message_string = 'output of "' // TRIM( var ) // '" requi' // & |
---|
| 721 | 'res land_surface = .TRUE.' |
---|
| 722 | CALL message( 'check_parameters', 'PA0404', 1, 2, 0, 6, 0 ) |
---|
| 723 | ENDIF |
---|
| 724 | IF ( TRIM( var ) == 'm_liq_eb*' .AND. .NOT. land_surface ) THEN |
---|
| 725 | message_string = 'output of "' // TRIM( var ) // '" requi' // & |
---|
| 726 | 'res land_surface = .TRUE.' |
---|
| 727 | CALL message( 'check_parameters', 'PA0404', 1, 2, 0, 6, 0 ) |
---|
| 728 | ENDIF |
---|
| 729 | IF ( TRIM( var ) == 'qsws_eb*' .AND. .NOT. land_surface ) THEN |
---|
| 730 | message_string = 'output of "' // TRIM( var ) // '" requi' // & |
---|
| 731 | 'res land_surface = .TRUE.' |
---|
| 732 | CALL message( 'check_parameters', 'PA0404', 1, 2, 0, 6, 0 ) |
---|
| 733 | ENDIF |
---|
| 734 | IF ( TRIM( var ) == 'qsws_liq_eb*' .AND. .NOT. land_surface ) & |
---|
| 735 | THEN |
---|
| 736 | message_string = 'output of "' // TRIM( var ) // '" requi' // & |
---|
| 737 | 'res land_surface = .TRUE.' |
---|
| 738 | CALL message( 'check_parameters', 'PA0404', 1, 2, 0, 6, 0 ) |
---|
| 739 | ENDIF |
---|
| 740 | IF ( TRIM( var ) == 'qsws_soil_eb*' .AND. .NOT. land_surface ) & |
---|
| 741 | THEN |
---|
| 742 | message_string = 'output of "' // TRIM( var ) // '" requi' // & |
---|
| 743 | 'res land_surface = .TRUE.' |
---|
| 744 | CALL message( 'check_parameters', 'PA0404', 1, 2, 0, 6, 0 ) |
---|
| 745 | ENDIF |
---|
| 746 | IF ( TRIM( var ) == 'qsws_veg_eb*' .AND. .NOT. land_surface ) & |
---|
| 747 | THEN |
---|
| 748 | message_string = 'output of "' // TRIM( var ) // '" requi' // & |
---|
| 749 | 'res land_surface = .TRUE.' |
---|
| 750 | CALL message( 'check_parameters', 'PA0404', 1, 2, 0, 6, 0 ) |
---|
| 751 | ENDIF |
---|
| 752 | IF ( TRIM( var ) == 'r_a*' .AND. .NOT. land_surface ) & |
---|
| 753 | THEN |
---|
| 754 | message_string = 'output of "' // TRIM( var ) // '" requi' // & |
---|
| 755 | 'res land_surface = .TRUE.' |
---|
| 756 | CALL message( 'check_parameters', 'PA0404', 1, 2, 0, 6, 0 ) |
---|
| 757 | ENDIF |
---|
| 758 | IF ( TRIM( var ) == 'r_s*' .AND. .NOT. land_surface ) & |
---|
| 759 | THEN |
---|
| 760 | message_string = 'output of "' // TRIM( var ) // '" requi' // & |
---|
| 761 | 'res land_surface = .TRUE.' |
---|
| 762 | CALL message( 'check_parameters', 'PA0404', 1, 2, 0, 6, 0 ) |
---|
| 763 | ENDIF |
---|
[1496] | 764 | |
---|
[1976] | 765 | IF ( TRIM( var ) == 'lai*' ) unit = 'none' |
---|
| 766 | IF ( TRIM( var ) == 'c_liq*' ) unit = 'none' |
---|
| 767 | IF ( TRIM( var ) == 'c_soil*') unit = 'none' |
---|
| 768 | IF ( TRIM( var ) == 'c_veg*' ) unit = 'none' |
---|
| 769 | IF ( TRIM( var ) == 'ghf_eb*') unit = 'W/m2' |
---|
| 770 | IF ( TRIM( var ) == 'm_liq_eb*' ) unit = 'm' |
---|
| 771 | IF ( TRIM( var ) == 'qsws_eb*' ) unit = 'W/m2' |
---|
| 772 | IF ( TRIM( var ) == 'qsws_liq_eb*' ) unit = 'W/m2' |
---|
| 773 | IF ( TRIM( var ) == 'qsws_soil_eb*' ) unit = 'W/m2' |
---|
| 774 | IF ( TRIM( var ) == 'qsws_veg_eb*' ) unit = 'W/m2' |
---|
| 775 | IF ( TRIM( var ) == 'r_a*') unit = 's/m' |
---|
| 776 | IF ( TRIM( var ) == 'r_s*') unit = 's/m' |
---|
| 777 | IF ( TRIM( var ) == 'shf_eb*') unit = 'W/m2' |
---|
[1817] | 778 | |
---|
[1976] | 779 | CASE DEFAULT |
---|
| 780 | unit = 'illegal' |
---|
[1500] | 781 | |
---|
[1976] | 782 | END SELECT |
---|
[1551] | 783 | |
---|
| 784 | |
---|
[1976] | 785 | END SUBROUTINE lsm_check_data_output |
---|
[1551] | 786 | |
---|
[1976] | 787 | |
---|
| 788 | !------------------------------------------------------------------------------! |
---|
| 789 | ! Description: |
---|
| 790 | ! ------------ |
---|
| 791 | !> Check data output of profiles for land surface model |
---|
| 792 | !------------------------------------------------------------------------------! |
---|
[1817] | 793 | SUBROUTINE lsm_check_data_output_pr( variable, var_count, unit, dopr_unit ) |
---|
| 794 | |
---|
[1976] | 795 | USE control_parameters, & |
---|
| 796 | ONLY: data_output_pr, message_string |
---|
[1551] | 797 | |
---|
[1976] | 798 | USE indices |
---|
[1691] | 799 | |
---|
[1976] | 800 | USE profil_parameter |
---|
[1551] | 801 | |
---|
[1976] | 802 | USE statistics |
---|
[1551] | 803 | |
---|
[1976] | 804 | IMPLICIT NONE |
---|
[1817] | 805 | |
---|
[1976] | 806 | CHARACTER (LEN=*) :: unit !< |
---|
| 807 | CHARACTER (LEN=*) :: variable !< |
---|
| 808 | CHARACTER (LEN=*) :: dopr_unit !< local value of dopr_unit |
---|
[1817] | 809 | |
---|
[1976] | 810 | INTEGER(iwp) :: user_pr_index !< |
---|
| 811 | INTEGER(iwp) :: var_count !< |
---|
[1551] | 812 | |
---|
[1976] | 813 | SELECT CASE ( TRIM( variable ) ) |
---|
[1817] | 814 | |
---|
[1976] | 815 | CASE ( 't_soil', '#t_soil' ) |
---|
| 816 | IF ( .NOT. land_surface ) THEN |
---|
| 817 | message_string = 'data_output_pr = ' // & |
---|
| 818 | TRIM( data_output_pr(var_count) ) // ' is' // & |
---|
| 819 | 'not implemented for land_surface = .FALSE.' |
---|
| 820 | CALL message( 'check_parameters', 'PA0402', 1, 2, 0, 6, 0 ) |
---|
| 821 | ELSE |
---|
| 822 | dopr_index(var_count) = 89 |
---|
| 823 | dopr_unit = 'K' |
---|
| 824 | hom(0:nzs-1,2,89,:) = SPREAD( - zs, 2, statistic_regions+1 ) |
---|
| 825 | IF ( data_output_pr(var_count)(1:1) == '#' ) THEN |
---|
| 826 | dopr_initial_index(var_count) = 90 |
---|
| 827 | hom(0:nzs-1,2,90,:) = SPREAD( - zs, 2, statistic_regions+1 ) |
---|
| 828 | data_output_pr(var_count) = data_output_pr(var_count)(2:) |
---|
[1817] | 829 | ENDIF |
---|
[1976] | 830 | unit = dopr_unit |
---|
| 831 | ENDIF |
---|
[1551] | 832 | |
---|
[1976] | 833 | CASE ( 'm_soil', '#m_soil' ) |
---|
| 834 | IF ( .NOT. land_surface ) THEN |
---|
| 835 | message_string = 'data_output_pr = ' // & |
---|
| 836 | TRIM( data_output_pr(var_count) ) // ' is' // & |
---|
| 837 | ' not implemented for land_surface = .FALSE.' |
---|
| 838 | CALL message( 'check_parameters', 'PA0402', 1, 2, 0, 6, 0 ) |
---|
| 839 | ELSE |
---|
| 840 | dopr_index(var_count) = 91 |
---|
| 841 | dopr_unit = 'm3/m3' |
---|
| 842 | hom(0:nzs-1,2,91,:) = SPREAD( - zs, 2, statistic_regions+1 ) |
---|
| 843 | IF ( data_output_pr(var_count)(1:1) == '#' ) THEN |
---|
| 844 | dopr_initial_index(var_count) = 92 |
---|
| 845 | hom(0:nzs-1,2,92,:) = SPREAD( - zs, 2, statistic_regions+1 ) |
---|
| 846 | data_output_pr(var_count) = data_output_pr(var_count)(2:) |
---|
[1817] | 847 | ENDIF |
---|
[1976] | 848 | unit = dopr_unit |
---|
| 849 | ENDIF |
---|
[1551] | 850 | |
---|
| 851 | |
---|
[1976] | 852 | CASE DEFAULT |
---|
| 853 | unit = 'illegal' |
---|
[1500] | 854 | |
---|
[1976] | 855 | END SELECT |
---|
[1496] | 856 | |
---|
| 857 | |
---|
[1976] | 858 | END SUBROUTINE lsm_check_data_output_pr |
---|
[1817] | 859 | |
---|
| 860 | |
---|
| 861 | !------------------------------------------------------------------------------! |
---|
| 862 | ! Description: |
---|
| 863 | ! ------------ |
---|
| 864 | !> Check parameters routine for land surface model |
---|
| 865 | !------------------------------------------------------------------------------! |
---|
[1976] | 866 | SUBROUTINE lsm_check_parameters |
---|
[1496] | 867 | |
---|
[1976] | 868 | USE control_parameters, & |
---|
| 869 | ONLY: bc_pt_b, bc_q_b, constant_flux_layer, message_string, & |
---|
| 870 | most_method, topography |
---|
[1817] | 871 | |
---|
[1976] | 872 | USE radiation_model_mod, & |
---|
| 873 | ONLY: radiation |
---|
[1817] | 874 | |
---|
| 875 | |
---|
[1976] | 876 | IMPLICIT NONE |
---|
[1817] | 877 | |
---|
| 878 | |
---|
[1496] | 879 | ! |
---|
[1976] | 880 | !-- Dirichlet boundary conditions are required as the surface fluxes are |
---|
| 881 | !-- calculated from the temperature/humidity gradients in the land surface |
---|
| 882 | !-- model |
---|
| 883 | IF ( bc_pt_b == 'neumann' .OR. bc_q_b == 'neumann' ) THEN |
---|
| 884 | message_string = 'lsm requires setting of'// & |
---|
| 885 | 'bc_pt_b = "dirichlet" and '// & |
---|
| 886 | 'bc_q_b = "dirichlet"' |
---|
| 887 | CALL message( 'check_parameters', 'PA0399', 1, 2, 0, 6, 0 ) |
---|
| 888 | ENDIF |
---|
[1496] | 889 | |
---|
[1976] | 890 | IF ( .NOT. constant_flux_layer ) THEN |
---|
| 891 | message_string = 'lsm requires '// & |
---|
| 892 | 'constant_flux_layer = .T.' |
---|
| 893 | CALL message( 'check_parameters', 'PA0400', 1, 2, 0, 6, 0 ) |
---|
| 894 | ENDIF |
---|
[1496] | 895 | |
---|
[1976] | 896 | IF ( topography /= 'flat' ) THEN |
---|
| 897 | message_string = 'lsm cannot be used ' // & |
---|
| 898 | 'in combination with topography /= "flat"' |
---|
| 899 | CALL message( 'check_parameters', 'PA0415', 1, 2, 0, 6, 0 ) |
---|
| 900 | ENDIF |
---|
[1496] | 901 | |
---|
[1976] | 902 | IF ( ( veg_type == 14 .OR. veg_type == 15 ) .AND. & |
---|
| 903 | most_method == 'lookup' ) THEN |
---|
| 904 | WRITE( message_string, * ) 'veg_type = ', veg_type, ' is not ', & |
---|
| 905 | 'allowed in combination with ', & |
---|
| 906 | 'most_method = ', most_method |
---|
| 907 | CALL message( 'check_parameters', 'PA0417', 1, 2, 0, 6, 0 ) |
---|
| 908 | ENDIF |
---|
| 909 | |
---|
| 910 | IF ( veg_type == 0 ) THEN |
---|
| 911 | IF ( SUM( root_fraction ) /= 1.0_wp ) THEN |
---|
| 912 | message_string = 'veg_type = 0 (user_defined)'// & |
---|
| 913 | 'requires setting of root_fraction(0:3)'// & |
---|
| 914 | '/= 9999999.9 and SUM(root_fraction) = 1' |
---|
| 915 | CALL message( 'check_parameters', 'PA0401', 1, 2, 0, 6, 0 ) |
---|
[1817] | 916 | ENDIF |
---|
[1551] | 917 | |
---|
[1976] | 918 | IF ( min_canopy_resistance == 9999999.9_wp ) THEN |
---|
| 919 | message_string = 'veg_type = 0 (user defined)'// & |
---|
| 920 | 'requires setting of min_canopy_resistance'// & |
---|
| 921 | '/= 9999999.9' |
---|
| 922 | CALL message( 'check_parameters', 'PA0401', 1, 2, 0, 6, 0 ) |
---|
| 923 | ENDIF |
---|
[1551] | 924 | |
---|
[1976] | 925 | IF ( leaf_area_index == 9999999.9_wp ) THEN |
---|
| 926 | message_string = 'veg_type = 0 (user_defined)'// & |
---|
| 927 | 'requires setting of leaf_area_index'// & |
---|
| 928 | '/= 9999999.9' |
---|
| 929 | CALL message( 'check_parameters', 'PA0401', 1, 2, 0, 6, 0 ) |
---|
| 930 | ENDIF |
---|
[1551] | 931 | |
---|
[1976] | 932 | IF ( vegetation_coverage == 9999999.9_wp ) THEN |
---|
| 933 | message_string = 'veg_type = 0 (user_defined)'// & |
---|
| 934 | 'requires setting of vegetation_coverage'// & |
---|
| 935 | '/= 9999999.9' |
---|
[1817] | 936 | CALL message( 'check_parameters', 'PA0401', 1, 2, 0, 6, 0 ) |
---|
[1976] | 937 | ENDIF |
---|
[1551] | 938 | |
---|
[1976] | 939 | IF ( canopy_resistance_coefficient == 9999999.9_wp) THEN |
---|
| 940 | message_string = 'veg_type = 0 (user_defined)'// & |
---|
| 941 | 'requires setting of'// & |
---|
| 942 | 'canopy_resistance_coefficient /= 9999999.9' |
---|
| 943 | CALL message( 'check_parameters', 'PA0401', 1, 2, 0, 6, 0 ) |
---|
| 944 | ENDIF |
---|
[1551] | 945 | |
---|
[1976] | 946 | IF ( lambda_surface_stable == 9999999.9_wp ) THEN |
---|
| 947 | message_string = 'veg_type = 0 (user_defined)'// & |
---|
| 948 | 'requires setting of lambda_surface_stable'// & |
---|
| 949 | '/= 9999999.9' |
---|
| 950 | CALL message( 'check_parameters', 'PA0401', 1, 2, 0, 6, 0 ) |
---|
| 951 | ENDIF |
---|
[1551] | 952 | |
---|
[1976] | 953 | IF ( lambda_surface_unstable == 9999999.9_wp ) THEN |
---|
| 954 | message_string = 'veg_type = 0 (user_defined)'// & |
---|
| 955 | 'requires setting of lambda_surface_unstable'// & |
---|
| 956 | '/= 9999999.9' |
---|
| 957 | CALL message( 'check_parameters', 'PA0401', 1, 2, 0, 6, 0 ) |
---|
| 958 | ENDIF |
---|
[1551] | 959 | |
---|
[1976] | 960 | IF ( f_shortwave_incoming == 9999999.9_wp ) THEN |
---|
| 961 | message_string = 'veg_type = 0 (user_defined)'// & |
---|
| 962 | 'requires setting of f_shortwave_incoming'// & |
---|
| 963 | '/= 9999999.9' |
---|
| 964 | CALL message( 'check_parameters', 'PA0401', 1, 2, 0, 6, 0 ) |
---|
| 965 | ENDIF |
---|
[1551] | 966 | |
---|
[1976] | 967 | IF ( z0_eb == 9999999.9_wp ) THEN |
---|
| 968 | message_string = 'veg_type = 0 (user_defined)'// & |
---|
| 969 | 'requires setting of z0_eb'// & |
---|
| 970 | '/= 9999999.9' |
---|
| 971 | CALL message( 'check_parameters', 'PA0401', 1, 2, 0, 6, 0 ) |
---|
| 972 | ENDIF |
---|
[1551] | 973 | |
---|
[1976] | 974 | IF ( z0h_eb == 9999999.9_wp ) THEN |
---|
| 975 | message_string = 'veg_type = 0 (user_defined)'// & |
---|
| 976 | 'requires setting of z0h_eb'// & |
---|
| 977 | '/= 9999999.9' |
---|
| 978 | CALL message( 'check_parameters', 'PA0401', 1, 2, 0, 6, 0 ) |
---|
[1551] | 979 | ENDIF |
---|
| 980 | |
---|
[1496] | 981 | |
---|
[1976] | 982 | ENDIF |
---|
[1817] | 983 | |
---|
[1976] | 984 | IF ( soil_type == 0 ) THEN |
---|
[1817] | 985 | |
---|
[1976] | 986 | IF ( alpha_vangenuchten == 9999999.9_wp ) THEN |
---|
| 987 | message_string = 'soil_type = 0 (user_defined)'// & |
---|
| 988 | 'requires setting of alpha_vangenuchten'// & |
---|
| 989 | '/= 9999999.9' |
---|
| 990 | CALL message( 'check_parameters', 'PA0403', 1, 2, 0, 6, 0 ) |
---|
| 991 | ENDIF |
---|
[1817] | 992 | |
---|
[1976] | 993 | IF ( l_vangenuchten == 9999999.9_wp ) THEN |
---|
| 994 | message_string = 'soil_type = 0 (user_defined)'// & |
---|
| 995 | 'requires setting of l_vangenuchten'// & |
---|
| 996 | '/= 9999999.9' |
---|
| 997 | CALL message( 'check_parameters', 'PA0403', 1, 2, 0, 6, 0 ) |
---|
| 998 | ENDIF |
---|
[1817] | 999 | |
---|
[1976] | 1000 | IF ( n_vangenuchten == 9999999.9_wp ) THEN |
---|
| 1001 | message_string = 'soil_type = 0 (user_defined)'// & |
---|
| 1002 | 'requires setting of n_vangenuchten'// & |
---|
| 1003 | '/= 9999999.9' |
---|
| 1004 | CALL message( 'check_parameters', 'PA0403', 1, 2, 0, 6, 0 ) |
---|
| 1005 | ENDIF |
---|
[1817] | 1006 | |
---|
[1976] | 1007 | IF ( hydraulic_conductivity == 9999999.9_wp ) THEN |
---|
| 1008 | message_string = 'soil_type = 0 (user_defined)'// & |
---|
| 1009 | 'requires setting of hydraulic_conductivity'// & |
---|
| 1010 | '/= 9999999.9' |
---|
| 1011 | CALL message( 'check_parameters', 'PA0403', 1, 2, 0, 6, 0 ) |
---|
| 1012 | ENDIF |
---|
[1496] | 1013 | |
---|
[1976] | 1014 | IF ( saturation_moisture == 9999999.9_wp ) THEN |
---|
| 1015 | message_string = 'soil_type = 0 (user_defined)'// & |
---|
| 1016 | 'requires setting of saturation_moisture'// & |
---|
| 1017 | '/= 9999999.9' |
---|
| 1018 | CALL message( 'check_parameters', 'PA0403', 1, 2, 0, 6, 0 ) |
---|
| 1019 | ENDIF |
---|
[1496] | 1020 | |
---|
[1976] | 1021 | IF ( field_capacity == 9999999.9_wp ) THEN |
---|
| 1022 | message_string = 'soil_type = 0 (user_defined)'// & |
---|
| 1023 | 'requires setting of field_capacity'// & |
---|
| 1024 | '/= 9999999.9' |
---|
| 1025 | CALL message( 'check_parameters', 'PA0403', 1, 2, 0, 6, 0 ) |
---|
| 1026 | ENDIF |
---|
[1553] | 1027 | |
---|
[1976] | 1028 | IF ( wilting_point == 9999999.9_wp ) THEN |
---|
| 1029 | message_string = 'soil_type = 0 (user_defined)'// & |
---|
| 1030 | 'requires setting of wilting_point'// & |
---|
| 1031 | '/= 9999999.9' |
---|
| 1032 | CALL message( 'check_parameters', 'PA0403', 1, 2, 0, 6, 0 ) |
---|
[1496] | 1033 | ENDIF |
---|
| 1034 | |
---|
[1976] | 1035 | IF ( residual_moisture == 9999999.9_wp ) THEN |
---|
| 1036 | message_string = 'soil_type = 0 (user_defined)'// & |
---|
| 1037 | 'requires setting of residual_moisture'// & |
---|
| 1038 | '/= 9999999.9' |
---|
| 1039 | CALL message( 'check_parameters', 'PA0403', 1, 2, 0, 6, 0 ) |
---|
[1788] | 1040 | ENDIF |
---|
[1976] | 1041 | |
---|
| 1042 | ENDIF |
---|
| 1043 | |
---|
| 1044 | IF ( .NOT. radiation ) THEN |
---|
| 1045 | message_string = 'lsm requires '// & |
---|
| 1046 | 'radiation = .T.' |
---|
| 1047 | CALL message( 'check_parameters', 'PA0400', 1, 2, 0, 6, 0 ) |
---|
| 1048 | ENDIF |
---|
[1817] | 1049 | |
---|
| 1050 | |
---|
[1976] | 1051 | END SUBROUTINE lsm_check_parameters |
---|
[1817] | 1052 | |
---|
[1496] | 1053 | !------------------------------------------------------------------------------! |
---|
| 1054 | ! Description: |
---|
| 1055 | ! ------------ |
---|
[1682] | 1056 | !> Solver for the energy balance at the surface. |
---|
[1496] | 1057 | !------------------------------------------------------------------------------! |
---|
[1976] | 1058 | SUBROUTINE lsm_energy_balance |
---|
[1496] | 1059 | |
---|
| 1060 | |
---|
[1976] | 1061 | IMPLICIT NONE |
---|
[1496] | 1062 | |
---|
[1976] | 1063 | INTEGER(iwp) :: i !< running index |
---|
| 1064 | INTEGER(iwp) :: j !< running index |
---|
| 1065 | INTEGER(iwp) :: k, ks !< running index |
---|
[1496] | 1066 | |
---|
[1976] | 1067 | REAL(wp) :: c_surface_tmp,& !< temporary variable for storing the volumetric heat capacity of the surface |
---|
| 1068 | f1, & !< resistance correction term 1 |
---|
| 1069 | f2, & !< resistance correction term 2 |
---|
| 1070 | f3, & !< resistance correction term 3 |
---|
| 1071 | m_min, & !< minimum soil moisture |
---|
| 1072 | e, & !< water vapour pressure |
---|
| 1073 | e_s, & !< water vapour saturation pressure |
---|
| 1074 | e_s_dt, & !< derivate of e_s with respect to T |
---|
| 1075 | tend, & !< tendency |
---|
| 1076 | dq_s_dt, & !< derivate of q_s with respect to T |
---|
| 1077 | coef_1, & !< coef. for prognostic equation |
---|
| 1078 | coef_2, & !< coef. for prognostic equation |
---|
| 1079 | f_qsws, & !< factor for qsws_eb |
---|
| 1080 | f_qsws_veg, & !< factor for qsws_veg_eb |
---|
| 1081 | f_qsws_soil, & !< factor for qsws_soil_eb |
---|
| 1082 | f_qsws_liq, & !< factor for qsws_liq_eb |
---|
| 1083 | f_shf, & !< factor for shf_eb |
---|
| 1084 | lambda_surface, & !< Current value of lambda_surface |
---|
| 1085 | m_liq_eb_max, & !< maxmimum value of the liq. water reservoir |
---|
| 1086 | pt1, & !< potential temperature at first grid level |
---|
| 1087 | qv1 !< specific humidity at first grid level |
---|
[1496] | 1088 | |
---|
| 1089 | ! |
---|
[1976] | 1090 | !-- Calculate the exner function for the current time step |
---|
| 1091 | exn = ( surface_pressure / 1000.0_wp )**0.286_wp |
---|
[1496] | 1092 | |
---|
[1976] | 1093 | DO i = nxlg, nxrg |
---|
| 1094 | DO j = nysg, nyng |
---|
| 1095 | k = nzb_s_inner(j,i) |
---|
[1496] | 1096 | |
---|
| 1097 | ! |
---|
[1976] | 1098 | !-- Set lambda_surface according to stratification between skin layer and soil |
---|
| 1099 | IF ( .NOT. pave_surface(j,i) ) THEN |
---|
[1788] | 1100 | |
---|
[1976] | 1101 | c_surface_tmp = c_surface |
---|
[1788] | 1102 | |
---|
[1976] | 1103 | IF ( t_surface(j,i) >= t_soil(nzb_soil,j,i)) THEN |
---|
| 1104 | lambda_surface = lambda_surface_s(j,i) |
---|
[1496] | 1105 | ELSE |
---|
[1976] | 1106 | lambda_surface = lambda_surface_u(j,i) |
---|
| 1107 | ENDIF |
---|
| 1108 | ELSE |
---|
[1788] | 1109 | |
---|
[1976] | 1110 | c_surface_tmp = pave_heat_capacity * dz_soil(nzb_soil) * 0.5_wp |
---|
| 1111 | lambda_surface = pave_heat_conductivity * ddz_soil(nzb_soil) |
---|
[1788] | 1112 | |
---|
[1976] | 1113 | ENDIF |
---|
[1500] | 1114 | |
---|
[1496] | 1115 | ! |
---|
[1976] | 1116 | !-- First step: calculate aerodyamic resistance. As pt, us, ts |
---|
| 1117 | !-- are not available for the prognostic time step, data from the last |
---|
| 1118 | !-- time step is used here. Note that this formulation is the |
---|
| 1119 | !-- equivalent to the ECMWF formulation using drag coefficients |
---|
| 1120 | IF ( cloud_physics ) THEN |
---|
| 1121 | pt1 = pt(k+1,j,i) + l_d_cp * pt_d_t(k+1) * ql(k+1,j,i) |
---|
| 1122 | qv1 = q(k+1,j,i) - ql(k+1,j,i) |
---|
| 1123 | ELSE |
---|
| 1124 | pt1 = pt(k+1,j,i) |
---|
| 1125 | qv1 = q(k+1,j,i) |
---|
| 1126 | ENDIF |
---|
[1496] | 1127 | |
---|
[1976] | 1128 | r_a(j,i) = (pt1 - pt(k,j,i)) / (ts(j,i) * us(j,i) + 1.0E-20_wp) |
---|
[1691] | 1129 | |
---|
[1496] | 1130 | ! |
---|
[1976] | 1131 | !-- Make sure that the resistance does not drop to zero for neutral |
---|
| 1132 | !-- stratification |
---|
| 1133 | IF ( ABS(r_a(j,i)) < 1.0_wp ) r_a(j,i) = 1.0_wp |
---|
[1709] | 1134 | |
---|
| 1135 | ! |
---|
[1976] | 1136 | !-- Second step: calculate canopy resistance r_canopy |
---|
| 1137 | !-- f1-f3 here are defined as 1/f1-f3 as in ECMWF documentation |
---|
[1496] | 1138 | |
---|
[1976] | 1139 | !-- f1: correction for incoming shortwave radiation (stomata close at |
---|
| 1140 | !-- night) |
---|
| 1141 | f1 = MIN( 1.0_wp, ( 0.004_wp * rad_sw_in(k,j,i) + 0.05_wp ) / & |
---|
| 1142 | (0.81_wp * (0.004_wp * rad_sw_in(k,j,i) & |
---|
| 1143 | + 1.0_wp)) ) |
---|
[1496] | 1144 | |
---|
[1709] | 1145 | |
---|
[1976] | 1146 | |
---|
[1496] | 1147 | ! |
---|
[1976] | 1148 | !-- f2: correction for soil moisture availability to plants (the |
---|
| 1149 | !-- integrated soil moisture must thus be considered here) |
---|
| 1150 | !-- f2 = 0 for very dry soils |
---|
| 1151 | m_total = 0.0_wp |
---|
| 1152 | DO ks = nzb_soil, nzt_soil |
---|
| 1153 | m_total = m_total + root_fr(ks,j,i) & |
---|
| 1154 | * MAX(m_soil(ks,j,i),m_wilt(j,i)) |
---|
| 1155 | ENDDO |
---|
[1496] | 1156 | |
---|
[1976] | 1157 | IF ( m_total > m_wilt(j,i) .AND. m_total < m_fc(j,i) ) THEN |
---|
| 1158 | f2 = ( m_total - m_wilt(j,i) ) / (m_fc(j,i) - m_wilt(j,i) ) |
---|
| 1159 | ELSEIF ( m_total >= m_fc(j,i) ) THEN |
---|
| 1160 | f2 = 1.0_wp |
---|
| 1161 | ELSE |
---|
| 1162 | f2 = 1.0E-20_wp |
---|
| 1163 | ENDIF |
---|
[1496] | 1164 | |
---|
| 1165 | ! |
---|
[1976] | 1166 | !-- Calculate water vapour pressure at saturation |
---|
| 1167 | e_s = 0.01_wp * 610.78_wp * EXP( 17.269_wp * ( t_surface(j,i) & |
---|
| 1168 | - 273.16_wp ) / ( t_surface(j,i) - 35.86_wp ) ) |
---|
[1496] | 1169 | |
---|
| 1170 | ! |
---|
[1976] | 1171 | !-- f3: correction for vapour pressure deficit |
---|
| 1172 | IF ( g_d(j,i) /= 0.0_wp ) THEN |
---|
[1496] | 1173 | ! |
---|
[1976] | 1174 | !-- Calculate vapour pressure |
---|
| 1175 | e = qv1 * surface_pressure / 0.622_wp |
---|
| 1176 | f3 = EXP ( -g_d(j,i) * (e_s - e) ) |
---|
| 1177 | ELSE |
---|
| 1178 | f3 = 1.0_wp |
---|
| 1179 | ENDIF |
---|
[1496] | 1180 | |
---|
| 1181 | ! |
---|
[1976] | 1182 | !-- Calculate canopy resistance. In case that c_veg is 0 (bare soils), |
---|
| 1183 | !-- this calculation is obsolete, as r_canopy is not used below. |
---|
| 1184 | !-- To do: check for very dry soil -> r_canopy goes to infinity |
---|
| 1185 | r_canopy(j,i) = r_canopy_min(j,i) / (lai(j,i) * f1 * f2 * f3 & |
---|
| 1186 | + 1.0E-20_wp) |
---|
[1496] | 1187 | |
---|
| 1188 | ! |
---|
[1976] | 1189 | !-- Third step: calculate bare soil resistance r_soil. The Clapp & |
---|
| 1190 | !-- Hornberger parametrization does not consider c_veg. |
---|
| 1191 | IF ( soil_type_2d(j,i) /= 7 ) THEN |
---|
| 1192 | m_min = c_veg(j,i) * m_wilt(j,i) + (1.0_wp - c_veg(j,i)) * & |
---|
| 1193 | m_res(j,i) |
---|
| 1194 | ELSE |
---|
| 1195 | m_min = m_wilt(j,i) |
---|
| 1196 | ENDIF |
---|
[1496] | 1197 | |
---|
[1976] | 1198 | f2 = ( m_soil(nzb_soil,j,i) - m_min ) / ( m_fc(j,i) - m_min ) |
---|
| 1199 | f2 = MAX(f2,1.0E-20_wp) |
---|
| 1200 | f2 = MIN(f2,1.0_wp) |
---|
[1496] | 1201 | |
---|
[1976] | 1202 | r_soil(j,i) = r_soil_min(j,i) / f2 |
---|
[1496] | 1203 | |
---|
| 1204 | ! |
---|
[1976] | 1205 | !-- Calculate the maximum possible liquid water amount on plants and |
---|
| 1206 | !-- bare surface. For vegetated surfaces, a maximum depth of 0.2 mm is |
---|
| 1207 | !-- assumed, while paved surfaces might hold up 1 mm of water. The |
---|
| 1208 | !-- liquid water fraction for paved surfaces is calculated after |
---|
| 1209 | !-- Noilhan & Planton (1989), while the ECMWF formulation is used for |
---|
| 1210 | !-- vegetated surfaces and bare soils. |
---|
| 1211 | IF ( pave_surface(j,i) ) THEN |
---|
| 1212 | m_liq_eb_max = m_max_depth * 5.0_wp |
---|
| 1213 | c_liq(j,i) = MIN( 1.0_wp, (m_liq_eb(j,i) / m_liq_eb_max)**0.67 ) |
---|
| 1214 | ELSE |
---|
| 1215 | m_liq_eb_max = m_max_depth * ( c_veg(j,i) * lai(j,i) & |
---|
| 1216 | + (1.0_wp - c_veg(j,i)) ) |
---|
| 1217 | c_liq(j,i) = MIN( 1.0_wp, m_liq_eb(j,i) / m_liq_eb_max ) |
---|
| 1218 | ENDIF |
---|
[1496] | 1219 | |
---|
[1551] | 1220 | ! |
---|
[1976] | 1221 | !-- Calculate saturation specific humidity |
---|
| 1222 | q_s = 0.622_wp * e_s / surface_pressure |
---|
[1500] | 1223 | |
---|
| 1224 | ! |
---|
[1976] | 1225 | !-- In case of dewfall, set evapotranspiration to zero |
---|
| 1226 | !-- All super-saturated water is then removed from the air |
---|
| 1227 | IF ( humidity .AND. q_s <= qv1 ) THEN |
---|
| 1228 | r_canopy(j,i) = 0.0_wp |
---|
| 1229 | r_soil(j,i) = 0.0_wp |
---|
| 1230 | ENDIF |
---|
[1496] | 1231 | |
---|
| 1232 | ! |
---|
[1976] | 1233 | !-- Calculate coefficients for the total evapotranspiration |
---|
| 1234 | !-- In case of water surface, set vegetation and soil fluxes to zero. |
---|
| 1235 | !-- For pavements, only evaporation of liquid water is possible. |
---|
| 1236 | IF ( water_surface(j,i) ) THEN |
---|
| 1237 | f_qsws_veg = 0.0_wp |
---|
| 1238 | f_qsws_soil = 0.0_wp |
---|
| 1239 | f_qsws_liq = rho_lv / r_a(j,i) |
---|
| 1240 | ELSEIF ( pave_surface (j,i) ) THEN |
---|
| 1241 | f_qsws_veg = 0.0_wp |
---|
| 1242 | f_qsws_soil = 0.0_wp |
---|
| 1243 | f_qsws_liq = rho_lv * c_liq(j,i) / r_a(j,i) |
---|
| 1244 | ELSE |
---|
| 1245 | f_qsws_veg = rho_lv * c_veg(j,i) * (1.0_wp - c_liq(j,i))/ & |
---|
| 1246 | (r_a(j,i) + r_canopy(j,i)) |
---|
| 1247 | f_qsws_soil = rho_lv * (1.0_wp - c_veg(j,i)) / (r_a(j,i) + & |
---|
| 1248 | r_soil(j,i)) |
---|
| 1249 | f_qsws_liq = rho_lv * c_veg(j,i) * c_liq(j,i) / r_a(j,i) |
---|
| 1250 | ENDIF |
---|
[1500] | 1251 | ! |
---|
[1976] | 1252 | !-- If soil moisture is below wilting point, plants do no longer |
---|
| 1253 | !-- transpirate. |
---|
| 1254 | ! IF ( m_soil(k,j,i) < m_wilt(j,i) ) THEN |
---|
| 1255 | ! f_qsws_veg = 0.0_wp |
---|
| 1256 | ! ENDIF |
---|
[1496] | 1257 | |
---|
[1976] | 1258 | f_shf = rho_cp / r_a(j,i) |
---|
| 1259 | f_qsws = f_qsws_veg + f_qsws_soil + f_qsws_liq |
---|
[1551] | 1260 | |
---|
| 1261 | ! |
---|
[1976] | 1262 | !-- Calculate derivative of q_s for Taylor series expansion |
---|
| 1263 | e_s_dt = e_s * ( 17.269_wp / (t_surface(j,i) - 35.86_wp) - & |
---|
| 1264 | 17.269_wp*(t_surface(j,i) - 273.16_wp) & |
---|
| 1265 | / (t_surface(j,i) - 35.86_wp)**2 ) |
---|
[1496] | 1266 | |
---|
[1976] | 1267 | dq_s_dt = 0.622_wp * e_s_dt / surface_pressure |
---|
[1496] | 1268 | |
---|
| 1269 | ! |
---|
[1976] | 1270 | !-- Add LW up so that it can be removed in prognostic equation |
---|
| 1271 | rad_net_l(j,i) = rad_net(j,i) + rad_lw_out(nzb,j,i) |
---|
[1496] | 1272 | |
---|
[1788] | 1273 | ! |
---|
[1976] | 1274 | !-- Calculate new skin temperature |
---|
| 1275 | IF ( humidity ) THEN |
---|
[1496] | 1276 | |
---|
| 1277 | ! |
---|
[1976] | 1278 | !-- Numerator of the prognostic equation |
---|
| 1279 | coef_1 = rad_net_l(j,i) + rad_lw_out_change_0(j,i) & |
---|
| 1280 | * t_surface(j,i) - rad_lw_out(nzb,j,i) & |
---|
| 1281 | + f_shf * pt1 + f_qsws * ( qv1 - q_s & |
---|
| 1282 | + dq_s_dt * t_surface(j,i) ) + lambda_surface & |
---|
| 1283 | * t_soil(nzb_soil,j,i) |
---|
[1496] | 1284 | |
---|
[1691] | 1285 | ! |
---|
[1976] | 1286 | !-- Denominator of the prognostic equation |
---|
| 1287 | coef_2 = rad_lw_out_change_0(j,i) + f_qsws * dq_s_dt & |
---|
| 1288 | + lambda_surface + f_shf / exn |
---|
| 1289 | ELSE |
---|
[1691] | 1290 | |
---|
| 1291 | ! |
---|
[1976] | 1292 | !-- Numerator of the prognostic equation |
---|
| 1293 | coef_1 = rad_net_l(j,i) + rad_lw_out_change_0(j,i) & |
---|
| 1294 | * t_surface(j,i) - rad_lw_out(nzb,j,i) & |
---|
| 1295 | + f_shf * pt1 + lambda_surface & |
---|
| 1296 | * t_soil(nzb_soil,j,i) |
---|
[1500] | 1297 | |
---|
| 1298 | ! |
---|
[1976] | 1299 | !-- Denominator of the prognostic equation |
---|
| 1300 | coef_2 = rad_lw_out_change_0(j,i) + lambda_surface + f_shf / exn |
---|
[1500] | 1301 | |
---|
[1976] | 1302 | ENDIF |
---|
[1691] | 1303 | |
---|
[1976] | 1304 | tend = 0.0_wp |
---|
[1691] | 1305 | |
---|
| 1306 | ! |
---|
[1976] | 1307 | !-- Implicit solution when the surface layer has no heat capacity, |
---|
| 1308 | !-- otherwise use RK3 scheme. |
---|
| 1309 | t_surface_p(j,i) = ( coef_1 * dt_3d * tsc(2) + c_surface_tmp * & |
---|
| 1310 | t_surface(j,i) ) / ( c_surface_tmp + coef_2 & |
---|
[1788] | 1311 | * dt_3d * tsc(2) ) |
---|
| 1312 | |
---|
[1496] | 1313 | ! |
---|
[1976] | 1314 | !-- Add RK3 term |
---|
| 1315 | IF ( c_surface_tmp /= 0.0_wp ) THEN |
---|
[1496] | 1316 | |
---|
[1976] | 1317 | t_surface_p(j,i) = t_surface_p(j,i) + dt_3d * tsc(3) & |
---|
| 1318 | * tt_surface_m(j,i) |
---|
[1691] | 1319 | |
---|
[1496] | 1320 | ! |
---|
[1976] | 1321 | !-- Calculate true tendency |
---|
| 1322 | tend = (t_surface_p(j,i) - t_surface(j,i) - dt_3d * tsc(3) & |
---|
| 1323 | * tt_surface_m(j,i)) / (dt_3d * tsc(2)) |
---|
[1496] | 1324 | ! |
---|
[1976] | 1325 | !-- Calculate t_surface tendencies for the next Runge-Kutta step |
---|
| 1326 | IF ( timestep_scheme(1:5) == 'runge' ) THEN |
---|
| 1327 | IF ( intermediate_timestep_count == 1 ) THEN |
---|
| 1328 | tt_surface_m(j,i) = tend |
---|
| 1329 | ELSEIF ( intermediate_timestep_count < & |
---|
| 1330 | intermediate_timestep_count_max ) THEN |
---|
| 1331 | tt_surface_m(j,i) = -9.5625_wp * tend + 5.3125_wp & |
---|
| 1332 | * tt_surface_m(j,i) |
---|
[1496] | 1333 | ENDIF |
---|
| 1334 | ENDIF |
---|
[1976] | 1335 | ENDIF |
---|
[1496] | 1336 | |
---|
| 1337 | ! |
---|
[1976] | 1338 | !-- In case of fast changes in the skin temperature, it is possible to |
---|
| 1339 | !-- update the radiative fluxes independently from the prescribed |
---|
| 1340 | !-- radiation call frequency. This effectively prevents oscillations, |
---|
| 1341 | !-- especially when setting skip_time_do_radiation /= 0. The threshold |
---|
| 1342 | !-- value of 0.2 used here is just a first guess. This method should be |
---|
| 1343 | !-- revised in the future as tests have shown that the threshold is |
---|
| 1344 | !-- often reached, when no oscillations would occur (causes immense |
---|
| 1345 | !-- computing time for the radiation code). |
---|
| 1346 | IF ( ABS( t_surface_p(j,i) - t_surface(j,i) ) > 0.2_wp .AND. & |
---|
| 1347 | unscheduled_radiation_calls ) THEN |
---|
| 1348 | force_radiation_call_l = .TRUE. |
---|
| 1349 | ENDIF |
---|
[1691] | 1350 | |
---|
[1976] | 1351 | pt(k,j,i) = t_surface_p(j,i) / exn |
---|
[1691] | 1352 | |
---|
| 1353 | ! |
---|
[1976] | 1354 | !-- Calculate fluxes |
---|
| 1355 | rad_net_l(j,i) = rad_net_l(j,i) + rad_lw_out_change_0(j,i) & |
---|
| 1356 | * t_surface(j,i) - rad_lw_out(nzb,j,i) & |
---|
| 1357 | - rad_lw_out_change_0(j,i) * t_surface_p(j,i) |
---|
[1691] | 1358 | |
---|
[1976] | 1359 | rad_net(j,i) = rad_net_l(j,i) |
---|
| 1360 | rad_lw_out(nzb,j,i) = rad_lw_out(nzb,j,i) + rad_lw_out_change_0(j,i) & |
---|
| 1361 | * ( t_surface_p(j,i) - t_surface(j,i) ) |
---|
[1691] | 1362 | |
---|
[1976] | 1363 | ghf_eb(j,i) = lambda_surface * (t_surface_p(j,i) & |
---|
| 1364 | - t_soil(nzb_soil,j,i)) |
---|
[1496] | 1365 | |
---|
[1976] | 1366 | shf_eb(j,i) = - f_shf * ( pt1 - pt(k,j,i) ) |
---|
[1691] | 1367 | |
---|
[1976] | 1368 | shf(j,i) = shf_eb(j,i) / rho_cp |
---|
[1691] | 1369 | |
---|
[1976] | 1370 | IF ( humidity ) THEN |
---|
| 1371 | qsws_eb(j,i) = - f_qsws * ( qv1 - q_s + dq_s_dt & |
---|
| 1372 | * t_surface(j,i) - dq_s_dt * t_surface_p(j,i) ) |
---|
[1496] | 1373 | |
---|
[1976] | 1374 | qsws(j,i) = qsws_eb(j,i) / rho_lv |
---|
[1691] | 1375 | |
---|
[1976] | 1376 | qsws_veg_eb(j,i) = - f_qsws_veg * ( qv1 - q_s & |
---|
| 1377 | + dq_s_dt * t_surface(j,i) - dq_s_dt & |
---|
| 1378 | * t_surface_p(j,i) ) |
---|
[1551] | 1379 | |
---|
[1976] | 1380 | qsws_soil_eb(j,i) = - f_qsws_soil * ( qv1 - q_s & |
---|
| 1381 | + dq_s_dt * t_surface(j,i) - dq_s_dt & |
---|
| 1382 | * t_surface_p(j,i) ) |
---|
[1551] | 1383 | |
---|
[1976] | 1384 | qsws_liq_eb(j,i) = - f_qsws_liq * ( qv1 - q_s & |
---|
| 1385 | + dq_s_dt * t_surface(j,i) - dq_s_dt & |
---|
| 1386 | * t_surface_p(j,i) ) |
---|
| 1387 | ENDIF |
---|
[1788] | 1388 | |
---|
[1500] | 1389 | ! |
---|
[1976] | 1390 | !-- Calculate the true surface resistance |
---|
| 1391 | IF ( qsws_eb(j,i) == 0.0_wp ) THEN |
---|
| 1392 | r_s(j,i) = 1.0E10_wp |
---|
| 1393 | ELSE |
---|
| 1394 | r_s(j,i) = - rho_lv * ( qv1 - q_s + dq_s_dt & |
---|
| 1395 | * t_surface(j,i) - dq_s_dt * t_surface_p(j,i) ) & |
---|
| 1396 | / qsws_eb(j,i) - r_a(j,i) |
---|
| 1397 | ENDIF |
---|
[1496] | 1398 | |
---|
| 1399 | ! |
---|
[1976] | 1400 | !-- Calculate change in liquid water reservoir due to dew fall or |
---|
| 1401 | !-- evaporation of liquid water |
---|
| 1402 | IF ( humidity ) THEN |
---|
[1496] | 1403 | ! |
---|
[1976] | 1404 | !-- If precipitation is activated, add rain water to qsws_liq_eb |
---|
| 1405 | !-- and qsws_soil_eb according the the vegetation coverage. |
---|
| 1406 | !-- precipitation_rate is given in mm. |
---|
| 1407 | IF ( precipitation ) THEN |
---|
[1571] | 1408 | |
---|
| 1409 | ! |
---|
[1976] | 1410 | !-- Add precipitation to liquid water reservoir, if possible. |
---|
| 1411 | !-- Otherwise, add the water to soil. In case of |
---|
| 1412 | !-- pavements, the exceeding water amount is implicitely removed |
---|
| 1413 | !-- as runoff as qsws_soil_eb is then not used in the soil model |
---|
| 1414 | IF ( m_liq_eb(j,i) /= m_liq_eb_max ) THEN |
---|
| 1415 | qsws_liq_eb(j,i) = qsws_liq_eb(j,i) & |
---|
| 1416 | + c_veg(j,i) * prr(k,j,i) * hyrho(k) & |
---|
| 1417 | * 0.001_wp * rho_l * l_v |
---|
| 1418 | ELSE |
---|
| 1419 | qsws_soil_eb(j,i) = qsws_soil_eb(j,i) & |
---|
| 1420 | + c_veg(j,i) * prr(k,j,i) * hyrho(k) & |
---|
| 1421 | * 0.001_wp * rho_l * l_v |
---|
[1496] | 1422 | ENDIF |
---|
[1691] | 1423 | |
---|
[1976] | 1424 | !-- Add precipitation to bare soil according to the bare soil |
---|
| 1425 | !-- coverage. |
---|
| 1426 | qsws_soil_eb(j,i) = qsws_soil_eb(j,i) + (1.0_wp & |
---|
| 1427 | - c_veg(j,i)) * prr(k,j,i) * hyrho(k) & |
---|
| 1428 | * 0.001_wp * rho_l * l_v |
---|
| 1429 | ENDIF |
---|
| 1430 | |
---|
[1500] | 1431 | ! |
---|
[1976] | 1432 | !-- If the air is saturated, check the reservoir water level |
---|
| 1433 | IF ( qsws_eb(j,i) < 0.0_wp ) THEN |
---|
[1691] | 1434 | |
---|
[1500] | 1435 | ! |
---|
[1976] | 1436 | !-- Check if reservoir is full (avoid values > m_liq_eb_max) |
---|
| 1437 | !-- In that case, qsws_liq_eb goes to qsws_soil_eb. In this |
---|
| 1438 | !-- case qsws_veg_eb is zero anyway (because c_liq = 1), |
---|
| 1439 | !-- so that tend is zero and no further check is needed |
---|
| 1440 | IF ( m_liq_eb(j,i) == m_liq_eb_max ) THEN |
---|
| 1441 | qsws_soil_eb(j,i) = qsws_soil_eb(j,i) & |
---|
| 1442 | + qsws_liq_eb(j,i) |
---|
[1949] | 1443 | |
---|
[1976] | 1444 | qsws_liq_eb(j,i) = 0.0_wp |
---|
| 1445 | ENDIF |
---|
[1496] | 1446 | |
---|
| 1447 | ! |
---|
[1976] | 1448 | !-- In case qsws_veg_eb becomes negative (unphysical behavior), |
---|
| 1449 | !-- let the water enter the liquid water reservoir as dew on the |
---|
| 1450 | !-- plant |
---|
| 1451 | IF ( qsws_veg_eb(j,i) < 0.0_wp ) THEN |
---|
| 1452 | qsws_liq_eb(j,i) = qsws_liq_eb(j,i) + qsws_veg_eb(j,i) |
---|
| 1453 | qsws_veg_eb(j,i) = 0.0_wp |
---|
| 1454 | ENDIF |
---|
| 1455 | ENDIF |
---|
[1496] | 1456 | |
---|
[1976] | 1457 | tend = - qsws_liq_eb(j,i) * drho_l_lv |
---|
| 1458 | m_liq_eb_p(j,i) = m_liq_eb(j,i) + dt_3d * ( tsc(2) * tend & |
---|
| 1459 | + tsc(3) * tm_liq_eb_m(j,i) ) |
---|
[1496] | 1460 | |
---|
| 1461 | ! |
---|
[1976] | 1462 | !-- Check if reservoir is overfull -> reduce to maximum |
---|
| 1463 | !-- (conservation of water is violated here) |
---|
| 1464 | m_liq_eb_p(j,i) = MIN(m_liq_eb_p(j,i),m_liq_eb_max) |
---|
[1496] | 1465 | |
---|
| 1466 | ! |
---|
[1976] | 1467 | !-- Check if reservoir is empty (avoid values < 0.0) |
---|
| 1468 | !-- (conservation of water is violated here) |
---|
| 1469 | m_liq_eb_p(j,i) = MAX(m_liq_eb_p(j,i),0.0_wp) |
---|
[1496] | 1470 | |
---|
| 1471 | |
---|
| 1472 | ! |
---|
[1976] | 1473 | !-- Calculate m_liq_eb tendencies for the next Runge-Kutta step |
---|
| 1474 | IF ( timestep_scheme(1:5) == 'runge' ) THEN |
---|
| 1475 | IF ( intermediate_timestep_count == 1 ) THEN |
---|
| 1476 | tm_liq_eb_m(j,i) = tend |
---|
| 1477 | ELSEIF ( intermediate_timestep_count < & |
---|
| 1478 | intermediate_timestep_count_max ) THEN |
---|
| 1479 | tm_liq_eb_m(j,i) = -9.5625_wp * tend + 5.3125_wp & |
---|
| 1480 | * tm_liq_eb_m(j,i) |
---|
[1496] | 1481 | ENDIF |
---|
[1500] | 1482 | ENDIF |
---|
| 1483 | |
---|
[1976] | 1484 | ENDIF |
---|
| 1485 | |
---|
[1500] | 1486 | ENDDO |
---|
[1976] | 1487 | ENDDO |
---|
[1496] | 1488 | |
---|
[1691] | 1489 | ! |
---|
[1976] | 1490 | !-- Make a logical OR for all processes. Force radiation call if at |
---|
| 1491 | !-- least one processor reached the threshold change in skin temperature |
---|
| 1492 | IF ( unscheduled_radiation_calls .AND. intermediate_timestep_count & |
---|
| 1493 | == intermediate_timestep_count_max-1 ) THEN |
---|
[1697] | 1494 | #if defined( __parallel ) |
---|
[1976] | 1495 | IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr ) |
---|
| 1496 | CALL MPI_ALLREDUCE( force_radiation_call_l, force_radiation_call, & |
---|
| 1497 | 1, MPI_LOGICAL, MPI_LOR, comm2d, ierr ) |
---|
[1691] | 1498 | #else |
---|
[1976] | 1499 | force_radiation_call = force_radiation_call_l |
---|
[1691] | 1500 | #endif |
---|
[1976] | 1501 | force_radiation_call_l = .FALSE. |
---|
| 1502 | ENDIF |
---|
[1691] | 1503 | |
---|
[1788] | 1504 | ! |
---|
[1976] | 1505 | !-- Calculate surface specific humidity |
---|
| 1506 | IF ( humidity ) THEN |
---|
| 1507 | CALL calc_q_surface |
---|
| 1508 | ENDIF |
---|
[1691] | 1509 | |
---|
[1788] | 1510 | ! |
---|
[1976] | 1511 | !-- Calculate new roughness lengths (for water surfaces only) |
---|
| 1512 | CALL calc_z0_water_surface |
---|
[1788] | 1513 | |
---|
| 1514 | |
---|
[1976] | 1515 | END SUBROUTINE lsm_energy_balance |
---|
[1496] | 1516 | |
---|
[1976] | 1517 | |
---|
[1817] | 1518 | !------------------------------------------------------------------------------! |
---|
| 1519 | ! Description: |
---|
| 1520 | ! ------------ |
---|
| 1521 | !> Header output for land surface model |
---|
| 1522 | !------------------------------------------------------------------------------! |
---|
| 1523 | SUBROUTINE lsm_header ( io ) |
---|
[1496] | 1524 | |
---|
[1817] | 1525 | |
---|
| 1526 | IMPLICIT NONE |
---|
| 1527 | |
---|
| 1528 | CHARACTER (LEN=86) :: t_soil_chr !< String for soil temperature profile |
---|
| 1529 | CHARACTER (LEN=86) :: roots_chr !< String for root profile |
---|
| 1530 | CHARACTER (LEN=86) :: vertical_index_chr !< String for the vertical index |
---|
| 1531 | CHARACTER (LEN=86) :: m_soil_chr !< String for soil moisture |
---|
| 1532 | CHARACTER (LEN=86) :: soil_depth_chr !< String for soil depth |
---|
| 1533 | CHARACTER (LEN=10) :: coor_chr !< Temporary string |
---|
| 1534 | |
---|
| 1535 | INTEGER(iwp) :: i !< Loop index over soil layers |
---|
| 1536 | |
---|
| 1537 | INTEGER(iwp), INTENT(IN) :: io !< Unit of the output file |
---|
| 1538 | |
---|
| 1539 | t_soil_chr = '' |
---|
| 1540 | m_soil_chr = '' |
---|
| 1541 | soil_depth_chr = '' |
---|
| 1542 | roots_chr = '' |
---|
| 1543 | vertical_index_chr = '' |
---|
| 1544 | |
---|
| 1545 | i = 1 |
---|
| 1546 | DO i = nzb_soil, nzt_soil |
---|
| 1547 | WRITE (coor_chr,'(F10.2,7X)') soil_temperature(i) |
---|
| 1548 | t_soil_chr = TRIM( t_soil_chr ) // ' ' // TRIM( coor_chr ) |
---|
| 1549 | |
---|
| 1550 | WRITE (coor_chr,'(F10.2,7X)') soil_moisture(i) |
---|
| 1551 | m_soil_chr = TRIM( m_soil_chr ) // ' ' // TRIM( coor_chr ) |
---|
| 1552 | |
---|
| 1553 | WRITE (coor_chr,'(F10.2,7X)') - zs(i) |
---|
| 1554 | soil_depth_chr = TRIM( soil_depth_chr ) // ' ' // TRIM( coor_chr ) |
---|
| 1555 | |
---|
| 1556 | WRITE (coor_chr,'(F10.2,7X)') root_fraction(i) |
---|
| 1557 | roots_chr = TRIM( roots_chr ) // ' ' // TRIM( coor_chr ) |
---|
| 1558 | |
---|
| 1559 | WRITE (coor_chr,'(I10,7X)') i |
---|
| 1560 | vertical_index_chr = TRIM( vertical_index_chr ) // ' ' // & |
---|
| 1561 | TRIM( coor_chr ) |
---|
| 1562 | ENDDO |
---|
| 1563 | |
---|
| 1564 | ! |
---|
| 1565 | !-- Write land surface model header |
---|
| 1566 | WRITE( io, 1 ) |
---|
| 1567 | IF ( conserve_water_content ) THEN |
---|
| 1568 | WRITE( io, 2 ) |
---|
| 1569 | ELSE |
---|
| 1570 | WRITE( io, 3 ) |
---|
| 1571 | ENDIF |
---|
| 1572 | |
---|
| 1573 | WRITE( io, 4 ) TRIM( veg_type_name(veg_type) ), & |
---|
| 1574 | TRIM (soil_type_name(soil_type) ) |
---|
| 1575 | WRITE( io, 5 ) TRIM( soil_depth_chr ), TRIM( t_soil_chr ), & |
---|
| 1576 | TRIM( m_soil_chr ), TRIM( roots_chr ), & |
---|
| 1577 | TRIM( vertical_index_chr ) |
---|
| 1578 | |
---|
| 1579 | 1 FORMAT (//' Land surface model information:'/ & |
---|
| 1580 | ' ------------------------------'/) |
---|
[1826] | 1581 | 2 FORMAT (' --> Soil bottom is closed (water content is conserved', & |
---|
| 1582 | ', default)') |
---|
| 1583 | 3 FORMAT (' --> Soil bottom is open (water content is not conserved)') |
---|
| 1584 | 4 FORMAT (' --> Land surface type : ',A,/ & |
---|
| 1585 | ' --> Soil porosity type : ',A) |
---|
[1817] | 1586 | 5 FORMAT (/' Initial soil temperature and moisture profile:'// & |
---|
| 1587 | ' Height: ',A,' m'/ & |
---|
| 1588 | ' Temperature: ',A,' K'/ & |
---|
| 1589 | ' Moisture: ',A,' m**3/m**3'/ & |
---|
| 1590 | ' Root fraction: ',A,' '/ & |
---|
[1826] | 1591 | ' Grid point: ',A) |
---|
[1817] | 1592 | |
---|
| 1593 | END SUBROUTINE lsm_header |
---|
| 1594 | |
---|
| 1595 | |
---|
[1496] | 1596 | !------------------------------------------------------------------------------! |
---|
| 1597 | ! Description: |
---|
| 1598 | ! ------------ |
---|
[1817] | 1599 | !> Initialization of the land surface model |
---|
| 1600 | !------------------------------------------------------------------------------! |
---|
| 1601 | SUBROUTINE lsm_init |
---|
| 1602 | |
---|
| 1603 | |
---|
| 1604 | IMPLICIT NONE |
---|
| 1605 | |
---|
| 1606 | INTEGER(iwp) :: i !< running index |
---|
| 1607 | INTEGER(iwp) :: j !< running index |
---|
| 1608 | INTEGER(iwp) :: k !< running index |
---|
| 1609 | |
---|
| 1610 | REAL(wp) :: pt1 !< potential temperature at first grid level |
---|
| 1611 | |
---|
| 1612 | |
---|
| 1613 | ! |
---|
| 1614 | !-- Calculate Exner function |
---|
| 1615 | exn = ( surface_pressure / 1000.0_wp )**0.286_wp |
---|
| 1616 | |
---|
| 1617 | |
---|
| 1618 | ! |
---|
| 1619 | !-- If no cloud physics is used, rho_surface has not been calculated before |
---|
| 1620 | IF ( .NOT. cloud_physics ) THEN |
---|
| 1621 | rho_surface = surface_pressure * 100.0_wp / ( r_d * pt_surface * exn ) |
---|
| 1622 | ENDIF |
---|
| 1623 | |
---|
| 1624 | ! |
---|
| 1625 | !-- Calculate frequently used parameters |
---|
| 1626 | rho_cp = cp * rho_surface |
---|
| 1627 | rd_d_rv = r_d / r_v |
---|
| 1628 | rho_lv = rho_surface * l_v |
---|
| 1629 | drho_l_lv = 1.0_wp / (rho_l * l_v) |
---|
| 1630 | |
---|
| 1631 | ! |
---|
| 1632 | !-- Set inital values for prognostic quantities |
---|
| 1633 | tt_surface_m = 0.0_wp |
---|
| 1634 | tt_soil_m = 0.0_wp |
---|
| 1635 | tm_soil_m = 0.0_wp |
---|
| 1636 | tm_liq_eb_m = 0.0_wp |
---|
| 1637 | c_liq = 0.0_wp |
---|
| 1638 | |
---|
| 1639 | ghf_eb = 0.0_wp |
---|
| 1640 | shf_eb = rho_cp * shf |
---|
| 1641 | |
---|
| 1642 | IF ( humidity ) THEN |
---|
[1826] | 1643 | qsws_eb = rho_lv * qsws |
---|
[1817] | 1644 | ELSE |
---|
| 1645 | qsws_eb = 0.0_wp |
---|
| 1646 | ENDIF |
---|
| 1647 | |
---|
| 1648 | qsws_liq_eb = 0.0_wp |
---|
| 1649 | qsws_soil_eb = 0.0_wp |
---|
| 1650 | qsws_veg_eb = 0.0_wp |
---|
| 1651 | |
---|
| 1652 | r_a = 50.0_wp |
---|
| 1653 | r_s = 50.0_wp |
---|
| 1654 | r_canopy = 0.0_wp |
---|
| 1655 | r_soil = 0.0_wp |
---|
| 1656 | |
---|
| 1657 | ! |
---|
| 1658 | !-- Allocate 3D soil model arrays |
---|
| 1659 | ALLOCATE ( root_fr(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) ) |
---|
| 1660 | ALLOCATE ( lambda_h(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) ) |
---|
| 1661 | ALLOCATE ( rho_c_total(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) ) |
---|
| 1662 | |
---|
| 1663 | lambda_h = 0.0_wp |
---|
| 1664 | ! |
---|
| 1665 | !-- If required, allocate humidity-related variables for the soil model |
---|
| 1666 | IF ( humidity ) THEN |
---|
| 1667 | ALLOCATE ( lambda_w(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) ) |
---|
| 1668 | ALLOCATE ( gamma_w(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) ) |
---|
| 1669 | |
---|
| 1670 | lambda_w = 0.0_wp |
---|
| 1671 | ENDIF |
---|
| 1672 | |
---|
| 1673 | ! |
---|
| 1674 | !-- Calculate grid spacings. Temperature and moisture are defined at |
---|
| 1675 | !-- the edges of the soil layers (_stag), whereas gradients/fluxes are defined |
---|
| 1676 | !-- at the centers |
---|
| 1677 | dz_soil(nzb_soil) = zs(nzb_soil) |
---|
| 1678 | |
---|
| 1679 | DO k = nzb_soil+1, nzt_soil |
---|
| 1680 | dz_soil(k) = zs(k) - zs(k-1) |
---|
| 1681 | ENDDO |
---|
| 1682 | dz_soil(nzt_soil+1) = dz_soil(nzt_soil) |
---|
| 1683 | |
---|
| 1684 | DO k = nzb_soil, nzt_soil-1 |
---|
| 1685 | dz_soil_stag(k) = 0.5_wp * (dz_soil(k+1) + dz_soil(k)) |
---|
| 1686 | ENDDO |
---|
| 1687 | dz_soil_stag(nzt_soil) = dz_soil(nzt_soil) |
---|
| 1688 | |
---|
| 1689 | ddz_soil = 1.0_wp / dz_soil |
---|
| 1690 | ddz_soil_stag = 1.0_wp / dz_soil_stag |
---|
| 1691 | |
---|
| 1692 | ! |
---|
| 1693 | !-- Initialize standard soil types. It is possible to overwrite each |
---|
| 1694 | !-- parameter by setting the respecticy NAMELIST variable to a |
---|
| 1695 | !-- value /= 9999999.9. |
---|
| 1696 | IF ( soil_type /= 0 ) THEN |
---|
| 1697 | |
---|
| 1698 | IF ( alpha_vangenuchten == 9999999.9_wp ) THEN |
---|
| 1699 | alpha_vangenuchten = soil_pars(0,soil_type) |
---|
| 1700 | ENDIF |
---|
| 1701 | |
---|
| 1702 | IF ( l_vangenuchten == 9999999.9_wp ) THEN |
---|
| 1703 | l_vangenuchten = soil_pars(1,soil_type) |
---|
| 1704 | ENDIF |
---|
| 1705 | |
---|
| 1706 | IF ( n_vangenuchten == 9999999.9_wp ) THEN |
---|
| 1707 | n_vangenuchten = soil_pars(2,soil_type) |
---|
| 1708 | ENDIF |
---|
| 1709 | |
---|
| 1710 | IF ( hydraulic_conductivity == 9999999.9_wp ) THEN |
---|
| 1711 | hydraulic_conductivity = soil_pars(3,soil_type) |
---|
| 1712 | ENDIF |
---|
| 1713 | |
---|
| 1714 | IF ( saturation_moisture == 9999999.9_wp ) THEN |
---|
| 1715 | saturation_moisture = m_soil_pars(0,soil_type) |
---|
| 1716 | ENDIF |
---|
| 1717 | |
---|
| 1718 | IF ( field_capacity == 9999999.9_wp ) THEN |
---|
| 1719 | field_capacity = m_soil_pars(1,soil_type) |
---|
| 1720 | ENDIF |
---|
| 1721 | |
---|
| 1722 | IF ( wilting_point == 9999999.9_wp ) THEN |
---|
| 1723 | wilting_point = m_soil_pars(2,soil_type) |
---|
| 1724 | ENDIF |
---|
| 1725 | |
---|
| 1726 | IF ( residual_moisture == 9999999.9_wp ) THEN |
---|
| 1727 | residual_moisture = m_soil_pars(3,soil_type) |
---|
| 1728 | ENDIF |
---|
| 1729 | |
---|
| 1730 | ENDIF |
---|
| 1731 | |
---|
| 1732 | ! |
---|
| 1733 | !-- Map values to the respective 2D arrays |
---|
| 1734 | alpha_vg = alpha_vangenuchten |
---|
| 1735 | l_vg = l_vangenuchten |
---|
| 1736 | n_vg = n_vangenuchten |
---|
| 1737 | gamma_w_sat = hydraulic_conductivity |
---|
| 1738 | m_sat = saturation_moisture |
---|
| 1739 | m_fc = field_capacity |
---|
| 1740 | m_wilt = wilting_point |
---|
| 1741 | m_res = residual_moisture |
---|
| 1742 | r_soil_min = min_soil_resistance |
---|
| 1743 | |
---|
| 1744 | ! |
---|
| 1745 | !-- Initial run actions |
---|
| 1746 | IF ( TRIM( initializing_actions ) /= 'read_restart_data' ) THEN |
---|
| 1747 | |
---|
| 1748 | t_soil = 0.0_wp |
---|
| 1749 | m_liq_eb = 0.0_wp |
---|
| 1750 | m_soil = 0.0_wp |
---|
| 1751 | |
---|
| 1752 | ! |
---|
| 1753 | !-- Map user settings of T and q for each soil layer |
---|
| 1754 | !-- (make sure that the soil moisture does not drop below the permanent |
---|
| 1755 | !-- wilting point) -> problems with devision by zero) |
---|
| 1756 | DO k = nzb_soil, nzt_soil |
---|
| 1757 | t_soil(k,:,:) = soil_temperature(k) |
---|
| 1758 | m_soil(k,:,:) = MAX(soil_moisture(k),m_wilt(:,:)) |
---|
| 1759 | soil_moisture(k) = MAX(soil_moisture(k),wilting_point) |
---|
| 1760 | ENDDO |
---|
| 1761 | t_soil(nzt_soil+1,:,:) = soil_temperature(nzt_soil+1) |
---|
| 1762 | |
---|
| 1763 | ! |
---|
| 1764 | !-- Calculate surface temperature |
---|
| 1765 | t_surface = pt_surface * exn |
---|
| 1766 | |
---|
| 1767 | ! |
---|
| 1768 | !-- Set artifical values for ts and us so that r_a has its initial value |
---|
| 1769 | !-- for the first time step |
---|
| 1770 | DO i = nxlg, nxrg |
---|
| 1771 | DO j = nysg, nyng |
---|
| 1772 | k = nzb_s_inner(j,i) |
---|
| 1773 | |
---|
| 1774 | IF ( cloud_physics ) THEN |
---|
| 1775 | pt1 = pt(k+1,j,i) + l_d_cp * pt_d_t(k+1) * ql(k+1,j,i) |
---|
| 1776 | ELSE |
---|
| 1777 | pt1 = pt(k+1,j,i) |
---|
| 1778 | ENDIF |
---|
| 1779 | |
---|
| 1780 | ! |
---|
| 1781 | !-- Assure that r_a cannot be zero at model start |
---|
| 1782 | IF ( pt1 == pt(k,j,i) ) pt1 = pt1 + 1.0E-10_wp |
---|
| 1783 | |
---|
| 1784 | us(j,i) = 0.1_wp |
---|
| 1785 | ts(j,i) = (pt1 - pt(k,j,i)) / r_a(j,i) |
---|
| 1786 | shf(j,i) = - us(j,i) * ts(j,i) |
---|
| 1787 | ENDDO |
---|
| 1788 | ENDDO |
---|
| 1789 | |
---|
| 1790 | ! |
---|
| 1791 | !-- Actions for restart runs |
---|
| 1792 | ELSE |
---|
| 1793 | |
---|
| 1794 | DO i = nxlg, nxrg |
---|
| 1795 | DO j = nysg, nyng |
---|
| 1796 | k = nzb_s_inner(j,i) |
---|
| 1797 | t_surface(j,i) = pt(k,j,i) * exn |
---|
| 1798 | ENDDO |
---|
| 1799 | ENDDO |
---|
| 1800 | |
---|
| 1801 | ENDIF |
---|
| 1802 | |
---|
| 1803 | DO k = nzb_soil, nzt_soil |
---|
| 1804 | root_fr(k,:,:) = root_fraction(k) |
---|
| 1805 | ENDDO |
---|
| 1806 | |
---|
| 1807 | IF ( veg_type /= 0 ) THEN |
---|
| 1808 | IF ( min_canopy_resistance == 9999999.9_wp ) THEN |
---|
| 1809 | min_canopy_resistance = veg_pars(0,veg_type) |
---|
| 1810 | ENDIF |
---|
| 1811 | IF ( leaf_area_index == 9999999.9_wp ) THEN |
---|
| 1812 | leaf_area_index = veg_pars(1,veg_type) |
---|
| 1813 | ENDIF |
---|
| 1814 | IF ( vegetation_coverage == 9999999.9_wp ) THEN |
---|
| 1815 | vegetation_coverage = veg_pars(2,veg_type) |
---|
| 1816 | ENDIF |
---|
| 1817 | IF ( canopy_resistance_coefficient == 9999999.9_wp ) THEN |
---|
| 1818 | canopy_resistance_coefficient= veg_pars(3,veg_type) |
---|
| 1819 | ENDIF |
---|
| 1820 | IF ( lambda_surface_stable == 9999999.9_wp ) THEN |
---|
| 1821 | lambda_surface_stable = surface_pars(0,veg_type) |
---|
| 1822 | ENDIF |
---|
| 1823 | IF ( lambda_surface_unstable == 9999999.9_wp ) THEN |
---|
| 1824 | lambda_surface_unstable = surface_pars(1,veg_type) |
---|
| 1825 | ENDIF |
---|
| 1826 | IF ( f_shortwave_incoming == 9999999.9_wp ) THEN |
---|
| 1827 | f_shortwave_incoming = surface_pars(2,veg_type) |
---|
| 1828 | ENDIF |
---|
| 1829 | IF ( z0_eb == 9999999.9_wp ) THEN |
---|
| 1830 | roughness_length = roughness_par(0,veg_type) |
---|
| 1831 | z0_eb = roughness_par(0,veg_type) |
---|
| 1832 | ENDIF |
---|
| 1833 | IF ( z0h_eb == 9999999.9_wp ) THEN |
---|
| 1834 | z0h_eb = roughness_par(1,veg_type) |
---|
| 1835 | ENDIF |
---|
| 1836 | IF ( z0q_eb == 9999999.9_wp ) THEN |
---|
| 1837 | z0q_eb = roughness_par(2,veg_type) |
---|
| 1838 | ENDIF |
---|
| 1839 | z0h_factor = z0h_eb / ( z0_eb + 1.0E-20_wp ) |
---|
| 1840 | |
---|
| 1841 | IF ( ANY( root_fraction == 9999999.9_wp ) ) THEN |
---|
| 1842 | DO k = nzb_soil, nzt_soil |
---|
| 1843 | root_fr(k,:,:) = root_distribution(k,veg_type) |
---|
| 1844 | root_fraction(k) = root_distribution(k,veg_type) |
---|
| 1845 | ENDDO |
---|
| 1846 | ENDIF |
---|
| 1847 | |
---|
| 1848 | ELSE |
---|
| 1849 | |
---|
| 1850 | IF ( z0_eb == 9999999.9_wp ) THEN |
---|
| 1851 | z0_eb = roughness_length |
---|
| 1852 | ENDIF |
---|
| 1853 | IF ( z0h_eb == 9999999.9_wp ) THEN |
---|
| 1854 | z0h_eb = z0_eb * z0h_factor |
---|
| 1855 | ENDIF |
---|
| 1856 | IF ( z0q_eb == 9999999.9_wp ) THEN |
---|
| 1857 | z0q_eb = z0_eb * z0h_factor |
---|
| 1858 | ENDIF |
---|
| 1859 | |
---|
| 1860 | ENDIF |
---|
| 1861 | |
---|
| 1862 | ! |
---|
| 1863 | !-- For surfaces covered with pavement, set depth of the pavement (with dry |
---|
| 1864 | !-- soil below). The depth must be greater than the first soil layer depth |
---|
| 1865 | IF ( veg_type == 20 ) THEN |
---|
| 1866 | IF ( pave_depth == 9999999.9_wp ) THEN |
---|
| 1867 | pave_depth = zs(nzb_soil) |
---|
| 1868 | ELSE |
---|
| 1869 | pave_depth = MAX( zs(nzb_soil), pave_depth ) |
---|
| 1870 | ENDIF |
---|
| 1871 | ENDIF |
---|
| 1872 | |
---|
| 1873 | ! |
---|
| 1874 | !-- Map vegetation and soil types to 2D array to allow for heterogeneous |
---|
| 1875 | !-- surfaces via user interface see below |
---|
| 1876 | veg_type_2d = veg_type |
---|
| 1877 | soil_type_2d = soil_type |
---|
| 1878 | |
---|
| 1879 | ! |
---|
| 1880 | !-- Map vegetation parameters to the respective 2D arrays |
---|
| 1881 | r_canopy_min = min_canopy_resistance |
---|
| 1882 | lai = leaf_area_index |
---|
| 1883 | c_veg = vegetation_coverage |
---|
| 1884 | g_d = canopy_resistance_coefficient |
---|
| 1885 | lambda_surface_s = lambda_surface_stable |
---|
| 1886 | lambda_surface_u = lambda_surface_unstable |
---|
| 1887 | f_sw_in = f_shortwave_incoming |
---|
| 1888 | z0 = z0_eb |
---|
| 1889 | z0h = z0h_eb |
---|
| 1890 | z0q = z0q_eb |
---|
| 1891 | |
---|
| 1892 | ! |
---|
| 1893 | !-- Possibly do user-defined actions (e.g. define heterogeneous land surface) |
---|
| 1894 | CALL user_init_land_surface |
---|
| 1895 | |
---|
| 1896 | ! |
---|
| 1897 | !-- Set flag parameter if vegetation type was set to a water surface. Also |
---|
| 1898 | !-- set temperature to a constant value in all "soil" layers. |
---|
| 1899 | DO i = nxlg, nxrg |
---|
| 1900 | DO j = nysg, nyng |
---|
| 1901 | IF ( veg_type_2d(j,i) == 14 .OR. veg_type_2d(j,i) == 15 ) THEN |
---|
| 1902 | water_surface(j,i) = .TRUE. |
---|
[1856] | 1903 | t_soil(:,j,i) = t_surface(j,i) |
---|
[1817] | 1904 | ELSEIF ( veg_type_2d(j,i) == 20 ) THEN |
---|
| 1905 | pave_surface(j,i) = .TRUE. |
---|
| 1906 | m_soil(:,j,i) = 0.0_wp |
---|
| 1907 | ENDIF |
---|
| 1908 | |
---|
| 1909 | ENDDO |
---|
| 1910 | ENDDO |
---|
| 1911 | |
---|
| 1912 | ! |
---|
| 1913 | !-- Calculate new roughness lengths (for water surfaces only) |
---|
| 1914 | CALL calc_z0_water_surface |
---|
| 1915 | |
---|
| 1916 | t_soil_p = t_soil |
---|
| 1917 | m_soil_p = m_soil |
---|
| 1918 | m_liq_eb_p = m_liq_eb |
---|
| 1919 | t_surface_p = t_surface |
---|
| 1920 | |
---|
| 1921 | |
---|
| 1922 | |
---|
| 1923 | !-- Store initial profiles of t_soil and m_soil (assuming they are |
---|
| 1924 | !-- horizontally homogeneous on this PE) |
---|
| 1925 | hom(nzb_soil:nzt_soil,1,90,:) = SPREAD( t_soil(nzb_soil:nzt_soil, & |
---|
| 1926 | nysg,nxlg), 2, & |
---|
| 1927 | statistic_regions+1 ) |
---|
| 1928 | hom(nzb_soil:nzt_soil,1,92,:) = SPREAD( m_soil(nzb_soil:nzt_soil, & |
---|
| 1929 | nysg,nxlg), 2, & |
---|
| 1930 | statistic_regions+1 ) |
---|
| 1931 | |
---|
| 1932 | END SUBROUTINE lsm_init |
---|
| 1933 | |
---|
| 1934 | |
---|
| 1935 | !------------------------------------------------------------------------------! |
---|
| 1936 | ! Description: |
---|
| 1937 | ! ------------ |
---|
| 1938 | !> Allocate land surface model arrays and define pointers |
---|
| 1939 | !------------------------------------------------------------------------------! |
---|
| 1940 | SUBROUTINE lsm_init_arrays |
---|
| 1941 | |
---|
| 1942 | |
---|
| 1943 | IMPLICIT NONE |
---|
| 1944 | |
---|
| 1945 | ! |
---|
| 1946 | !-- Allocate surface and soil temperature / humidity |
---|
| 1947 | #if defined( __nopointer ) |
---|
| 1948 | ALLOCATE ( m_liq_eb(nysg:nyng,nxlg:nxrg) ) |
---|
| 1949 | ALLOCATE ( m_liq_eb_p(nysg:nyng,nxlg:nxrg) ) |
---|
| 1950 | ALLOCATE ( m_soil(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) ) |
---|
| 1951 | ALLOCATE ( m_soil_p(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) ) |
---|
| 1952 | ALLOCATE ( t_surface(nysg:nyng,nxlg:nxrg) ) |
---|
| 1953 | ALLOCATE ( t_surface_p(nysg:nyng,nxlg:nxrg) ) |
---|
| 1954 | ALLOCATE ( t_soil(nzb_soil:nzt_soil+1,nysg:nyng,nxlg:nxrg) ) |
---|
| 1955 | ALLOCATE ( t_soil_p(nzb_soil:nzt_soil+1,nysg:nyng,nxlg:nxrg) ) |
---|
| 1956 | #else |
---|
| 1957 | ALLOCATE ( m_liq_eb_1(nysg:nyng,nxlg:nxrg) ) |
---|
| 1958 | ALLOCATE ( m_liq_eb_2(nysg:nyng,nxlg:nxrg) ) |
---|
| 1959 | ALLOCATE ( m_soil_1(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) ) |
---|
| 1960 | ALLOCATE ( m_soil_2(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) ) |
---|
| 1961 | ALLOCATE ( t_surface_1(nysg:nyng,nxlg:nxrg) ) |
---|
| 1962 | ALLOCATE ( t_surface_2(nysg:nyng,nxlg:nxrg) ) |
---|
| 1963 | ALLOCATE ( t_soil_1(nzb_soil:nzt_soil+1,nysg:nyng,nxlg:nxrg) ) |
---|
| 1964 | ALLOCATE ( t_soil_2(nzb_soil:nzt_soil+1,nysg:nyng,nxlg:nxrg) ) |
---|
| 1965 | #endif |
---|
| 1966 | |
---|
| 1967 | ! |
---|
| 1968 | !-- Allocate intermediate timestep arrays |
---|
| 1969 | ALLOCATE ( tm_liq_eb_m(nysg:nyng,nxlg:nxrg) ) |
---|
| 1970 | ALLOCATE ( tm_soil_m(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) ) |
---|
| 1971 | ALLOCATE ( tt_surface_m(nysg:nyng,nxlg:nxrg) ) |
---|
| 1972 | ALLOCATE ( tt_soil_m(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) ) |
---|
| 1973 | |
---|
| 1974 | ! |
---|
| 1975 | !-- Allocate 2D vegetation model arrays |
---|
| 1976 | ALLOCATE ( alpha_vg(nysg:nyng,nxlg:nxrg) ) |
---|
| 1977 | ALLOCATE ( building_surface(nysg:nyng,nxlg:nxrg) ) |
---|
| 1978 | ALLOCATE ( c_liq(nysg:nyng,nxlg:nxrg) ) |
---|
| 1979 | ALLOCATE ( c_veg(nysg:nyng,nxlg:nxrg) ) |
---|
| 1980 | ALLOCATE ( f_sw_in(nysg:nyng,nxlg:nxrg) ) |
---|
| 1981 | ALLOCATE ( ghf_eb(nysg:nyng,nxlg:nxrg) ) |
---|
| 1982 | ALLOCATE ( gamma_w_sat(nysg:nyng,nxlg:nxrg) ) |
---|
| 1983 | ALLOCATE ( g_d(nysg:nyng,nxlg:nxrg) ) |
---|
| 1984 | ALLOCATE ( lai(nysg:nyng,nxlg:nxrg) ) |
---|
| 1985 | ALLOCATE ( l_vg(nysg:nyng,nxlg:nxrg) ) |
---|
| 1986 | ALLOCATE ( lambda_surface_u(nysg:nyng,nxlg:nxrg) ) |
---|
| 1987 | ALLOCATE ( lambda_surface_s(nysg:nyng,nxlg:nxrg) ) |
---|
| 1988 | ALLOCATE ( m_fc(nysg:nyng,nxlg:nxrg) ) |
---|
| 1989 | ALLOCATE ( m_res(nysg:nyng,nxlg:nxrg) ) |
---|
| 1990 | ALLOCATE ( m_sat(nysg:nyng,nxlg:nxrg) ) |
---|
| 1991 | ALLOCATE ( m_wilt(nysg:nyng,nxlg:nxrg) ) |
---|
| 1992 | ALLOCATE ( n_vg(nysg:nyng,nxlg:nxrg) ) |
---|
| 1993 | ALLOCATE ( pave_surface(nysg:nyng,nxlg:nxrg) ) |
---|
| 1994 | ALLOCATE ( qsws_eb(nysg:nyng,nxlg:nxrg) ) |
---|
| 1995 | ALLOCATE ( qsws_soil_eb(nysg:nyng,nxlg:nxrg) ) |
---|
| 1996 | ALLOCATE ( qsws_liq_eb(nysg:nyng,nxlg:nxrg) ) |
---|
| 1997 | ALLOCATE ( qsws_veg_eb(nysg:nyng,nxlg:nxrg) ) |
---|
| 1998 | ALLOCATE ( rad_net_l(nysg:nyng,nxlg:nxrg) ) |
---|
| 1999 | ALLOCATE ( r_a(nysg:nyng,nxlg:nxrg) ) |
---|
| 2000 | ALLOCATE ( r_canopy(nysg:nyng,nxlg:nxrg) ) |
---|
| 2001 | ALLOCATE ( r_soil(nysg:nyng,nxlg:nxrg) ) |
---|
| 2002 | ALLOCATE ( r_soil_min(nysg:nyng,nxlg:nxrg) ) |
---|
| 2003 | ALLOCATE ( r_s(nysg:nyng,nxlg:nxrg) ) |
---|
| 2004 | ALLOCATE ( r_canopy_min(nysg:nyng,nxlg:nxrg) ) |
---|
| 2005 | ALLOCATE ( shf_eb(nysg:nyng,nxlg:nxrg) ) |
---|
| 2006 | ALLOCATE ( soil_type_2d(nysg:nyng,nxlg:nxrg) ) |
---|
| 2007 | ALLOCATE ( veg_type_2d(nysg:nyng,nxlg:nxrg) ) |
---|
| 2008 | ALLOCATE ( water_surface(nysg:nyng,nxlg:nxrg) ) |
---|
| 2009 | |
---|
[1978] | 2010 | water_surface = .FALSE. |
---|
| 2011 | pave_surface = .FALSE. |
---|
| 2012 | |
---|
[1817] | 2013 | #if ! defined( __nopointer ) |
---|
| 2014 | ! |
---|
| 2015 | !-- Initial assignment of the pointers |
---|
| 2016 | t_soil => t_soil_1; t_soil_p => t_soil_2 |
---|
| 2017 | t_surface => t_surface_1; t_surface_p => t_surface_2 |
---|
| 2018 | m_soil => m_soil_1; m_soil_p => m_soil_2 |
---|
| 2019 | m_liq_eb => m_liq_eb_1; m_liq_eb_p => m_liq_eb_2 |
---|
| 2020 | #endif |
---|
| 2021 | |
---|
| 2022 | |
---|
| 2023 | END SUBROUTINE lsm_init_arrays |
---|
| 2024 | |
---|
| 2025 | |
---|
| 2026 | !------------------------------------------------------------------------------! |
---|
| 2027 | ! Description: |
---|
| 2028 | ! ------------ |
---|
| 2029 | !> Parin for &lsmpar for land surface model |
---|
| 2030 | !------------------------------------------------------------------------------! |
---|
| 2031 | SUBROUTINE lsm_parin |
---|
| 2032 | |
---|
| 2033 | |
---|
| 2034 | IMPLICIT NONE |
---|
| 2035 | |
---|
| 2036 | CHARACTER (LEN=80) :: line !< dummy string that contains the current line of the parameter file |
---|
| 2037 | |
---|
| 2038 | NAMELIST /lsm_par/ alpha_vangenuchten, c_surface, & |
---|
| 2039 | canopy_resistance_coefficient, & |
---|
| 2040 | conserve_water_content, & |
---|
| 2041 | f_shortwave_incoming, field_capacity, & |
---|
| 2042 | hydraulic_conductivity, & |
---|
| 2043 | lambda_surface_stable, & |
---|
| 2044 | lambda_surface_unstable, leaf_area_index, & |
---|
| 2045 | l_vangenuchten, min_canopy_resistance, & |
---|
| 2046 | min_soil_resistance, n_vangenuchten, & |
---|
| 2047 | pave_depth, pave_heat_capacity, & |
---|
| 2048 | pave_heat_conductivity, & |
---|
| 2049 | residual_moisture, root_fraction, & |
---|
| 2050 | saturation_moisture, skip_time_do_lsm, & |
---|
| 2051 | soil_moisture, soil_temperature, soil_type, & |
---|
| 2052 | vegetation_coverage, veg_type, wilting_point,& |
---|
| 2053 | zs, z0_eb, z0h_eb, z0q_eb |
---|
| 2054 | |
---|
| 2055 | line = ' ' |
---|
| 2056 | |
---|
| 2057 | ! |
---|
| 2058 | !-- Try to find land surface model package |
---|
| 2059 | REWIND ( 11 ) |
---|
| 2060 | line = ' ' |
---|
| 2061 | DO WHILE ( INDEX( line, '&lsm_par' ) == 0 ) |
---|
| 2062 | READ ( 11, '(A)', END=10 ) line |
---|
| 2063 | ENDDO |
---|
| 2064 | BACKSPACE ( 11 ) |
---|
| 2065 | |
---|
| 2066 | ! |
---|
| 2067 | !-- Read user-defined namelist |
---|
| 2068 | READ ( 11, lsm_par ) |
---|
| 2069 | |
---|
| 2070 | ! |
---|
| 2071 | !-- Set flag that indicates that the land surface model is switched on |
---|
| 2072 | land_surface = .TRUE. |
---|
| 2073 | |
---|
| 2074 | 10 CONTINUE |
---|
| 2075 | |
---|
| 2076 | |
---|
| 2077 | END SUBROUTINE lsm_parin |
---|
| 2078 | |
---|
| 2079 | |
---|
| 2080 | !------------------------------------------------------------------------------! |
---|
| 2081 | ! Description: |
---|
| 2082 | ! ------------ |
---|
[1682] | 2083 | !> Soil model as part of the land surface model. The model predicts soil |
---|
| 2084 | !> temperature and water content. |
---|
[1496] | 2085 | !------------------------------------------------------------------------------! |
---|
| 2086 | SUBROUTINE lsm_soil_model |
---|
| 2087 | |
---|
| 2088 | |
---|
| 2089 | IMPLICIT NONE |
---|
| 2090 | |
---|
[1682] | 2091 | INTEGER(iwp) :: i !< running index |
---|
| 2092 | INTEGER(iwp) :: j !< running index |
---|
| 2093 | INTEGER(iwp) :: k !< running index |
---|
[1496] | 2094 | |
---|
[1682] | 2095 | REAL(wp) :: h_vg !< Van Genuchten coef. h |
---|
[1496] | 2096 | |
---|
[1682] | 2097 | REAL(wp), DIMENSION(nzb_soil:nzt_soil) :: gamma_temp, & !< temp. gamma |
---|
[1691] | 2098 | lambda_temp, & !< temp. lambda |
---|
| 2099 | tend !< tendency |
---|
[1496] | 2100 | |
---|
[1691] | 2101 | DO i = nxlg, nxrg |
---|
| 2102 | DO j = nysg, nyng |
---|
[1788] | 2103 | |
---|
| 2104 | IF ( pave_surface(j,i) ) THEN |
---|
| 2105 | rho_c_total(nzb_soil,j,i) = pave_heat_capacity |
---|
| 2106 | lambda_temp(nzb_soil) = pave_heat_conductivity |
---|
| 2107 | ENDIF |
---|
| 2108 | |
---|
| 2109 | IF ( .NOT. water_surface(j,i) ) THEN |
---|
| 2110 | DO k = nzb_soil, nzt_soil |
---|
| 2111 | |
---|
| 2112 | |
---|
| 2113 | IF ( pave_surface(j,i) .AND. zs(k) <= pave_depth ) THEN |
---|
| 2114 | |
---|
| 2115 | rho_c_total(k,j,i) = pave_heat_capacity |
---|
| 2116 | lambda_temp(k) = pave_heat_conductivity |
---|
| 2117 | |
---|
| 2118 | ELSE |
---|
[1496] | 2119 | ! |
---|
[1788] | 2120 | !-- Calculate volumetric heat capacity of the soil, taking |
---|
| 2121 | !-- into account water content |
---|
| 2122 | rho_c_total(k,j,i) = (rho_c_soil * (1.0_wp - m_sat(j,i)) & |
---|
| 2123 | + rho_c_water * m_soil(k,j,i)) |
---|
[1496] | 2124 | |
---|
| 2125 | ! |
---|
[1788] | 2126 | !-- Calculate soil heat conductivity at the center of the soil |
---|
| 2127 | !-- layers |
---|
| 2128 | lambda_h_sat = lambda_h_sm ** (1.0_wp - m_sat(j,i)) * & |
---|
| 2129 | lambda_h_water ** m_soil(k,j,i) |
---|
[1691] | 2130 | |
---|
[1817] | 2131 | ke = 1.0_wp + LOG10(MAX(0.1_wp,m_soil(k,j,i) & |
---|
| 2132 | / m_sat(j,i))) |
---|
[1691] | 2133 | |
---|
[1788] | 2134 | lambda_temp(k) = ke * (lambda_h_sat - lambda_h_dry) + & |
---|
| 2135 | lambda_h_dry |
---|
| 2136 | ENDIF |
---|
[1496] | 2137 | |
---|
[1788] | 2138 | ENDDO |
---|
[1496] | 2139 | |
---|
| 2140 | ! |
---|
[1788] | 2141 | !-- Calculate soil heat conductivity (lambda_h) at the _stag level |
---|
| 2142 | !-- using linear interpolation. For pavement surface, the |
---|
| 2143 | !-- true pavement depth is considered |
---|
| 2144 | DO k = nzb_soil, nzt_soil-1 |
---|
| 2145 | IF ( pave_surface(j,i) .AND. zs(k) < pave_depth & |
---|
| 2146 | .AND. zs(k+1) > pave_depth ) THEN |
---|
| 2147 | lambda_h(k,j,i) = ( pave_depth - zs(k) ) / dz_soil(k+1) & |
---|
| 2148 | * lambda_temp(k) & |
---|
| 2149 | + ( 1.0_wp - ( pave_depth - zs(k) ) & |
---|
| 2150 | / dz_soil(k+1) ) * lambda_temp(k+1) |
---|
| 2151 | ELSE |
---|
| 2152 | lambda_h(k,j,i) = ( lambda_temp(k+1) + lambda_temp(k) ) & |
---|
| 2153 | * 0.5_wp |
---|
| 2154 | ENDIF |
---|
| 2155 | ENDDO |
---|
| 2156 | lambda_h(nzt_soil,j,i) = lambda_temp(nzt_soil) |
---|
[1496] | 2157 | |
---|
[1788] | 2158 | |
---|
| 2159 | |
---|
| 2160 | |
---|
[1496] | 2161 | ! |
---|
[1788] | 2162 | !-- Prognostic equation for soil temperature t_soil |
---|
| 2163 | tend(:) = 0.0_wp |
---|
[1496] | 2164 | |
---|
[1788] | 2165 | tend(nzb_soil) = (1.0_wp/rho_c_total(nzb_soil,j,i)) * & |
---|
| 2166 | ( lambda_h(nzb_soil,j,i) * ( t_soil(nzb_soil+1,j,i) & |
---|
| 2167 | - t_soil(nzb_soil,j,i) ) * ddz_soil(nzb_soil+1) & |
---|
| 2168 | + ghf_eb(j,i) ) * ddz_soil_stag(nzb_soil) |
---|
[1691] | 2169 | |
---|
[1788] | 2170 | DO k = nzb_soil+1, nzt_soil |
---|
| 2171 | tend(k) = (1.0_wp/rho_c_total(k,j,i)) & |
---|
| 2172 | * ( lambda_h(k,j,i) & |
---|
| 2173 | * ( t_soil(k+1,j,i) - t_soil(k,j,i) ) & |
---|
| 2174 | * ddz_soil(k+1) & |
---|
| 2175 | - lambda_h(k-1,j,i) & |
---|
| 2176 | * ( t_soil(k,j,i) - t_soil(k-1,j,i) ) & |
---|
| 2177 | * ddz_soil(k) & |
---|
| 2178 | ) * ddz_soil_stag(k) |
---|
[1691] | 2179 | |
---|
[1788] | 2180 | ENDDO |
---|
[1691] | 2181 | |
---|
[1788] | 2182 | t_soil_p(nzb_soil:nzt_soil,j,i) = t_soil(nzb_soil:nzt_soil,j,i)& |
---|
| 2183 | + dt_3d * ( tsc(2) & |
---|
| 2184 | * tend(nzb_soil:nzt_soil) & |
---|
| 2185 | + tsc(3) & |
---|
| 2186 | * tt_soil_m(:,j,i) ) |
---|
[1496] | 2187 | |
---|
| 2188 | ! |
---|
[1788] | 2189 | !-- Calculate t_soil tendencies for the next Runge-Kutta step |
---|
| 2190 | IF ( timestep_scheme(1:5) == 'runge' ) THEN |
---|
| 2191 | IF ( intermediate_timestep_count == 1 ) THEN |
---|
| 2192 | DO k = nzb_soil, nzt_soil |
---|
| 2193 | tt_soil_m(k,j,i) = tend(k) |
---|
| 2194 | ENDDO |
---|
| 2195 | ELSEIF ( intermediate_timestep_count < & |
---|
| 2196 | intermediate_timestep_count_max ) THEN |
---|
| 2197 | DO k = nzb_soil, nzt_soil |
---|
| 2198 | tt_soil_m(k,j,i) = -9.5625_wp * tend(k) + 5.3125_wp & |
---|
[1551] | 2199 | * tt_soil_m(k,j,i) |
---|
[1788] | 2200 | ENDDO |
---|
| 2201 | ENDIF |
---|
[1496] | 2202 | ENDIF |
---|
| 2203 | |
---|
| 2204 | |
---|
[1788] | 2205 | DO k = nzb_soil, nzt_soil |
---|
[1551] | 2206 | |
---|
[1496] | 2207 | ! |
---|
[1788] | 2208 | !-- Calculate soil diffusivity at the center of the soil layers |
---|
| 2209 | lambda_temp(k) = (- b_ch * gamma_w_sat(j,i) * psi_sat & |
---|
| 2210 | / m_sat(j,i) ) * ( MAX( m_soil(k,j,i), & |
---|
| 2211 | m_wilt(j,i) ) / m_sat(j,i) )**( & |
---|
| 2212 | b_ch + 2.0_wp ) |
---|
[1496] | 2213 | |
---|
| 2214 | ! |
---|
[1788] | 2215 | !-- Parametrization of Van Genuchten |
---|
| 2216 | IF ( soil_type /= 7 ) THEN |
---|
[1551] | 2217 | ! |
---|
[1788] | 2218 | !-- Calculate the hydraulic conductivity after Van Genuchten |
---|
| 2219 | !-- (1980) |
---|
| 2220 | h_vg = ( ( (m_res(j,i) - m_sat(j,i)) / ( m_res(j,i) - & |
---|
| 2221 | MAX( m_soil(k,j,i), m_wilt(j,i) ) ) )**( & |
---|
| 2222 | n_vg(j,i) / (n_vg(j,i) - 1.0_wp ) ) - 1.0_wp & |
---|
| 2223 | )**( 1.0_wp / n_vg(j,i) ) / alpha_vg(j,i) |
---|
[1496] | 2224 | |
---|
[1691] | 2225 | |
---|
[1788] | 2226 | gamma_temp(k) = gamma_w_sat(j,i) * ( ( (1.0_wp + & |
---|
| 2227 | ( alpha_vg(j,i) * h_vg )**n_vg(j,i))**( & |
---|
| 2228 | 1.0_wp - 1.0_wp / n_vg(j,i) ) - ( & |
---|
| 2229 | alpha_vg(j,i) * h_vg )**( n_vg(j,i) & |
---|
| 2230 | - 1.0_wp) )**2 ) & |
---|
| 2231 | / ( ( 1.0_wp + ( alpha_vg(j,i) * h_vg & |
---|
| 2232 | )**n_vg(j,i) )**( ( 1.0_wp - 1.0_wp & |
---|
| 2233 | / n_vg(j,i) ) *( l_vg(j,i) + 2.0_wp) ) ) |
---|
[1496] | 2234 | |
---|
[1551] | 2235 | ! |
---|
[1788] | 2236 | !-- Parametrization of Clapp & Hornberger |
---|
| 2237 | ELSE |
---|
| 2238 | gamma_temp(k) = gamma_w_sat(j,i) * ( m_soil(k,j,i) & |
---|
| 2239 | / m_sat(j,i) )**(2.0_wp * b_ch + 3.0_wp) |
---|
| 2240 | ENDIF |
---|
[1551] | 2241 | |
---|
[1788] | 2242 | ENDDO |
---|
[1496] | 2243 | |
---|
| 2244 | ! |
---|
[1788] | 2245 | !-- Prognostic equation for soil moisture content. Only performed, |
---|
| 2246 | !-- when humidity is enabled in the atmosphere and the surface type |
---|
| 2247 | !-- is not pavement (implies dry soil below). |
---|
| 2248 | IF ( humidity .AND. .NOT. pave_surface(j,i) ) THEN |
---|
| 2249 | ! |
---|
| 2250 | !-- Calculate soil diffusivity (lambda_w) at the _stag level |
---|
| 2251 | !-- using linear interpolation. To do: replace this with |
---|
| 2252 | !-- ECMWF-IFS Eq. 8.81 |
---|
| 2253 | DO k = nzb_soil, nzt_soil-1 |
---|
[1496] | 2254 | |
---|
[1788] | 2255 | lambda_w(k,j,i) = ( lambda_temp(k+1) + lambda_temp(k) ) & |
---|
| 2256 | * 0.5_wp |
---|
| 2257 | gamma_w(k,j,i) = ( gamma_temp(k+1) + gamma_temp(k) ) & |
---|
| 2258 | * 0.5_wp |
---|
[1496] | 2259 | |
---|
[1788] | 2260 | ENDDO |
---|
[1496] | 2261 | |
---|
| 2262 | ! |
---|
| 2263 | ! |
---|
[1817] | 2264 | !-- In case of a closed bottom (= water content is conserved), |
---|
| 2265 | !-- set hydraulic conductivity to zero to that no water will be |
---|
| 2266 | !-- lost in the bottom layer. |
---|
[1788] | 2267 | IF ( conserve_water_content ) THEN |
---|
| 2268 | gamma_w(nzt_soil,j,i) = 0.0_wp |
---|
| 2269 | ELSE |
---|
| 2270 | gamma_w(nzt_soil,j,i) = gamma_temp(nzt_soil) |
---|
| 2271 | ENDIF |
---|
[1496] | 2272 | |
---|
[1817] | 2273 | !-- The root extraction (= root_extr * qsws_veg_eb / (rho_l |
---|
| 2274 | !-- * l_v)) ensures the mass conservation for water. The |
---|
| 2275 | !-- transpiration of plants equals the cumulative withdrawals by |
---|
| 2276 | !-- the roots in the soil. The scheme takes into account the |
---|
| 2277 | !-- availability of water in the soil layers as well as the root |
---|
| 2278 | !-- fraction in the respective layer. Layer with moisture below |
---|
| 2279 | !-- wilting point will not contribute, which reflects the |
---|
| 2280 | !-- preference of plants to take water from moister layers. |
---|
[1496] | 2281 | |
---|
| 2282 | ! |
---|
[1817] | 2283 | !-- Calculate the root extraction (ECMWF 7.69, the sum of |
---|
| 2284 | !-- root_extr = 1). The energy balance solver guarantees a |
---|
| 2285 | !-- positive transpiration, so that there is no need for an |
---|
| 2286 | !-- additional check. |
---|
[1966] | 2287 | m_total = 0.0_wp |
---|
[1691] | 2288 | DO k = nzb_soil, nzt_soil |
---|
[1788] | 2289 | IF ( m_soil(k,j,i) > m_wilt(j,i) ) THEN |
---|
| 2290 | m_total = m_total + root_fr(k,j,i) * m_soil(k,j,i) |
---|
| 2291 | ENDIF |
---|
| 2292 | ENDDO |
---|
[1496] | 2293 | |
---|
[1788] | 2294 | IF ( m_total > 0.0_wp ) THEN |
---|
| 2295 | DO k = nzb_soil, nzt_soil |
---|
| 2296 | IF ( m_soil(k,j,i) > m_wilt(j,i) ) THEN |
---|
[1817] | 2297 | root_extr(k) = root_fr(k,j,i) * m_soil(k,j,i) & |
---|
| 2298 | / m_total |
---|
[1788] | 2299 | ELSE |
---|
| 2300 | root_extr(k) = 0.0_wp |
---|
| 2301 | ENDIF |
---|
| 2302 | ENDDO |
---|
| 2303 | ENDIF |
---|
| 2304 | |
---|
[1496] | 2305 | ! |
---|
[1788] | 2306 | !-- Prognostic equation for soil water content m_soil. |
---|
| 2307 | tend(:) = 0.0_wp |
---|
| 2308 | |
---|
| 2309 | tend(nzb_soil) = ( lambda_w(nzb_soil,j,i) * ( & |
---|
[1551] | 2310 | m_soil(nzb_soil+1,j,i) - m_soil(nzb_soil,j,i) ) & |
---|
[1691] | 2311 | * ddz_soil(nzb_soil+1) - gamma_w(nzb_soil,j,i) - ( & |
---|
[1788] | 2312 | root_extr(nzb_soil) * qsws_veg_eb(j,i) & |
---|
| 2313 | + qsws_soil_eb(j,i) ) * drho_l_lv ) & |
---|
| 2314 | * ddz_soil_stag(nzb_soil) |
---|
[1496] | 2315 | |
---|
[1788] | 2316 | DO k = nzb_soil+1, nzt_soil-1 |
---|
| 2317 | tend(k) = ( lambda_w(k,j,i) * ( m_soil(k+1,j,i) & |
---|
| 2318 | - m_soil(k,j,i) ) * ddz_soil(k+1) & |
---|
| 2319 | - gamma_w(k,j,i) & |
---|
| 2320 | - lambda_w(k-1,j,i) * (m_soil(k,j,i) - & |
---|
| 2321 | m_soil(k-1,j,i)) * ddz_soil(k) & |
---|
| 2322 | + gamma_w(k-1,j,i) - (root_extr(k) & |
---|
| 2323 | * qsws_veg_eb(j,i) * drho_l_lv) & |
---|
| 2324 | ) * ddz_soil_stag(k) |
---|
[1496] | 2325 | |
---|
[1788] | 2326 | ENDDO |
---|
| 2327 | tend(nzt_soil) = ( - gamma_w(nzt_soil,j,i) & |
---|
| 2328 | - lambda_w(nzt_soil-1,j,i) & |
---|
| 2329 | * (m_soil(nzt_soil,j,i) & |
---|
| 2330 | - m_soil(nzt_soil-1,j,i)) & |
---|
| 2331 | * ddz_soil(nzt_soil) & |
---|
| 2332 | + gamma_w(nzt_soil-1,j,i) - ( & |
---|
| 2333 | root_extr(nzt_soil) & |
---|
| 2334 | * qsws_veg_eb(j,i) * drho_l_lv ) & |
---|
| 2335 | ) * ddz_soil_stag(nzt_soil) |
---|
[1496] | 2336 | |
---|
[1788] | 2337 | m_soil_p(nzb_soil:nzt_soil,j,i) = m_soil(nzb_soil:nzt_soil,j,i)& |
---|
| 2338 | + dt_3d * ( tsc(2) * tend(:) & |
---|
| 2339 | + tsc(3) * tm_soil_m(:,j,i) ) |
---|
| 2340 | |
---|
[1496] | 2341 | ! |
---|
[1788] | 2342 | !-- Account for dry soils (find a better solution here!) |
---|
| 2343 | DO k = nzb_soil, nzt_soil |
---|
| 2344 | IF ( m_soil_p(k,j,i) < 0.0_wp ) m_soil_p(k,j,i) = 0.0_wp |
---|
| 2345 | ENDDO |
---|
[1496] | 2346 | |
---|
| 2347 | ! |
---|
[1788] | 2348 | !-- Calculate m_soil tendencies for the next Runge-Kutta step |
---|
| 2349 | IF ( timestep_scheme(1:5) == 'runge' ) THEN |
---|
| 2350 | IF ( intermediate_timestep_count == 1 ) THEN |
---|
| 2351 | DO k = nzb_soil, nzt_soil |
---|
| 2352 | tm_soil_m(k,j,i) = tend(k) |
---|
| 2353 | ENDDO |
---|
| 2354 | ELSEIF ( intermediate_timestep_count < & |
---|
| 2355 | intermediate_timestep_count_max ) THEN |
---|
| 2356 | DO k = nzb_soil, nzt_soil |
---|
| 2357 | tm_soil_m(k,j,i) = -9.5625_wp * tend(k) + 5.3125_wp& |
---|
[1496] | 2358 | * tm_soil_m(k,j,i) |
---|
[1788] | 2359 | ENDDO |
---|
| 2360 | ENDIF |
---|
[1496] | 2361 | ENDIF |
---|
[1788] | 2362 | |
---|
[1496] | 2363 | ENDIF |
---|
| 2364 | |
---|
| 2365 | ENDIF |
---|
| 2366 | |
---|
| 2367 | ENDDO |
---|
| 2368 | ENDDO |
---|
| 2369 | |
---|
[1788] | 2370 | END SUBROUTINE lsm_soil_model |
---|
| 2371 | |
---|
[1817] | 2372 | |
---|
| 2373 | !------------------------------------------------------------------------------! |
---|
| 2374 | ! Description: |
---|
| 2375 | ! ------------ |
---|
| 2376 | !> Swapping of timelevels |
---|
| 2377 | !------------------------------------------------------------------------------! |
---|
| 2378 | SUBROUTINE lsm_swap_timelevel ( mod_count ) |
---|
[1788] | 2379 | |
---|
[1817] | 2380 | IMPLICIT NONE |
---|
| 2381 | |
---|
| 2382 | INTEGER, INTENT(IN) :: mod_count |
---|
| 2383 | |
---|
| 2384 | #if defined( __nopointer ) |
---|
| 2385 | |
---|
| 2386 | t_surface = t_surface_p |
---|
| 2387 | t_soil = t_soil_p |
---|
| 2388 | IF ( humidity ) THEN |
---|
| 2389 | m_soil = m_soil_p |
---|
| 2390 | m_liq_eb = m_liq_eb_p |
---|
| 2391 | ENDIF |
---|
| 2392 | |
---|
| 2393 | #else |
---|
| 2394 | |
---|
| 2395 | SELECT CASE ( mod_count ) |
---|
| 2396 | |
---|
| 2397 | CASE ( 0 ) |
---|
| 2398 | |
---|
| 2399 | t_surface => t_surface_1; t_surface_p => t_surface_2 |
---|
| 2400 | t_soil => t_soil_1; t_soil_p => t_soil_2 |
---|
| 2401 | IF ( humidity ) THEN |
---|
| 2402 | m_soil => m_soil_1; m_soil_p => m_soil_2 |
---|
| 2403 | m_liq_eb => m_liq_eb_1; m_liq_eb_p => m_liq_eb_2 |
---|
| 2404 | ENDIF |
---|
| 2405 | |
---|
| 2406 | |
---|
| 2407 | CASE ( 1 ) |
---|
| 2408 | |
---|
| 2409 | t_surface => t_surface_2; t_surface_p => t_surface_1 |
---|
| 2410 | t_soil => t_soil_2; t_soil_p => t_soil_1 |
---|
| 2411 | IF ( humidity ) THEN |
---|
| 2412 | m_soil => m_soil_2; m_soil_p => m_soil_1 |
---|
| 2413 | m_liq_eb => m_liq_eb_2; m_liq_eb_p => m_liq_eb_1 |
---|
| 2414 | ENDIF |
---|
| 2415 | |
---|
| 2416 | END SELECT |
---|
| 2417 | #endif |
---|
| 2418 | |
---|
| 2419 | END SUBROUTINE lsm_swap_timelevel |
---|
| 2420 | |
---|
| 2421 | |
---|
[1972] | 2422 | |
---|
| 2423 | |
---|
[1788] | 2424 | !------------------------------------------------------------------------------! |
---|
[1972] | 2425 | ! |
---|
[1788] | 2426 | ! Description: |
---|
| 2427 | ! ------------ |
---|
[1976] | 2428 | !> Subroutine for averaging 3D data |
---|
[1972] | 2429 | !------------------------------------------------------------------------------! |
---|
| 2430 | SUBROUTINE lsm_3d_data_averaging( mode, variable ) |
---|
| 2431 | |
---|
| 2432 | |
---|
| 2433 | USE control_parameters |
---|
| 2434 | |
---|
| 2435 | USE indices |
---|
| 2436 | |
---|
| 2437 | USE kinds |
---|
| 2438 | |
---|
| 2439 | IMPLICIT NONE |
---|
| 2440 | |
---|
| 2441 | CHARACTER (LEN=*) :: mode !< |
---|
| 2442 | CHARACTER (LEN=*) :: variable !< |
---|
| 2443 | |
---|
| 2444 | INTEGER(iwp) :: i !< |
---|
| 2445 | INTEGER(iwp) :: j !< |
---|
| 2446 | INTEGER(iwp) :: k !< |
---|
| 2447 | |
---|
| 2448 | IF ( mode == 'allocate' ) THEN |
---|
| 2449 | |
---|
| 2450 | SELECT CASE ( TRIM( variable ) ) |
---|
| 2451 | |
---|
| 2452 | CASE ( 'c_liq*' ) |
---|
| 2453 | IF ( .NOT. ALLOCATED( c_liq_av ) ) THEN |
---|
| 2454 | ALLOCATE( c_liq_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 2455 | ENDIF |
---|
| 2456 | c_liq_av = 0.0_wp |
---|
| 2457 | |
---|
| 2458 | CASE ( 'c_soil*' ) |
---|
| 2459 | IF ( .NOT. ALLOCATED( c_soil_av ) ) THEN |
---|
| 2460 | ALLOCATE( c_soil_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 2461 | ENDIF |
---|
| 2462 | c_soil_av = 0.0_wp |
---|
| 2463 | |
---|
| 2464 | CASE ( 'c_veg*' ) |
---|
| 2465 | IF ( .NOT. ALLOCATED( c_veg_av ) ) THEN |
---|
| 2466 | ALLOCATE( c_veg_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 2467 | ENDIF |
---|
| 2468 | c_veg_av = 0.0_wp |
---|
| 2469 | |
---|
| 2470 | CASE ( 'ghf_eb*' ) |
---|
| 2471 | IF ( .NOT. ALLOCATED( ghf_eb_av ) ) THEN |
---|
| 2472 | ALLOCATE( ghf_eb_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 2473 | ENDIF |
---|
| 2474 | ghf_eb_av = 0.0_wp |
---|
| 2475 | |
---|
| 2476 | CASE ( 'lai*' ) |
---|
| 2477 | IF ( .NOT. ALLOCATED( lai_av ) ) THEN |
---|
| 2478 | ALLOCATE( lai_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 2479 | ENDIF |
---|
| 2480 | lai_av = 0.0_wp |
---|
| 2481 | |
---|
| 2482 | CASE ( 'm_liq_eb*' ) |
---|
| 2483 | IF ( .NOT. ALLOCATED( m_liq_eb_av ) ) THEN |
---|
| 2484 | ALLOCATE( m_liq_eb_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 2485 | ENDIF |
---|
| 2486 | m_liq_eb_av = 0.0_wp |
---|
| 2487 | |
---|
| 2488 | CASE ( 'm_soil' ) |
---|
| 2489 | IF ( .NOT. ALLOCATED( m_soil_av ) ) THEN |
---|
| 2490 | ALLOCATE( m_soil_av(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) ) |
---|
| 2491 | ENDIF |
---|
| 2492 | m_soil_av = 0.0_wp |
---|
| 2493 | |
---|
| 2494 | CASE ( 'qsws_eb*' ) |
---|
| 2495 | IF ( .NOT. ALLOCATED( qsws_eb_av ) ) THEN |
---|
| 2496 | ALLOCATE( qsws_eb_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 2497 | ENDIF |
---|
| 2498 | qsws_eb_av = 0.0_wp |
---|
| 2499 | |
---|
| 2500 | CASE ( 'qsws_liq_eb*' ) |
---|
| 2501 | IF ( .NOT. ALLOCATED( qsws_liq_eb_av ) ) THEN |
---|
| 2502 | ALLOCATE( qsws_liq_eb_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 2503 | ENDIF |
---|
| 2504 | qsws_liq_eb_av = 0.0_wp |
---|
| 2505 | |
---|
| 2506 | CASE ( 'qsws_soil_eb*' ) |
---|
| 2507 | IF ( .NOT. ALLOCATED( qsws_soil_eb_av ) ) THEN |
---|
| 2508 | ALLOCATE( qsws_soil_eb_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 2509 | ENDIF |
---|
| 2510 | qsws_soil_eb_av = 0.0_wp |
---|
| 2511 | |
---|
| 2512 | CASE ( 'qsws_veg_eb*' ) |
---|
| 2513 | IF ( .NOT. ALLOCATED( qsws_veg_eb_av ) ) THEN |
---|
| 2514 | ALLOCATE( qsws_veg_eb_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 2515 | ENDIF |
---|
| 2516 | qsws_veg_eb_av = 0.0_wp |
---|
| 2517 | |
---|
| 2518 | CASE ( 'r_a*' ) |
---|
| 2519 | IF ( .NOT. ALLOCATED( r_a_av ) ) THEN |
---|
| 2520 | ALLOCATE( r_a_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 2521 | ENDIF |
---|
| 2522 | r_a_av = 0.0_wp |
---|
| 2523 | |
---|
| 2524 | CASE ( 'r_s*' ) |
---|
| 2525 | IF ( .NOT. ALLOCATED( r_s_av ) ) THEN |
---|
| 2526 | ALLOCATE( r_s_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 2527 | ENDIF |
---|
| 2528 | r_s_av = 0.0_wp |
---|
| 2529 | |
---|
| 2530 | CASE ( 'shf_eb*' ) |
---|
| 2531 | IF ( .NOT. ALLOCATED( shf_eb_av ) ) THEN |
---|
| 2532 | ALLOCATE( shf_eb_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 2533 | ENDIF |
---|
| 2534 | shf_eb_av = 0.0_wp |
---|
| 2535 | |
---|
| 2536 | CASE ( 't_soil' ) |
---|
| 2537 | IF ( .NOT. ALLOCATED( t_soil_av ) ) THEN |
---|
| 2538 | ALLOCATE( t_soil_av(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) ) |
---|
| 2539 | ENDIF |
---|
| 2540 | t_soil_av = 0.0_wp |
---|
| 2541 | |
---|
| 2542 | CASE DEFAULT |
---|
| 2543 | CONTINUE |
---|
| 2544 | |
---|
| 2545 | END SELECT |
---|
| 2546 | |
---|
| 2547 | ELSEIF ( mode == 'sum' ) THEN |
---|
| 2548 | |
---|
| 2549 | SELECT CASE ( TRIM( variable ) ) |
---|
| 2550 | |
---|
| 2551 | CASE ( 'c_liq*' ) |
---|
| 2552 | DO i = nxlg, nxrg |
---|
| 2553 | DO j = nysg, nyng |
---|
| 2554 | c_liq_av(j,i) = c_liq_av(j,i) + c_liq(j,i) |
---|
| 2555 | ENDDO |
---|
| 2556 | ENDDO |
---|
| 2557 | |
---|
| 2558 | CASE ( 'c_soil*' ) |
---|
| 2559 | DO i = nxlg, nxrg |
---|
| 2560 | DO j = nysg, nyng |
---|
| 2561 | c_soil_av(j,i) = c_soil_av(j,i) + (1.0 - c_veg(j,i)) |
---|
| 2562 | ENDDO |
---|
| 2563 | ENDDO |
---|
| 2564 | |
---|
| 2565 | CASE ( 'c_veg*' ) |
---|
| 2566 | DO i = nxlg, nxrg |
---|
| 2567 | DO j = nysg, nyng |
---|
| 2568 | c_veg_av(j,i) = c_veg_av(j,i) + c_veg(j,i) |
---|
| 2569 | ENDDO |
---|
| 2570 | ENDDO |
---|
| 2571 | |
---|
| 2572 | CASE ( 'ghf_eb*' ) |
---|
| 2573 | DO i = nxlg, nxrg |
---|
| 2574 | DO j = nysg, nyng |
---|
| 2575 | ghf_eb_av(j,i) = ghf_eb_av(j,i) + ghf_eb(j,i) |
---|
| 2576 | ENDDO |
---|
| 2577 | ENDDO |
---|
| 2578 | |
---|
| 2579 | CASE ( 'lai*' ) |
---|
| 2580 | DO i = nxlg, nxrg |
---|
| 2581 | DO j = nysg, nyng |
---|
| 2582 | lai_av(j,i) = lai_av(j,i) + lai(j,i) |
---|
| 2583 | ENDDO |
---|
| 2584 | ENDDO |
---|
| 2585 | |
---|
| 2586 | CASE ( 'm_liq_eb*' ) |
---|
| 2587 | DO i = nxlg, nxrg |
---|
| 2588 | DO j = nysg, nyng |
---|
| 2589 | m_liq_eb_av(j,i) = m_liq_eb_av(j,i) + m_liq_eb(j,i) |
---|
| 2590 | ENDDO |
---|
| 2591 | ENDDO |
---|
| 2592 | |
---|
| 2593 | CASE ( 'm_soil' ) |
---|
| 2594 | DO i = nxlg, nxrg |
---|
| 2595 | DO j = nysg, nyng |
---|
| 2596 | DO k = nzb_soil, nzt_soil |
---|
| 2597 | m_soil_av(k,j,i) = m_soil_av(k,j,i) + m_soil(k,j,i) |
---|
| 2598 | ENDDO |
---|
| 2599 | ENDDO |
---|
| 2600 | ENDDO |
---|
| 2601 | |
---|
| 2602 | CASE ( 'qsws_eb*' ) |
---|
| 2603 | DO i = nxlg, nxrg |
---|
| 2604 | DO j = nysg, nyng |
---|
| 2605 | qsws_eb_av(j,i) = qsws_eb_av(j,i) + qsws_eb(j,i) |
---|
| 2606 | ENDDO |
---|
| 2607 | ENDDO |
---|
| 2608 | |
---|
| 2609 | CASE ( 'qsws_liq_eb*' ) |
---|
| 2610 | DO i = nxlg, nxrg |
---|
| 2611 | DO j = nysg, nyng |
---|
| 2612 | qsws_liq_eb_av(j,i) = qsws_liq_eb_av(j,i) + qsws_liq_eb(j,i) |
---|
| 2613 | ENDDO |
---|
| 2614 | ENDDO |
---|
| 2615 | |
---|
| 2616 | CASE ( 'qsws_soil_eb*' ) |
---|
| 2617 | DO i = nxlg, nxrg |
---|
| 2618 | DO j = nysg, nyng |
---|
| 2619 | qsws_soil_eb_av(j,i) = qsws_soil_eb_av(j,i) + qsws_soil_eb(j,i) |
---|
| 2620 | ENDDO |
---|
| 2621 | ENDDO |
---|
| 2622 | |
---|
| 2623 | CASE ( 'qsws_veg_eb*' ) |
---|
| 2624 | DO i = nxlg, nxrg |
---|
| 2625 | DO j = nysg, nyng |
---|
| 2626 | qsws_veg_eb_av(j,i) = qsws_veg_eb_av(j,i) + qsws_veg_eb(j,i) |
---|
| 2627 | ENDDO |
---|
| 2628 | ENDDO |
---|
| 2629 | |
---|
| 2630 | CASE ( 'r_a*' ) |
---|
| 2631 | DO i = nxlg, nxrg |
---|
| 2632 | DO j = nysg, nyng |
---|
| 2633 | r_a_av(j,i) = r_a_av(j,i) + r_a(j,i) |
---|
| 2634 | ENDDO |
---|
| 2635 | ENDDO |
---|
| 2636 | |
---|
| 2637 | CASE ( 'r_s*' ) |
---|
| 2638 | DO i = nxlg, nxrg |
---|
| 2639 | DO j = nysg, nyng |
---|
| 2640 | r_s_av(j,i) = r_s_av(j,i) + r_s(j,i) |
---|
| 2641 | ENDDO |
---|
| 2642 | ENDDO |
---|
| 2643 | |
---|
| 2644 | CASE ( 'shf_eb*' ) |
---|
| 2645 | DO i = nxlg, nxrg |
---|
| 2646 | DO j = nysg, nyng |
---|
| 2647 | shf_eb_av(j,i) = shf_eb_av(j,i) + shf_eb(j,i) |
---|
| 2648 | ENDDO |
---|
| 2649 | ENDDO |
---|
| 2650 | |
---|
| 2651 | CASE ( 't_soil' ) |
---|
| 2652 | DO i = nxlg, nxrg |
---|
| 2653 | DO j = nysg, nyng |
---|
| 2654 | DO k = nzb_soil, nzt_soil |
---|
| 2655 | t_soil_av(k,j,i) = t_soil_av(k,j,i) + t_soil(k,j,i) |
---|
| 2656 | ENDDO |
---|
| 2657 | ENDDO |
---|
| 2658 | ENDDO |
---|
| 2659 | |
---|
| 2660 | CASE DEFAULT |
---|
| 2661 | CONTINUE |
---|
| 2662 | |
---|
| 2663 | END SELECT |
---|
| 2664 | |
---|
| 2665 | ELSEIF ( mode == 'average' ) THEN |
---|
| 2666 | |
---|
| 2667 | SELECT CASE ( TRIM( variable ) ) |
---|
| 2668 | |
---|
| 2669 | CASE ( 'c_liq*' ) |
---|
| 2670 | DO i = nxlg, nxrg |
---|
| 2671 | DO j = nysg, nyng |
---|
| 2672 | c_liq_av(j,i) = c_liq_av(j,i) / REAL( average_count_3d, KIND=wp ) |
---|
| 2673 | ENDDO |
---|
| 2674 | ENDDO |
---|
| 2675 | |
---|
| 2676 | CASE ( 'c_soil*' ) |
---|
| 2677 | DO i = nxlg, nxrg |
---|
| 2678 | DO j = nysg, nyng |
---|
| 2679 | c_soil_av(j,i) = c_soil_av(j,i) / REAL( average_count_3d, KIND=wp ) |
---|
| 2680 | ENDDO |
---|
| 2681 | ENDDO |
---|
| 2682 | |
---|
| 2683 | CASE ( 'c_veg*' ) |
---|
| 2684 | DO i = nxlg, nxrg |
---|
| 2685 | DO j = nysg, nyng |
---|
| 2686 | c_veg_av(j,i) = c_veg_av(j,i) / REAL( average_count_3d, KIND=wp ) |
---|
| 2687 | ENDDO |
---|
| 2688 | ENDDO |
---|
| 2689 | |
---|
| 2690 | CASE ( 'ghf_eb*' ) |
---|
| 2691 | DO i = nxlg, nxrg |
---|
| 2692 | DO j = nysg, nyng |
---|
| 2693 | ghf_eb_av(j,i) = ghf_eb_av(j,i) / REAL( average_count_3d, KIND=wp ) |
---|
| 2694 | ENDDO |
---|
| 2695 | ENDDO |
---|
| 2696 | |
---|
| 2697 | CASE ( 'lai*' ) |
---|
| 2698 | DO i = nxlg, nxrg |
---|
| 2699 | DO j = nysg, nyng |
---|
| 2700 | lai_av(j,i) = lai_av(j,i) / REAL( average_count_3d, KIND=wp ) |
---|
| 2701 | ENDDO |
---|
| 2702 | ENDDO |
---|
| 2703 | |
---|
| 2704 | CASE ( 'm_liq_eb*' ) |
---|
| 2705 | DO i = nxlg, nxrg |
---|
| 2706 | DO j = nysg, nyng |
---|
| 2707 | m_liq_eb_av(j,i) = m_liq_eb_av(j,i) / REAL( average_count_3d, KIND=wp ) |
---|
| 2708 | ENDDO |
---|
| 2709 | ENDDO |
---|
| 2710 | |
---|
| 2711 | CASE ( 'm_soil' ) |
---|
| 2712 | DO i = nxlg, nxrg |
---|
| 2713 | DO j = nysg, nyng |
---|
| 2714 | DO k = nzb_soil, nzt_soil |
---|
| 2715 | m_soil_av(k,j,i) = m_soil_av(k,j,i) / REAL( average_count_3d, KIND=wp ) |
---|
| 2716 | ENDDO |
---|
| 2717 | ENDDO |
---|
| 2718 | ENDDO |
---|
| 2719 | |
---|
| 2720 | CASE ( 'qsws_eb*' ) |
---|
| 2721 | DO i = nxlg, nxrg |
---|
| 2722 | DO j = nysg, nyng |
---|
| 2723 | qsws_eb_av(j,i) = qsws_eb_av(j,i) / REAL( average_count_3d, KIND=wp ) |
---|
| 2724 | ENDDO |
---|
| 2725 | ENDDO |
---|
| 2726 | |
---|
| 2727 | CASE ( 'qsws_liq_eb*' ) |
---|
| 2728 | DO i = nxlg, nxrg |
---|
| 2729 | DO j = nysg, nyng |
---|
| 2730 | qsws_liq_eb_av(j,i) = qsws_liq_eb_av(j,i) / REAL( average_count_3d, KIND=wp ) |
---|
| 2731 | ENDDO |
---|
| 2732 | ENDDO |
---|
| 2733 | |
---|
| 2734 | CASE ( 'qsws_soil_eb*' ) |
---|
| 2735 | DO i = nxlg, nxrg |
---|
| 2736 | DO j = nysg, nyng |
---|
| 2737 | qsws_soil_eb_av(j,i) = qsws_soil_eb_av(j,i) / REAL( average_count_3d, KIND=wp ) |
---|
| 2738 | ENDDO |
---|
| 2739 | ENDDO |
---|
| 2740 | |
---|
| 2741 | CASE ( 'qsws_veg_eb*' ) |
---|
| 2742 | DO i = nxlg, nxrg |
---|
| 2743 | DO j = nysg, nyng |
---|
| 2744 | qsws_veg_eb_av(j,i) = qsws_veg_eb_av(j,i) / REAL( average_count_3d, KIND=wp ) |
---|
| 2745 | ENDDO |
---|
| 2746 | ENDDO |
---|
| 2747 | |
---|
| 2748 | CASE ( 'r_a*' ) |
---|
| 2749 | DO i = nxlg, nxrg |
---|
| 2750 | DO j = nysg, nyng |
---|
| 2751 | r_a_av(j,i) = r_a_av(j,i) / REAL( average_count_3d, KIND=wp ) |
---|
| 2752 | ENDDO |
---|
| 2753 | ENDDO |
---|
| 2754 | |
---|
| 2755 | CASE ( 'r_s*' ) |
---|
| 2756 | DO i = nxlg, nxrg |
---|
| 2757 | DO j = nysg, nyng |
---|
| 2758 | r_s_av(j,i) = r_s_av(j,i) / REAL( average_count_3d, KIND=wp ) |
---|
| 2759 | ENDDO |
---|
| 2760 | ENDDO |
---|
| 2761 | |
---|
| 2762 | CASE ( 't_soil' ) |
---|
| 2763 | DO i = nxlg, nxrg |
---|
| 2764 | DO j = nysg, nyng |
---|
| 2765 | DO k = nzb_soil, nzt_soil |
---|
| 2766 | t_soil_av(k,j,i) = t_soil_av(k,j,i) / REAL( average_count_3d, KIND=wp ) |
---|
| 2767 | ENDDO |
---|
| 2768 | ENDDO |
---|
| 2769 | ENDDO |
---|
| 2770 | |
---|
| 2771 | END SELECT |
---|
| 2772 | |
---|
| 2773 | ENDIF |
---|
| 2774 | |
---|
| 2775 | END SUBROUTINE lsm_3d_data_averaging |
---|
| 2776 | |
---|
| 2777 | |
---|
| 2778 | !------------------------------------------------------------------------------! |
---|
| 2779 | ! |
---|
| 2780 | ! Description: |
---|
| 2781 | ! ------------ |
---|
[1976] | 2782 | !> Subroutine defining appropriate grid for netcdf variables. |
---|
[1972] | 2783 | !> It is called out from subroutine netcdf. |
---|
| 2784 | !------------------------------------------------------------------------------! |
---|
[1976] | 2785 | SUBROUTINE lsm_define_netcdf_grid( var, found, grid_x, grid_y, grid_z ) |
---|
[1972] | 2786 | |
---|
[1976] | 2787 | IMPLICIT NONE |
---|
[1972] | 2788 | |
---|
[1976] | 2789 | CHARACTER (LEN=*), INTENT(IN) :: var !< |
---|
| 2790 | LOGICAL, INTENT(OUT) :: found !< |
---|
| 2791 | CHARACTER (LEN=*), INTENT(OUT) :: grid_x !< |
---|
| 2792 | CHARACTER (LEN=*), INTENT(OUT) :: grid_y !< |
---|
| 2793 | CHARACTER (LEN=*), INTENT(OUT) :: grid_z !< |
---|
[1972] | 2794 | |
---|
[1976] | 2795 | found = .TRUE. |
---|
[1972] | 2796 | |
---|
| 2797 | ! |
---|
[1976] | 2798 | !-- Check for the grid |
---|
| 2799 | SELECT CASE ( TRIM( var ) ) |
---|
[1972] | 2800 | |
---|
[1976] | 2801 | CASE ( 'm_soil', 't_soil', 'm_soil_xy', 't_soil_xy', 'm_soil_xz', & |
---|
| 2802 | 't_soil_xz', 'm_soil_yz', 't_soil_yz' ) |
---|
| 2803 | grid_x = 'x' |
---|
| 2804 | grid_y = 'y' |
---|
| 2805 | grid_z = 'zs' |
---|
[1972] | 2806 | |
---|
[1976] | 2807 | CASE DEFAULT |
---|
| 2808 | found = .FALSE. |
---|
| 2809 | grid_x = 'none' |
---|
| 2810 | grid_y = 'none' |
---|
| 2811 | grid_z = 'none' |
---|
| 2812 | END SELECT |
---|
[1972] | 2813 | |
---|
[1976] | 2814 | END SUBROUTINE lsm_define_netcdf_grid |
---|
[1972] | 2815 | |
---|
| 2816 | !------------------------------------------------------------------------------! |
---|
| 2817 | ! |
---|
| 2818 | ! Description: |
---|
| 2819 | ! ------------ |
---|
[1976] | 2820 | !> Subroutine defining 3D output variables |
---|
[1972] | 2821 | !------------------------------------------------------------------------------! |
---|
| 2822 | SUBROUTINE lsm_data_output_2d( av, variable, found, grid, mode, local_pf, & |
---|
| 2823 | two_d, nzb_do, nzt_do ) |
---|
| 2824 | |
---|
| 2825 | USE indices |
---|
| 2826 | |
---|
| 2827 | USE kinds |
---|
| 2828 | |
---|
| 2829 | |
---|
| 2830 | IMPLICIT NONE |
---|
| 2831 | |
---|
| 2832 | CHARACTER (LEN=*) :: grid !< |
---|
| 2833 | CHARACTER (LEN=*) :: mode !< |
---|
| 2834 | CHARACTER (LEN=*) :: variable !< |
---|
| 2835 | |
---|
| 2836 | INTEGER(iwp) :: av !< |
---|
| 2837 | INTEGER(iwp) :: i !< |
---|
| 2838 | INTEGER(iwp) :: j !< |
---|
| 2839 | INTEGER(iwp) :: k !< |
---|
| 2840 | INTEGER(iwp) :: nzb_do !< |
---|
| 2841 | INTEGER(iwp) :: nzt_do !< |
---|
| 2842 | |
---|
| 2843 | LOGICAL :: found !< |
---|
| 2844 | LOGICAL :: two_d !< flag parameter that indicates 2D variables (horizontal cross sections) |
---|
| 2845 | |
---|
| 2846 | REAL(wp), DIMENSION(nxlg:nxrg,nysg:nyng,nzb:nzt+1) :: local_pf !< |
---|
| 2847 | |
---|
| 2848 | found = .TRUE. |
---|
| 2849 | |
---|
| 2850 | SELECT CASE ( TRIM( variable ) ) |
---|
| 2851 | |
---|
| 2852 | |
---|
| 2853 | CASE ( 'c_liq*_xy' ) ! 2d-array |
---|
| 2854 | IF ( av == 0 ) THEN |
---|
| 2855 | DO i = nxlg, nxrg |
---|
| 2856 | DO j = nysg, nyng |
---|
| 2857 | local_pf(i,j,nzb+1) = c_liq(j,i) * c_veg(j,i) |
---|
| 2858 | ENDDO |
---|
| 2859 | ENDDO |
---|
| 2860 | ELSE |
---|
| 2861 | DO i = nxlg, nxrg |
---|
| 2862 | DO j = nysg, nyng |
---|
| 2863 | local_pf(i,j,nzb+1) = c_liq_av(j,i) |
---|
| 2864 | ENDDO |
---|
| 2865 | ENDDO |
---|
| 2866 | ENDIF |
---|
| 2867 | |
---|
| 2868 | two_d = .TRUE. |
---|
| 2869 | grid = 'zu1' |
---|
| 2870 | |
---|
| 2871 | CASE ( 'c_soil*_xy' ) ! 2d-array |
---|
| 2872 | IF ( av == 0 ) THEN |
---|
| 2873 | DO i = nxlg, nxrg |
---|
| 2874 | DO j = nysg, nyng |
---|
| 2875 | local_pf(i,j,nzb+1) = 1.0_wp - c_veg(j,i) |
---|
| 2876 | ENDDO |
---|
| 2877 | ENDDO |
---|
| 2878 | ELSE |
---|
| 2879 | DO i = nxlg, nxrg |
---|
| 2880 | DO j = nysg, nyng |
---|
| 2881 | local_pf(i,j,nzb+1) = c_soil_av(j,i) |
---|
| 2882 | ENDDO |
---|
| 2883 | ENDDO |
---|
| 2884 | ENDIF |
---|
| 2885 | |
---|
| 2886 | two_d = .TRUE. |
---|
| 2887 | grid = 'zu1' |
---|
| 2888 | |
---|
| 2889 | CASE ( 'c_veg*_xy' ) ! 2d-array |
---|
| 2890 | IF ( av == 0 ) THEN |
---|
| 2891 | DO i = nxlg, nxrg |
---|
| 2892 | DO j = nysg, nyng |
---|
| 2893 | local_pf(i,j,nzb+1) = c_veg(j,i) |
---|
| 2894 | ENDDO |
---|
| 2895 | ENDDO |
---|
| 2896 | ELSE |
---|
| 2897 | DO i = nxlg, nxrg |
---|
| 2898 | DO j = nysg, nyng |
---|
| 2899 | local_pf(i,j,nzb+1) = c_veg_av(j,i) |
---|
| 2900 | ENDDO |
---|
| 2901 | ENDDO |
---|
| 2902 | ENDIF |
---|
| 2903 | |
---|
| 2904 | two_d = .TRUE. |
---|
| 2905 | grid = 'zu1' |
---|
| 2906 | |
---|
| 2907 | CASE ( 'ghf_eb*_xy' ) ! 2d-array |
---|
| 2908 | IF ( av == 0 ) THEN |
---|
| 2909 | DO i = nxlg, nxrg |
---|
| 2910 | DO j = nysg, nyng |
---|
| 2911 | local_pf(i,j,nzb+1) = ghf_eb(j,i) |
---|
| 2912 | ENDDO |
---|
| 2913 | ENDDO |
---|
| 2914 | ELSE |
---|
| 2915 | DO i = nxlg, nxrg |
---|
| 2916 | DO j = nysg, nyng |
---|
| 2917 | local_pf(i,j,nzb+1) = ghf_eb_av(j,i) |
---|
| 2918 | ENDDO |
---|
| 2919 | ENDDO |
---|
| 2920 | ENDIF |
---|
| 2921 | |
---|
| 2922 | two_d = .TRUE. |
---|
| 2923 | grid = 'zu1' |
---|
| 2924 | |
---|
| 2925 | CASE ( 'lai*_xy' ) ! 2d-array |
---|
| 2926 | IF ( av == 0 ) THEN |
---|
| 2927 | DO i = nxlg, nxrg |
---|
| 2928 | DO j = nysg, nyng |
---|
| 2929 | local_pf(i,j,nzb+1) = lai(j,i) |
---|
| 2930 | ENDDO |
---|
| 2931 | ENDDO |
---|
| 2932 | ELSE |
---|
| 2933 | DO i = nxlg, nxrg |
---|
| 2934 | DO j = nysg, nyng |
---|
| 2935 | local_pf(i,j,nzb+1) = lai_av(j,i) |
---|
| 2936 | ENDDO |
---|
| 2937 | ENDDO |
---|
| 2938 | ENDIF |
---|
| 2939 | |
---|
| 2940 | two_d = .TRUE. |
---|
| 2941 | grid = 'zu1' |
---|
| 2942 | |
---|
| 2943 | CASE ( 'm_liq_eb*_xy' ) ! 2d-array |
---|
| 2944 | IF ( av == 0 ) THEN |
---|
| 2945 | DO i = nxlg, nxrg |
---|
| 2946 | DO j = nysg, nyng |
---|
| 2947 | local_pf(i,j,nzb+1) = m_liq_eb(j,i) |
---|
| 2948 | ENDDO |
---|
| 2949 | ENDDO |
---|
| 2950 | ELSE |
---|
| 2951 | DO i = nxlg, nxrg |
---|
| 2952 | DO j = nysg, nyng |
---|
| 2953 | local_pf(i,j,nzb+1) = m_liq_eb_av(j,i) |
---|
| 2954 | ENDDO |
---|
| 2955 | ENDDO |
---|
| 2956 | ENDIF |
---|
| 2957 | |
---|
| 2958 | two_d = .TRUE. |
---|
| 2959 | grid = 'zu1' |
---|
| 2960 | |
---|
| 2961 | CASE ( 'm_soil_xy', 'm_soil_xz', 'm_soil_yz' ) |
---|
| 2962 | IF ( av == 0 ) THEN |
---|
| 2963 | DO i = nxlg, nxrg |
---|
| 2964 | DO j = nysg, nyng |
---|
| 2965 | DO k = nzb_soil, nzt_soil |
---|
| 2966 | local_pf(i,j,k) = m_soil(k,j,i) |
---|
| 2967 | ENDDO |
---|
| 2968 | ENDDO |
---|
| 2969 | ENDDO |
---|
| 2970 | ELSE |
---|
| 2971 | DO i = nxlg, nxrg |
---|
| 2972 | DO j = nysg, nyng |
---|
| 2973 | DO k = nzb_soil, nzt_soil |
---|
| 2974 | local_pf(i,j,k) = m_soil_av(k,j,i) |
---|
| 2975 | ENDDO |
---|
| 2976 | ENDDO |
---|
| 2977 | ENDDO |
---|
| 2978 | ENDIF |
---|
| 2979 | |
---|
| 2980 | nzb_do = nzb_soil |
---|
| 2981 | nzt_do = nzt_soil |
---|
| 2982 | |
---|
| 2983 | IF ( mode == 'xy' ) grid = 'zs' |
---|
| 2984 | |
---|
| 2985 | CASE ( 'qsws_eb*_xy' ) ! 2d-array |
---|
| 2986 | IF ( av == 0 ) THEN |
---|
| 2987 | DO i = nxlg, nxrg |
---|
| 2988 | DO j = nysg, nyng |
---|
| 2989 | local_pf(i,j,nzb+1) = qsws_eb(j,i) |
---|
| 2990 | ENDDO |
---|
| 2991 | ENDDO |
---|
| 2992 | ELSE |
---|
| 2993 | DO i = nxlg, nxrg |
---|
| 2994 | DO j = nysg, nyng |
---|
| 2995 | local_pf(i,j,nzb+1) = qsws_eb_av(j,i) |
---|
| 2996 | ENDDO |
---|
| 2997 | ENDDO |
---|
| 2998 | ENDIF |
---|
| 2999 | |
---|
| 3000 | two_d = .TRUE. |
---|
| 3001 | grid = 'zu1' |
---|
| 3002 | |
---|
| 3003 | CASE ( 'qsws_liq_eb*_xy' ) ! 2d-array |
---|
| 3004 | IF ( av == 0 ) THEN |
---|
| 3005 | DO i = nxlg, nxrg |
---|
| 3006 | DO j = nysg, nyng |
---|
| 3007 | local_pf(i,j,nzb+1) = qsws_liq_eb(j,i) |
---|
| 3008 | ENDDO |
---|
| 3009 | ENDDO |
---|
| 3010 | ELSE |
---|
| 3011 | DO i = nxlg, nxrg |
---|
| 3012 | DO j = nysg, nyng |
---|
| 3013 | local_pf(i,j,nzb+1) = qsws_liq_eb_av(j,i) |
---|
| 3014 | ENDDO |
---|
| 3015 | ENDDO |
---|
| 3016 | ENDIF |
---|
| 3017 | |
---|
| 3018 | two_d = .TRUE. |
---|
| 3019 | grid = 'zu1' |
---|
| 3020 | |
---|
| 3021 | CASE ( 'qsws_soil_eb*_xy' ) ! 2d-array |
---|
| 3022 | IF ( av == 0 ) THEN |
---|
| 3023 | DO i = nxlg, nxrg |
---|
| 3024 | DO j = nysg, nyng |
---|
| 3025 | local_pf(i,j,nzb+1) = qsws_soil_eb(j,i) |
---|
| 3026 | ENDDO |
---|
| 3027 | ENDDO |
---|
| 3028 | ELSE |
---|
| 3029 | DO i = nxlg, nxrg |
---|
| 3030 | DO j = nysg, nyng |
---|
| 3031 | local_pf(i,j,nzb+1) = qsws_soil_eb_av(j,i) |
---|
| 3032 | ENDDO |
---|
| 3033 | ENDDO |
---|
| 3034 | ENDIF |
---|
| 3035 | |
---|
| 3036 | two_d = .TRUE. |
---|
| 3037 | grid = 'zu1' |
---|
| 3038 | |
---|
| 3039 | CASE ( 'qsws_veg_eb*_xy' ) ! 2d-array |
---|
| 3040 | IF ( av == 0 ) THEN |
---|
| 3041 | DO i = nxlg, nxrg |
---|
| 3042 | DO j = nysg, nyng |
---|
| 3043 | local_pf(i,j,nzb+1) = qsws_veg_eb(j,i) |
---|
| 3044 | ENDDO |
---|
| 3045 | ENDDO |
---|
| 3046 | ELSE |
---|
| 3047 | DO i = nxlg, nxrg |
---|
| 3048 | DO j = nysg, nyng |
---|
| 3049 | local_pf(i,j,nzb+1) = qsws_veg_eb_av(j,i) |
---|
| 3050 | ENDDO |
---|
| 3051 | ENDDO |
---|
| 3052 | ENDIF |
---|
| 3053 | |
---|
| 3054 | two_d = .TRUE. |
---|
| 3055 | grid = 'zu1' |
---|
| 3056 | |
---|
| 3057 | |
---|
| 3058 | CASE ( 'r_a*_xy' ) ! 2d-array |
---|
| 3059 | IF ( av == 0 ) THEN |
---|
| 3060 | DO i = nxlg, nxrg |
---|
| 3061 | DO j = nysg, nyng |
---|
| 3062 | local_pf(i,j,nzb+1) = r_a(j,i) |
---|
| 3063 | ENDDO |
---|
| 3064 | ENDDO |
---|
| 3065 | ELSE |
---|
| 3066 | DO i = nxlg, nxrg |
---|
| 3067 | DO j = nysg, nyng |
---|
| 3068 | local_pf(i,j,nzb+1) = r_a_av(j,i) |
---|
| 3069 | ENDDO |
---|
| 3070 | ENDDO |
---|
| 3071 | ENDIF |
---|
| 3072 | |
---|
| 3073 | two_d = .TRUE. |
---|
| 3074 | grid = 'zu1' |
---|
| 3075 | |
---|
| 3076 | CASE ( 'r_s*_xy' ) ! 2d-array |
---|
| 3077 | IF ( av == 0 ) THEN |
---|
| 3078 | DO i = nxlg, nxrg |
---|
| 3079 | DO j = nysg, nyng |
---|
| 3080 | local_pf(i,j,nzb+1) = r_s(j,i) |
---|
| 3081 | ENDDO |
---|
| 3082 | ENDDO |
---|
| 3083 | ELSE |
---|
| 3084 | DO i = nxlg, nxrg |
---|
| 3085 | DO j = nysg, nyng |
---|
| 3086 | local_pf(i,j,nzb+1) = r_s_av(j,i) |
---|
| 3087 | ENDDO |
---|
| 3088 | ENDDO |
---|
| 3089 | ENDIF |
---|
| 3090 | |
---|
| 3091 | two_d = .TRUE. |
---|
| 3092 | grid = 'zu1' |
---|
| 3093 | |
---|
| 3094 | CASE ( 'shf_eb*_xy' ) ! 2d-array |
---|
| 3095 | IF ( av == 0 ) THEN |
---|
| 3096 | DO i = nxlg, nxrg |
---|
| 3097 | DO j = nysg, nyng |
---|
| 3098 | local_pf(i,j,nzb+1) = shf_eb(j,i) |
---|
| 3099 | ENDDO |
---|
| 3100 | ENDDO |
---|
| 3101 | ELSE |
---|
| 3102 | DO i = nxlg, nxrg |
---|
| 3103 | DO j = nysg, nyng |
---|
| 3104 | local_pf(i,j,nzb+1) = shf_eb_av(j,i) |
---|
| 3105 | ENDDO |
---|
| 3106 | ENDDO |
---|
| 3107 | ENDIF |
---|
| 3108 | |
---|
| 3109 | two_d = .TRUE. |
---|
| 3110 | grid = 'zu1' |
---|
| 3111 | |
---|
| 3112 | CASE ( 't_soil_xy', 't_soil_xz', 't_soil_yz' ) |
---|
| 3113 | IF ( av == 0 ) THEN |
---|
| 3114 | DO i = nxlg, nxrg |
---|
| 3115 | DO j = nysg, nyng |
---|
| 3116 | DO k = nzb_soil, nzt_soil |
---|
| 3117 | local_pf(i,j,k) = t_soil(k,j,i) |
---|
| 3118 | ENDDO |
---|
| 3119 | ENDDO |
---|
| 3120 | ENDDO |
---|
| 3121 | ELSE |
---|
| 3122 | DO i = nxlg, nxrg |
---|
| 3123 | DO j = nysg, nyng |
---|
| 3124 | DO k = nzb_soil, nzt_soil |
---|
| 3125 | local_pf(i,j,k) = t_soil_av(k,j,i) |
---|
| 3126 | ENDDO |
---|
| 3127 | ENDDO |
---|
| 3128 | ENDDO |
---|
| 3129 | ENDIF |
---|
| 3130 | |
---|
| 3131 | nzb_do = nzb_soil |
---|
| 3132 | nzt_do = nzt_soil |
---|
| 3133 | |
---|
| 3134 | IF ( mode == 'xy' ) grid = 'zs' |
---|
| 3135 | |
---|
| 3136 | CASE DEFAULT |
---|
| 3137 | found = .FALSE. |
---|
| 3138 | grid = 'none' |
---|
| 3139 | |
---|
| 3140 | END SELECT |
---|
| 3141 | |
---|
| 3142 | END SUBROUTINE lsm_data_output_2d |
---|
| 3143 | |
---|
| 3144 | |
---|
| 3145 | !------------------------------------------------------------------------------! |
---|
| 3146 | ! |
---|
| 3147 | ! Description: |
---|
| 3148 | ! ------------ |
---|
[1976] | 3149 | !> Subroutine defining 3D output variables |
---|
[1972] | 3150 | !------------------------------------------------------------------------------! |
---|
| 3151 | SUBROUTINE lsm_data_output_3d( av, variable, found, local_pf ) |
---|
| 3152 | |
---|
| 3153 | |
---|
| 3154 | USE indices |
---|
| 3155 | |
---|
| 3156 | USE kinds |
---|
| 3157 | |
---|
| 3158 | |
---|
| 3159 | IMPLICIT NONE |
---|
| 3160 | |
---|
| 3161 | CHARACTER (LEN=*) :: variable !< |
---|
| 3162 | |
---|
| 3163 | INTEGER(iwp) :: av !< |
---|
| 3164 | INTEGER(iwp) :: i !< |
---|
| 3165 | INTEGER(iwp) :: j !< |
---|
| 3166 | INTEGER(iwp) :: k !< |
---|
| 3167 | |
---|
| 3168 | LOGICAL :: found !< |
---|
| 3169 | |
---|
| 3170 | REAL(sp), DIMENSION(nxlg:nxrg,nysg:nyng,nzb_soil:nzt_soil) :: local_pf !< |
---|
| 3171 | |
---|
| 3172 | |
---|
| 3173 | found = .TRUE. |
---|
| 3174 | |
---|
| 3175 | |
---|
| 3176 | SELECT CASE ( TRIM( variable ) ) |
---|
| 3177 | |
---|
| 3178 | |
---|
| 3179 | CASE ( 'm_soil' ) |
---|
| 3180 | |
---|
| 3181 | IF ( av == 0 ) THEN |
---|
| 3182 | DO i = nxlg, nxrg |
---|
| 3183 | DO j = nysg, nyng |
---|
| 3184 | DO k = nzb_soil, nzt_soil |
---|
| 3185 | local_pf(i,j,k) = m_soil(k,j,i) |
---|
| 3186 | ENDDO |
---|
| 3187 | ENDDO |
---|
| 3188 | ENDDO |
---|
| 3189 | ELSE |
---|
| 3190 | DO i = nxlg, nxrg |
---|
| 3191 | DO j = nysg, nyng |
---|
| 3192 | DO k = nzb_soil, nzt_soil |
---|
| 3193 | local_pf(i,j,k) = m_soil_av(k,j,i) |
---|
| 3194 | ENDDO |
---|
| 3195 | ENDDO |
---|
| 3196 | ENDDO |
---|
| 3197 | ENDIF |
---|
| 3198 | |
---|
| 3199 | CASE ( 't_soil' ) |
---|
| 3200 | |
---|
| 3201 | IF ( av == 0 ) THEN |
---|
| 3202 | DO i = nxlg, nxrg |
---|
| 3203 | DO j = nysg, nyng |
---|
| 3204 | DO k = nzb_soil, nzt_soil |
---|
| 3205 | local_pf(i,j,k) = t_soil(k,j,i) |
---|
| 3206 | ENDDO |
---|
| 3207 | ENDDO |
---|
| 3208 | ENDDO |
---|
| 3209 | ELSE |
---|
| 3210 | DO i = nxlg, nxrg |
---|
| 3211 | DO j = nysg, nyng |
---|
| 3212 | DO k = nzb_soil, nzt_soil |
---|
| 3213 | local_pf(i,j,k) = t_soil_av(k,j,i) |
---|
| 3214 | ENDDO |
---|
| 3215 | ENDDO |
---|
| 3216 | ENDDO |
---|
| 3217 | ENDIF |
---|
| 3218 | |
---|
| 3219 | |
---|
| 3220 | CASE DEFAULT |
---|
| 3221 | found = .FALSE. |
---|
| 3222 | |
---|
| 3223 | END SELECT |
---|
| 3224 | |
---|
| 3225 | |
---|
| 3226 | END SUBROUTINE lsm_data_output_3d |
---|
| 3227 | |
---|
| 3228 | |
---|
| 3229 | !------------------------------------------------------------------------------! |
---|
| 3230 | ! |
---|
| 3231 | ! Description: |
---|
| 3232 | ! ------------ |
---|
| 3233 | !> Write restart data for land surface model |
---|
| 3234 | !------------------------------------------------------------------------------! |
---|
| 3235 | SUBROUTINE lsm_last_actions |
---|
| 3236 | |
---|
| 3237 | |
---|
| 3238 | USE control_parameters |
---|
| 3239 | |
---|
| 3240 | USE kinds |
---|
| 3241 | |
---|
| 3242 | IMPLICIT NONE |
---|
| 3243 | |
---|
| 3244 | IF ( write_binary(1:4) == 'true' ) THEN |
---|
| 3245 | IF ( ALLOCATED( c_liq_av ) ) THEN |
---|
| 3246 | WRITE ( 14 ) 'c_liq_av '; WRITE ( 14 ) c_liq_av |
---|
| 3247 | ENDIF |
---|
| 3248 | IF ( ALLOCATED( c_soil_av ) ) THEN |
---|
| 3249 | WRITE ( 14 ) 'c_soil_av '; WRITE ( 14 ) c_soil_av |
---|
| 3250 | ENDIF |
---|
| 3251 | IF ( ALLOCATED( c_veg_av ) ) THEN |
---|
| 3252 | WRITE ( 14 ) 'c_veg_av '; WRITE ( 14 ) c_veg_av |
---|
| 3253 | ENDIF |
---|
| 3254 | IF ( ALLOCATED( ghf_eb_av ) ) THEN |
---|
| 3255 | WRITE ( 14 ) 'ghf_eb_av '; WRITE ( 14 ) ghf_eb_av |
---|
| 3256 | ENDIF |
---|
| 3257 | IF ( ALLOCATED( lai_av ) ) THEN |
---|
| 3258 | WRITE ( 14 ) 'lai_av '; WRITE ( 14 ) lai_av |
---|
| 3259 | ENDIF |
---|
| 3260 | WRITE ( 14 ) 'm_liq_eb '; WRITE ( 14 ) m_liq_eb |
---|
| 3261 | IF ( ALLOCATED( m_liq_eb_av ) ) THEN |
---|
| 3262 | WRITE ( 14 ) 'm_liq_eb_av '; WRITE ( 14 ) m_liq_eb_av |
---|
| 3263 | ENDIF |
---|
| 3264 | WRITE ( 14 ) 'm_soil '; WRITE ( 14 ) m_soil |
---|
| 3265 | IF ( ALLOCATED( m_soil_av ) ) THEN |
---|
| 3266 | WRITE ( 14 ) 'm_soil_av '; WRITE ( 14 ) m_soil_av |
---|
| 3267 | ENDIF |
---|
| 3268 | IF ( ALLOCATED( qsws_eb_av ) ) THEN |
---|
| 3269 | WRITE ( 14 ) 'qsws_eb_av '; WRITE ( 14 ) qsws_eb_av |
---|
| 3270 | ENDIF |
---|
| 3271 | IF ( ALLOCATED( qsws_liq_eb_av ) ) THEN |
---|
| 3272 | WRITE ( 14 ) 'qsws_liq_eb_av '; WRITE ( 14 ) qsws_liq_eb_av |
---|
| 3273 | ENDIF |
---|
| 3274 | IF ( ALLOCATED( qsws_soil_eb_av ) ) THEN |
---|
| 3275 | WRITE ( 14 ) 'qsws_soil_eb_av '; WRITE ( 14 ) qsws_soil_eb_av |
---|
| 3276 | ENDIF |
---|
| 3277 | IF ( ALLOCATED( qsws_veg_eb_av ) ) THEN |
---|
| 3278 | WRITE ( 14 ) 'qsws_veg_eb_av '; WRITE ( 14 ) qsws_veg_eb_av |
---|
| 3279 | ENDIF |
---|
| 3280 | IF ( ALLOCATED( shf_eb_av ) ) THEN |
---|
| 3281 | WRITE ( 14 ) 'shf_eb_av '; WRITE ( 14 ) shf_eb_av |
---|
| 3282 | ENDIF |
---|
| 3283 | WRITE ( 14 ) 't_soil '; WRITE ( 14 ) t_soil |
---|
| 3284 | IF ( ALLOCATED( t_soil_av ) ) THEN |
---|
| 3285 | WRITE ( 14 ) 't_soil_av '; WRITE ( 14 ) t_soil_av |
---|
| 3286 | ENDIF |
---|
| 3287 | |
---|
| 3288 | WRITE ( 14 ) '*** end lsm *** ' |
---|
| 3289 | |
---|
| 3290 | ENDIF |
---|
| 3291 | |
---|
| 3292 | END SUBROUTINE lsm_last_actions |
---|
| 3293 | |
---|
| 3294 | |
---|
| 3295 | SUBROUTINE lsm_read_restart_data( i, nxlfa, nxl_on_file, nxrfa, nxr_on_file, & |
---|
| 3296 | nynfa, nyn_on_file, nysfa, nys_on_file, & |
---|
| 3297 | offset_xa, offset_ya, overlap_count, & |
---|
| 3298 | tmp_2d ) |
---|
| 3299 | |
---|
| 3300 | |
---|
| 3301 | USE control_parameters |
---|
| 3302 | |
---|
| 3303 | USE indices |
---|
| 3304 | |
---|
| 3305 | USE kinds |
---|
| 3306 | |
---|
| 3307 | USE pegrid |
---|
| 3308 | |
---|
| 3309 | IMPLICIT NONE |
---|
| 3310 | |
---|
| 3311 | CHARACTER (LEN=20) :: field_char !< |
---|
| 3312 | |
---|
| 3313 | INTEGER(iwp) :: i !< |
---|
| 3314 | INTEGER(iwp) :: k !< |
---|
| 3315 | INTEGER(iwp) :: nxlc !< |
---|
| 3316 | INTEGER(iwp) :: nxlf !< |
---|
| 3317 | INTEGER(iwp) :: nxl_on_file !< |
---|
| 3318 | INTEGER(iwp) :: nxrc !< |
---|
| 3319 | INTEGER(iwp) :: nxrf !< |
---|
| 3320 | INTEGER(iwp) :: nxr_on_file !< |
---|
| 3321 | INTEGER(iwp) :: nync !< |
---|
| 3322 | INTEGER(iwp) :: nynf !< |
---|
| 3323 | INTEGER(iwp) :: nyn_on_file !< |
---|
| 3324 | INTEGER(iwp) :: nysc !< |
---|
| 3325 | INTEGER(iwp) :: nysf !< |
---|
| 3326 | INTEGER(iwp) :: nys_on_file !< |
---|
| 3327 | INTEGER(iwp) :: overlap_count !< |
---|
| 3328 | |
---|
| 3329 | INTEGER(iwp), DIMENSION(numprocs_previous_run,1000) :: nxlfa !< |
---|
| 3330 | INTEGER(iwp), DIMENSION(numprocs_previous_run,1000) :: nxrfa !< |
---|
| 3331 | INTEGER(iwp), DIMENSION(numprocs_previous_run,1000) :: nynfa !< |
---|
| 3332 | INTEGER(iwp), DIMENSION(numprocs_previous_run,1000) :: nysfa !< |
---|
| 3333 | INTEGER(iwp), DIMENSION(numprocs_previous_run,1000) :: offset_xa !< |
---|
| 3334 | INTEGER(iwp), DIMENSION(numprocs_previous_run,1000) :: offset_ya !< |
---|
| 3335 | |
---|
| 3336 | REAL(wp), & |
---|
| 3337 | DIMENSION(nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) ::& |
---|
| 3338 | tmp_2d !< |
---|
| 3339 | |
---|
| 3340 | REAL(wp), & |
---|
| 3341 | DIMENSION(nzb_soil:nzt_soil+1,nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) ::& |
---|
| 3342 | tmp_3d !< |
---|
| 3343 | |
---|
| 3344 | REAL(wp), & |
---|
| 3345 | DIMENSION(nzb_soil:nzt_soil,nys_on_file-nbgp:nyn_on_file+nbgp,nxl_on_file-nbgp:nxr_on_file+nbgp) ::& |
---|
| 3346 | tmp_3d2 !< |
---|
| 3347 | |
---|
| 3348 | |
---|
| 3349 | IF ( initializing_actions == 'read_restart_data' ) THEN |
---|
| 3350 | READ ( 13 ) field_char |
---|
| 3351 | |
---|
| 3352 | DO WHILE ( TRIM( field_char ) /= '*** end lsm ***' ) |
---|
| 3353 | |
---|
| 3354 | DO k = 1, overlap_count |
---|
| 3355 | |
---|
| 3356 | nxlf = nxlfa(i,k) |
---|
| 3357 | nxlc = nxlfa(i,k) + offset_xa(i,k) |
---|
| 3358 | nxrf = nxrfa(i,k) |
---|
| 3359 | nxrc = nxrfa(i,k) + offset_xa(i,k) |
---|
| 3360 | nysf = nysfa(i,k) |
---|
| 3361 | nysc = nysfa(i,k) + offset_ya(i,k) |
---|
| 3362 | nynf = nynfa(i,k) |
---|
| 3363 | nync = nynfa(i,k) + offset_ya(i,k) |
---|
| 3364 | |
---|
| 3365 | |
---|
| 3366 | SELECT CASE ( TRIM( field_char ) ) |
---|
| 3367 | |
---|
| 3368 | CASE ( 'c_liq_av' ) |
---|
| 3369 | IF ( .NOT. ALLOCATED( c_liq_av ) ) THEN |
---|
| 3370 | ALLOCATE( c_liq_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 3371 | ENDIF |
---|
| 3372 | IF ( k == 1 ) READ ( 13 ) tmp_2d |
---|
| 3373 | c_liq_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & |
---|
| 3374 | tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) |
---|
| 3375 | |
---|
| 3376 | CASE ( 'c_soil_av' ) |
---|
| 3377 | IF ( .NOT. ALLOCATED( c_soil_av ) ) THEN |
---|
| 3378 | ALLOCATE( c_soil_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 3379 | ENDIF |
---|
| 3380 | IF ( k == 1 ) READ ( 13 ) tmp_2d |
---|
| 3381 | c_soil_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & |
---|
| 3382 | tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) |
---|
| 3383 | |
---|
| 3384 | CASE ( 'c_veg_av' ) |
---|
| 3385 | IF ( .NOT. ALLOCATED( c_veg_av ) ) THEN |
---|
| 3386 | ALLOCATE( c_veg_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 3387 | ENDIF |
---|
| 3388 | IF ( k == 1 ) READ ( 13 ) tmp_2d |
---|
| 3389 | c_veg_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & |
---|
| 3390 | tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) |
---|
| 3391 | |
---|
| 3392 | CASE ( 'ghf_eb_av' ) |
---|
| 3393 | IF ( .NOT. ALLOCATED( ghf_eb_av ) ) THEN |
---|
| 3394 | ALLOCATE( ghf_eb_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 3395 | ENDIF |
---|
| 3396 | IF ( k == 1 ) READ ( 13 ) tmp_2d |
---|
| 3397 | ghf_eb_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & |
---|
| 3398 | tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) |
---|
| 3399 | |
---|
| 3400 | CASE ( 'm_liq_eb' ) |
---|
| 3401 | IF ( k == 1 ) READ ( 13 ) tmp_2d |
---|
| 3402 | m_liq_eb(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & |
---|
| 3403 | tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) |
---|
| 3404 | |
---|
| 3405 | CASE ( 'lai_av' ) |
---|
| 3406 | IF ( .NOT. ALLOCATED( lai_av ) ) THEN |
---|
| 3407 | ALLOCATE( lai_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 3408 | ENDIF |
---|
| 3409 | IF ( k == 1 ) READ ( 13 ) tmp_2d |
---|
| 3410 | lai_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & |
---|
| 3411 | tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) |
---|
| 3412 | |
---|
| 3413 | CASE ( 'm_liq_eb_av' ) |
---|
| 3414 | IF ( .NOT. ALLOCATED( m_liq_eb_av ) ) THEN |
---|
| 3415 | ALLOCATE( m_liq_eb_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 3416 | ENDIF |
---|
| 3417 | IF ( k == 1 ) READ ( 13 ) tmp_2d |
---|
| 3418 | m_liq_eb_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & |
---|
| 3419 | tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) |
---|
| 3420 | |
---|
| 3421 | CASE ( 'm_soil' ) |
---|
| 3422 | IF ( k == 1 ) READ ( 13 ) tmp_3d2(:,:,:) |
---|
| 3423 | m_soil(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & |
---|
| 3424 | tmp_3d2(nzb_soil:nzt_soil,nysf-nbgp:nynf & |
---|
| 3425 | +nbgp,nxlf-nbgp:nxrf+nbgp) |
---|
| 3426 | |
---|
| 3427 | CASE ( 'm_soil_av' ) |
---|
| 3428 | IF ( .NOT. ALLOCATED( m_soil_av ) ) THEN |
---|
| 3429 | ALLOCATE( m_soil_av(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) ) |
---|
| 3430 | ENDIF |
---|
| 3431 | IF ( k == 1 ) READ ( 13 ) tmp_3d2(:,:,:) |
---|
| 3432 | m_soil_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & |
---|
| 3433 | tmp_3d2(nzb_soil:nzt_soil,nysf & |
---|
| 3434 | -nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) |
---|
| 3435 | |
---|
| 3436 | CASE ( 'qsws_eb_av' ) |
---|
| 3437 | IF ( .NOT. ALLOCATED( qsws_eb_av ) ) THEN |
---|
| 3438 | ALLOCATE( qsws_eb_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 3439 | ENDIF |
---|
| 3440 | IF ( k == 1 ) READ ( 13 ) tmp_2d |
---|
| 3441 | qsws_eb_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & |
---|
| 3442 | tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) |
---|
| 3443 | |
---|
| 3444 | CASE ( 'qsws_liq_eb_av' ) |
---|
| 3445 | IF ( .NOT. ALLOCATED( qsws_liq_eb_av ) ) THEN |
---|
| 3446 | ALLOCATE( qsws_liq_eb_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 3447 | ENDIF |
---|
| 3448 | IF ( k == 1 ) READ ( 13 ) tmp_2d |
---|
| 3449 | qsws_liq_eb_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & |
---|
| 3450 | tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) |
---|
| 3451 | CASE ( 'qsws_soil_eb_av' ) |
---|
| 3452 | IF ( .NOT. ALLOCATED( qsws_soil_eb_av ) ) THEN |
---|
| 3453 | ALLOCATE( qsws_soil_eb_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 3454 | ENDIF |
---|
| 3455 | IF ( k == 1 ) READ ( 13 ) tmp_2d |
---|
| 3456 | qsws_soil_eb_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & |
---|
| 3457 | tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) |
---|
| 3458 | |
---|
| 3459 | CASE ( 'qsws_veg_eb_av' ) |
---|
| 3460 | IF ( .NOT. ALLOCATED( qsws_veg_eb_av ) ) THEN |
---|
| 3461 | ALLOCATE( qsws_veg_eb_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 3462 | ENDIF |
---|
| 3463 | IF ( k == 1 ) READ ( 13 ) tmp_2d |
---|
| 3464 | qsws_veg_eb_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & |
---|
| 3465 | tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) |
---|
| 3466 | |
---|
| 3467 | CASE ( 'shf_eb_av' ) |
---|
| 3468 | IF ( .NOT. ALLOCATED( shf_eb_av ) ) THEN |
---|
| 3469 | ALLOCATE( shf_eb_av(nysg:nyng,nxlg:nxrg) ) |
---|
| 3470 | ENDIF |
---|
| 3471 | IF ( k == 1 ) READ ( 13 ) tmp_2d |
---|
| 3472 | shf_eb_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & |
---|
| 3473 | tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) |
---|
| 3474 | |
---|
| 3475 | CASE ( 't_soil' ) |
---|
| 3476 | IF ( k == 1 ) READ ( 13 ) tmp_3d |
---|
| 3477 | t_soil(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & |
---|
| 3478 | tmp_3d(:,nysf-nbgp:nynf+nbgp, & |
---|
| 3479 | nxlf-nbgp:nxrf+nbgp) |
---|
| 3480 | |
---|
| 3481 | CASE ( 't_soil_av' ) |
---|
| 3482 | IF ( .NOT. ALLOCATED( t_soil_av ) ) THEN |
---|
| 3483 | ALLOCATE( t_soil_av(nzb_soil:nzt_soil,nysg:nyng,nxlg:nxrg) ) |
---|
| 3484 | ENDIF |
---|
| 3485 | IF ( k == 1 ) READ ( 13 ) tmp_3d2(:,:,:) |
---|
| 3486 | t_soil_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & |
---|
| 3487 | tmp_3d(:,nysf-nbgp:nynf+nbgp, & |
---|
| 3488 | nxlf-nbgp:nxrf+nbgp) |
---|
| 3489 | |
---|
| 3490 | |
---|
| 3491 | CASE DEFAULT |
---|
| 3492 | WRITE( message_string, * ) 'unknown variable named "', & |
---|
| 3493 | TRIM( field_char ), '" found in', & |
---|
| 3494 | '&data from prior run on PE ', myid |
---|
| 3495 | CALL message( 'lsm_read_restart_data', 'PA0441', 1, 2, 0, 6, & |
---|
| 3496 | 0 ) |
---|
| 3497 | |
---|
| 3498 | END SELECT |
---|
| 3499 | |
---|
| 3500 | ENDDO |
---|
| 3501 | |
---|
| 3502 | READ ( 13 ) field_char |
---|
| 3503 | |
---|
| 3504 | ENDDO |
---|
| 3505 | ENDIF |
---|
| 3506 | |
---|
| 3507 | END SUBROUTINE lsm_read_restart_data |
---|
| 3508 | |
---|
| 3509 | !------------------------------------------------------------------------------! |
---|
| 3510 | ! Description: |
---|
| 3511 | ! ------------ |
---|
[1788] | 3512 | !> Calculation of roughness length for open water (lakes, ocean). The |
---|
| 3513 | !> parameterization follows Charnock (1955). Two different implementations |
---|
| 3514 | !> are available: as in ECMWF-IFS (Beljaars 1994) or as in FLake (Subin et al. |
---|
| 3515 | !> 2012) |
---|
| 3516 | !------------------------------------------------------------------------------! |
---|
| 3517 | SUBROUTINE calc_z0_water_surface |
---|
| 3518 | |
---|
| 3519 | USE control_parameters, & |
---|
| 3520 | ONLY: g, kappa, molecular_viscosity |
---|
| 3521 | |
---|
| 3522 | IMPLICIT NONE |
---|
| 3523 | |
---|
| 3524 | INTEGER :: i !< running index |
---|
| 3525 | INTEGER :: j !< running index |
---|
| 3526 | |
---|
| 3527 | REAL(wp), PARAMETER :: alpha_ch = 0.018_wp !< Charnock constant (0.01-0.11). Use 0.01 for FLake and 0.018 for ECMWF |
---|
| 3528 | ! REAL(wp), PARAMETER :: pr_number = 0.71_wp !< molecular Prandtl number in the Charnock parameterization (differs from prandtl_number) |
---|
| 3529 | ! REAL(wp), PARAMETER :: sc_number = 0.66_wp !< molecular Schmidt number in the Charnock parameterization |
---|
| 3530 | ! REAL(wp) :: re_0 !< near-surface roughness Reynolds number |
---|
| 3531 | |
---|
| 3532 | |
---|
| 3533 | DO i = nxlg, nxrg |
---|
| 3534 | DO j = nysg, nyng |
---|
| 3535 | IF ( water_surface(j,i) ) THEN |
---|
| 3536 | |
---|
[1496] | 3537 | ! |
---|
[1788] | 3538 | !-- Disabled: FLake parameterization. Ideally, the Charnock |
---|
| 3539 | !-- coefficient should depend on the water depth and the fetch |
---|
| 3540 | !-- length |
---|
| 3541 | ! re_0 = z0(j,i) * us(j,i) / molecular_viscosity |
---|
| 3542 | ! |
---|
| 3543 | ! z0(j,i) = MAX( 0.1_wp * molecular_viscosity / us(j,i), & |
---|
| 3544 | ! alpha_ch * us(j,i) / g ) |
---|
| 3545 | ! |
---|
| 3546 | ! z0h(j,i) = z0(j,i) * EXP( - kappa / pr_number * ( 4.0_wp * SQRT( re_0 ) - 3.2_wp ) ) |
---|
| 3547 | ! z0q(j,i) = z0(j,i) * EXP( - kappa / pr_number * ( 4.0_wp * SQRT( re_0 ) - 4.2_wp ) ) |
---|
[1496] | 3548 | |
---|
[1788] | 3549 | ! |
---|
| 3550 | !-- Set minimum roughness length for u* > 0.2 |
---|
| 3551 | ! IF ( us(j,i) > 0.2_wp ) THEN |
---|
| 3552 | ! z0h(j,i) = MAX( 1.0E-5_wp, z0h(j,i) ) |
---|
| 3553 | ! z0q(j,i) = MAX( 1.0E-5_wp, z0q(j,i) ) |
---|
| 3554 | ! ENDIF |
---|
[1496] | 3555 | |
---|
[1788] | 3556 | ! |
---|
| 3557 | !-- ECMWF IFS model parameterization after Beljaars (1994). At low |
---|
| 3558 | !-- wind speed, the sea surface becomes aerodynamically smooth and |
---|
| 3559 | !-- the roughness scales with the viscosity. At high wind speed, the |
---|
| 3560 | !-- Charnock relation is used. |
---|
| 3561 | z0(j,i) = ( 0.11_wp * molecular_viscosity / us(j,i) ) & |
---|
| 3562 | + ( alpha_ch * us(j,i)**2 / g ) |
---|
[1496] | 3563 | |
---|
[1788] | 3564 | z0h(j,i) = 0.40_wp * molecular_viscosity / us(j,i) |
---|
| 3565 | z0q(j,i) = 0.62_wp * molecular_viscosity / us(j,i) |
---|
[1496] | 3566 | |
---|
[1788] | 3567 | ENDIF |
---|
| 3568 | ENDDO |
---|
| 3569 | ENDDO |
---|
| 3570 | |
---|
| 3571 | END SUBROUTINE calc_z0_water_surface |
---|
| 3572 | |
---|
| 3573 | |
---|
[1496] | 3574 | !------------------------------------------------------------------------------! |
---|
| 3575 | ! Description: |
---|
| 3576 | ! ------------ |
---|
[1788] | 3577 | !> Calculation of specific humidity of the skin layer (surface). It is assumend |
---|
| 3578 | !> that the skin is always saturated. |
---|
[1496] | 3579 | !------------------------------------------------------------------------------! |
---|
[1551] | 3580 | SUBROUTINE calc_q_surface |
---|
[1496] | 3581 | |
---|
| 3582 | IMPLICIT NONE |
---|
| 3583 | |
---|
[1682] | 3584 | INTEGER :: i !< running index |
---|
| 3585 | INTEGER :: j !< running index |
---|
| 3586 | INTEGER :: k !< running index |
---|
[1788] | 3587 | |
---|
[1682] | 3588 | REAL(wp) :: resistance !< aerodynamic and soil resistance term |
---|
[1496] | 3589 | |
---|
[1691] | 3590 | DO i = nxlg, nxrg |
---|
| 3591 | DO j = nysg, nyng |
---|
[1496] | 3592 | k = nzb_s_inner(j,i) |
---|
| 3593 | |
---|
| 3594 | ! |
---|
| 3595 | !-- Calculate water vapour pressure at saturation |
---|
[1691] | 3596 | e_s = 0.01_wp * 610.78_wp * EXP( 17.269_wp * ( t_surface_p(j,i) & |
---|
| 3597 | - 273.16_wp ) / ( t_surface_p(j,i) - 35.86_wp ) ) |
---|
[1496] | 3598 | |
---|
| 3599 | ! |
---|
| 3600 | !-- Calculate specific humidity at saturation |
---|
| 3601 | q_s = 0.622_wp * e_s / surface_pressure |
---|
| 3602 | |
---|
| 3603 | resistance = r_a(j,i) / (r_a(j,i) + r_s(j,i)) |
---|
| 3604 | |
---|
| 3605 | ! |
---|
| 3606 | !-- Calculate specific humidity at surface |
---|
[1691] | 3607 | IF ( cloud_physics ) THEN |
---|
| 3608 | q(k,j,i) = resistance * q_s + (1.0_wp - resistance) & |
---|
| 3609 | * ( q(k+1,j,i) - ql(k+1,j,i) ) |
---|
| 3610 | ELSE |
---|
| 3611 | q(k,j,i) = resistance * q_s + (1.0_wp - resistance) & |
---|
| 3612 | * q(k+1,j,i) |
---|
| 3613 | ENDIF |
---|
[1496] | 3614 | |
---|
[1691] | 3615 | ! |
---|
| 3616 | !-- Update virtual potential temperature |
---|
| 3617 | vpt(k,j,i) = pt(k,j,i) * ( 1.0_wp + 0.61_wp * q(k,j,i) ) |
---|
| 3618 | |
---|
[1496] | 3619 | ENDDO |
---|
| 3620 | ENDDO |
---|
| 3621 | |
---|
[1551] | 3622 | END SUBROUTINE calc_q_surface |
---|
[1496] | 3623 | |
---|
[1788] | 3624 | |
---|
[1496] | 3625 | END MODULE land_surface_model_mod |
---|