Changeset 4700 for palm/trunk
- Timestamp:
- Sep 25, 2020 1:08:49 PM (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/large_scale_forcing_nudging_mod.f90
r4671 r4700 1 1 !> @file large_scale_forcing_nudging_mod.f90 2 !------------------------------------------------------------------------------ !2 !--------------------------------------------------------------------------------------------------! 3 3 ! This file is part of the PALM model system. 4 4 ! 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. 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/>. 5 ! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General 6 ! Public License as published by the Free Software Foundation, either version 3 of the License, or 7 ! (at your option) any later version. 8 ! 9 ! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the 10 ! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General 11 ! Public License for more details. 12 ! 13 ! You should have received a copy of the GNU General Public License along with PALM. If not, see 14 ! <http://www.gnu.org/licenses/>. 16 15 ! 17 16 ! Copyright 1997-2020 Leibniz Universitaet Hannover 18 !------------------------------------------------------------------------------ !17 !--------------------------------------------------------------------------------------------------! 19 18 ! 20 19 ! Current revisions: … … 25 24 ! ----------------- 26 25 ! $Id$ 26 ! file re-formatted to follow the PALM coding standard 27 ! 28 ! 4671 2020-09-09 20:27:58Z pavelkrc 27 29 ! Implementation of downward facing USM and LSM surfaces 28 ! 30 ! 29 31 ! 4360 2020-01-07 11:25:50Z suehring 30 ! Introduction of wall_flags_total_0, which currently sets bits based on static 31 ! topographyinformation used in wall_flags_static_032 ! 32 ! Introduction of wall_flags_total_0, which currently sets bits based on static topography 33 ! information used in wall_flags_static_0 34 ! 33 35 ! 4329 2019-12-10 15:46:36Z motisi 34 36 ! Renamed wall_flags_0 to wall_flags_static_0 35 ! 37 ! 36 38 ! 4182 2019-08-22 15:20:23Z scharf 37 39 ! Corrected "Former revisions" section 38 ! 40 ! 39 41 ! 3719 2019-02-06 13:10:18Z kanani 40 42 ! Removed USE cpulog (unused) 41 ! 43 ! 42 44 ! 3655 2019-01-07 16:51:22Z knoop 43 45 ! unused variables removed 46 ! 44 47 ! 2320 2017-07-21 12:47:43Z suehring 45 48 ! initial revision … … 47 50 ! Description: 48 51 ! ------------ 49 !> Calculates large scale forcings (geostrophic wind and subsidence velocity) as 50 !> well as surfacesfluxes dependent on time given in an external file (LSF_DATA).51 !> Moreover, module contains nudging routines, where u, v, pt and q are nudged 52 !> to given profiles on a relaxation timescale tnudge.53 !> Profiles are read in from NUDGING_DATA. 52 !> Calculates large scale forcings (geostrophic wind and subsidence velocity) as well as surfaces 53 !> fluxes dependent on time given in an external file (LSF_DATA). 54 !> Moreover, module contains nudging routines, where u, v, pt and q are nudged to given profiles on 55 !> a relaxation timescale tnudge. 56 !> Profiles are read in from NUDGING_DATA. 54 57 !> Code is based on Neggers et al. (2012) and also in parts on DALES and UCLA-LES. 55 58 !> @todo: Revise reading of ASCII-files 56 59 !> @todo: Remove unused variables and control flags 57 60 !> @todo: Revise large-scale facing of surface variables 58 !> @todo: Revise control flags lsf_exception, lsf_surf, lsf_vert, etc. 59 !-------------------------------------------------------------------------------- !61 !> @todo: Revise control flags lsf_exception, lsf_surf, lsf_vert, etc. 62 !--------------------------------------------------------------------------------------------------! 60 63 MODULE lsf_nudging_mod 61 64 62 USE arrays_3d, & 63 ONLY: dzw, e, diss, heatflux_input_conversion, pt, pt_init, q, & 64 q_init, s, tend, u, u_init, ug, v, v_init, vg, w, w_subs, & 65 waterflux_input_conversion, zu, zw 66 67 USE control_parameters, & 68 ONLY: bc_lr, bc_ns, bc_pt_b, bc_q_b, constant_diffusion, & 69 constant_heatflux, constant_waterflux, & 70 data_output_pr, dt_3d, end_time, & 71 humidity, initializing_actions, intermediate_timestep_count, & 72 ibc_pt_b, ibc_q_b, & 73 large_scale_forcing, large_scale_subsidence, lsf_surf, lsf_vert,& 74 lsf_exception, message_string, neutral, & 75 nudging, passive_scalar, pt_surface, ocean_mode, q_surface, & 76 surface_heatflux, surface_pressure, surface_waterflux, & 77 topography, use_subsidence_tendencies 78 65 USE arrays_3d, & 66 ONLY: dzw, diss, e, heatflux_input_conversion, pt, pt_init, q, q_init, s, tend, u, u_init,& 67 ug, v, v_init, vg, w, w_subs, waterflux_input_conversion, zu, zw 68 69 USE control_parameters, & 70 ONLY: bc_lr, bc_ns, bc_pt_b, bc_q_b, constant_diffusion, constant_heatflux, & 71 constant_waterflux, data_output_pr, dt_3d, end_time, humidity, initializing_actions,& 72 intermediate_timestep_count, ibc_pt_b, ibc_q_b, & 73 large_scale_forcing, large_scale_subsidence, lsf_surf, lsf_vert, lsf_exception, & 74 message_string, neutral, nudging, passive_scalar, pt_surface, ocean_mode, q_surface,& 75 surface_heatflux, surface_pressure, surface_waterflux, topography, & 76 use_subsidence_tendencies 77 79 78 USE grid_variables 80 79 81 USE indices, &82 ONLY: nbgp, ngp_sums_ls, nx, nxl, nxlg, nxlu, nxr, nxrg, ny, nys, 83 n ysv, nysg, nyn, nyng, nzb, nz, nzt, wall_flags_total_080 USE indices, & 81 ONLY: nbgp, ngp_sums_ls, nx, nxl, nxlg, nxlu, nxr, nxrg, ny, nys, nysv, nysg, nyn, nyng, & 82 nzb, nz, nzt, wall_flags_total_0 84 83 85 84 USE kinds … … 87 86 USE pegrid 88 87 89 USE surface_mod, &88 USE surface_mod, & 90 89 ONLY: surf_def_h, surf_lsm_h, surf_usm_h 91 90 92 USE statistics, &91 USE statistics, & 93 92 ONLY: hom, statistic_regions, sums_ls_l, weight_substep 94 93 95 94 INTEGER(iwp) :: nlsf = 1000 !< maximum number of profiles in LSF_DATA (large scale forcing) 96 95 INTEGER(iwp) :: ntnudge = 1000 !< maximum number of profiles in NUDGING_DATA (nudging) 97 98 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ptnudge !< vertical profile of pot. temperature interpolated to vertical grid (nudging)99 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: qnudge !< vertical profile of water vapor mixing ratio interpolated to vertical grid (nudging)100 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tnudge !< vertical profile of nudging time scale interpolated to vertical grid (nudging)101 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: td_lsa_lpt !< temperature tendency due to large scale advection (large scale forcing)102 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: td_lsa_q !< water vapor mixing ratio tendency due to large scale advection (large scale forcing)103 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: td_sub_lpt !< temperature tendency due to subsidence/ascent (large scale forcing)104 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: td_sub_q !< water vapor mixing ratio tendency due to subsidence/ascent (large scale forcing)105 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ug_vert !< vertical profile of geostrophic wind component in x-direction interpolated to vertical grid (large scale forcing)106 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: unudge !< vertical profile of wind component in x-direction interpolated to vertical grid (nudging)107 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: vnudge !< vertical profile of wind component in y-direction interpolated to vertical grid (nudging)108 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: vg_vert !< vertical profile of geostrophic wind component in y-direction interpolated to vertical grid (large scale forcing)109 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: wnudge !< vertical profile of subsidence/ascent velocity interpolated to vertical grid (nudging) ???110 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: wsubs_vert !< vertical profile of wind component in z-direction interpolated to vertical grid (nudging) ???111 112 REAL(wp), DIMENSION(:), ALLOCATABLE :: shf_surf !< time-dependent surface sensible heat flux (large scale forcing)113 REAL(wp), DIMENSION(:), ALLOCATABLE :: timenudge !< times at which vertical profiles are defined in NUDGING_DATA (nudging)114 REAL(wp), DIMENSION(:), ALLOCATABLE :: time_surf !< times at which surface values/fluxes are defined in LSF_DATA (large scale forcing)115 REAL(wp), DIMENSION(:), ALLOCATABLE :: time_vert !< times at which vertical profiles are defined in LSF_DATA (large scale forcing)116 117 REAL(wp), DIMENSION(:), ALLOCATABLE :: tmp_tnudge !< current nudging time scale118 96 119 97 REAL(wp), DIMENSION(:), ALLOCATABLE :: p_surf !< time-dependent surface pressure (large scale forcing) … … 121 99 REAL(wp), DIMENSION(:), ALLOCATABLE :: qsws_surf !< time-dependent surface latent heat flux (large scale forcing) 122 100 REAL(wp), DIMENSION(:), ALLOCATABLE :: q_surf !< time-dependent surface water vapor mixing ratio (large scale forcing) 101 REAL(wp), DIMENSION(:), ALLOCATABLE :: shf_surf !< time-dependent surface sensible heat flux (large scale forcing) 102 REAL(wp), DIMENSION(:), ALLOCATABLE :: timenudge !< times at which vertical profiles are defined in NUDGING_DATA (nudging) 103 REAL(wp), DIMENSION(:), ALLOCATABLE :: time_surf !< times at which surface values/fluxes are defined in LSF_DATA (large 104 !< scale forcing) 105 REAL(wp), DIMENSION(:), ALLOCATABLE :: time_vert !< times at which vertical profiles are defined in LSF_DATA (large scale 106 !< forcing) 107 REAL(wp), DIMENSION(:), ALLOCATABLE :: tmp_tnudge !< current nudging time scale 108 109 110 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ptnudge !< vertical profile of pot. temperature interpolated to vertical grid 111 !< (nudging) 112 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: qnudge !< vertical profile of water vapor mixing ratio interpolated to vertical 113 !< grid (nudging) 114 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: tnudge !< vertical profile of nudging time scale interpolated to vertical grid 115 !< (nudging) 116 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: td_lsa_lpt !< temperature tendency due to large scale advection (large scale forcing) 117 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: td_lsa_q !< water vapor mixing ratio tendency due to large scale advection (large 118 !< scale forcing) 119 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: td_sub_lpt !< temperature tendency due to subsidence/ascent (large scale forcing) 120 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: td_sub_q !< water vapor mixing ratio tendency due to subsidence/ascent (large scale 121 !< forcing) 122 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: ug_vert !< vertical profile of geostrophic wind component in x-direction 123 !< interpolated to vertical grid (large scale forcing) 124 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: unudge !< vertical profile of wind component in x-direction interpolated to 125 !< vertical grid (nudging) 126 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: vnudge !< vertical profile of wind component in y-direction interpolated to 127 !< vertical grid (nudging) 128 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: vg_vert !< vertical profile of geostrophic wind component in y-direction 129 !< interpolated to vertical grid (large scale forcing) 130 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: wnudge !< vertical profile of subsidence/ascent velocity interpolated to vertical 131 !< grid (nudging) ??? 132 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: wsubs_vert !< vertical profile of wind component in z-direction interpolated to 133 !< vertical grid (nudging) ??? 123 134 124 135 SAVE … … 126 137 ! 127 138 !-- Public subroutines 128 PUBLIC calc_tnudge, ls_forcing_surf, ls_forcing_vert, ls_advec, lsf_init, & 129 lsf_nudging_check_parameters, nudge_init, & 130 lsf_nudging_check_data_output_pr, lsf_nudging_header, & 131 nudge, nudge_ref 132 139 PUBLIC calc_tnudge, ls_forcing_surf, ls_forcing_vert, ls_advec, lsf_init, & 140 lsf_nudging_check_parameters, nudge_init, lsf_nudging_check_data_output_pr, & 141 lsf_nudging_header, nudge, nudge_ref 142 133 143 ! 134 144 !-- Public variables 135 PUBLIC qsws_surf, shf_surf, td_lsa_lpt, td_lsa_q, td_sub_lpt, & 136 td_sub_q, time_vert 145 PUBLIC qsws_surf, shf_surf, td_lsa_lpt, td_lsa_q, td_sub_lpt, td_sub_q, time_vert 137 146 138 147 … … 150 159 151 160 152 !------------------------------------------------------------------------------ !161 !--------------------------------------------------------------------------------------------------! 153 162 ! Description: 154 163 ! ------------ 155 164 !> @todo Missing subroutine description. 156 !------------------------------------------------------------------------------ !165 !--------------------------------------------------------------------------------------------------! 157 166 SUBROUTINE lsf_nudging_check_parameters 158 167 … … 161 170 !-- Check nudging and large scale forcing from external file 162 171 IF ( nudging .AND. ( .NOT. large_scale_forcing ) ) THEN 163 message_string = 'Nudging requires large_scale_forcing = .T.. &'// &164 'Surface fluxes and geostrophic wind should be &'//&165 'prescribed in file LSF_DATA'172 message_string = 'Nudging requires large_scale_forcing = .T.. &'// & 173 'Surface fluxes and geostrophic wind should be &'// & 174 'prescribed in file LSF_DATA' 166 175 CALL message( 'check_parameters', 'PA0374', 1, 2, 0, 6, 0 ) 167 176 ENDIF 168 177 169 IF ( large_scale_forcing .AND. ( bc_lr /= 'cyclic' .OR. & 170 bc_ns /= 'cyclic' ) ) THEN 171 message_string = 'Non-cyclic lateral boundaries do not allow for &'//& 172 'the usage of large scale forcing from external file.' 178 IF ( large_scale_forcing .AND. ( bc_lr /= 'cyclic' .OR. bc_ns /= 'cyclic' ) ) THEN 179 message_string = 'Non-cyclic lateral boundaries do not allow for &'// & 180 'the usage of large scale forcing from external file.' 173 181 CALL message( 'check_parameters', 'PA0375', 1, 2, 0, 6, 0 ) 174 182 ENDIF 175 183 176 IF ( large_scale_forcing .AND. ( .NOT.humidity ) ) THEN177 message_string = 'The usage of large scale forcing from external &'// &178 'file LSF_DATA requires humidity = .T..'184 IF ( large_scale_forcing .AND. ( .NOT. humidity ) ) THEN 185 message_string = 'The usage of large scale forcing from external &'// & 186 'file LSF_DATA requires humidity = .T..' 179 187 CALL message( 'check_parameters', 'PA0376', 1, 2, 0, 6, 0 ) 180 188 ENDIF 181 189 182 190 IF ( large_scale_forcing .AND. passive_scalar ) THEN 183 message_string = 'The usage of large scale forcing from external &'// &184 'file LSF_DATA is not implemented for passive scalars'191 message_string = 'The usage of large scale forcing from external &'// & 192 'file LSF_DATA is not implemented for passive scalars' 185 193 CALL message( 'check_parameters', 'PA0440', 1, 2, 0, 6, 0 ) 186 194 ENDIF 187 195 188 IF ( large_scale_forcing .AND. topography /= 'flat' & 189 .AND. .NOT. lsf_exception ) THEN 190 message_string = 'The usage of large scale forcing from external &'//& 191 'file LSF_DATA is not implemented for non-flat topography' 196 IF ( large_scale_forcing .AND. topography /= 'flat' .AND. .NOT. lsf_exception ) THEN 197 message_string = 'The usage of large scale forcing from external &'// & 198 'file LSF_DATA is not implemented for non-flat topography' 192 199 CALL message( 'check_parameters', 'PA0377', 1, 2, 0, 6, 0 ) 193 200 ENDIF 194 201 195 202 IF ( large_scale_forcing .AND. ocean_mode ) THEN 196 message_string = 'The usage of large scale forcing from external &'// &197 'file LSF_DATA is not implemented for ocean mode'203 message_string = 'The usage of large scale forcing from external &'// & 204 'file LSF_DATA is not implemented for ocean mode' 198 205 CALL message( 'check_parameters', 'PA0378', 1, 2, 0, 6, 0 ) 199 206 ENDIF … … 201 208 END SUBROUTINE lsf_nudging_check_parameters 202 209 203 !------------------------------------------------------------------------------ !210 !--------------------------------------------------------------------------------------------------! 204 211 ! Description: 205 212 ! ------------ 206 213 !> Check data output of profiles for land surface model 207 !------------------------------------------------------------------------------! 208 SUBROUTINE lsf_nudging_check_data_output_pr( variable, var_count, unit, & 209 dopr_unit ) 210 214 !--------------------------------------------------------------------------------------------------! 215 SUBROUTINE lsf_nudging_check_data_output_pr( variable, var_count, unit, dopr_unit ) 216 211 217 USE profil_parameter 212 218 213 219 IMPLICIT NONE 214 215 CHARACTER (LEN=*) :: unit !< 216 CHARACTER (LEN=*) :: variable !< 220 221 CHARACTER (LEN=*) :: unit !< 222 CHARACTER (LEN=*) :: variable !< 217 223 CHARACTER (LEN=*) :: dopr_unit !< local value of dopr_unit 218 219 INTEGER(iwp) :: var_count !< 224 225 INTEGER(iwp) :: var_count !< 220 226 221 227 SELECT CASE ( TRIM( variable ) ) 222 228 223 229 224 230 CASE ( 'td_lsa_thetal' ) 225 IF ( .NOT.large_scale_forcing ) THEN226 message_string = 'data_output_pr = ' // &227 TRIM( data_output_pr(var_count) ) // &228 ' is not implemented for ' // &231 IF ( .NOT. large_scale_forcing ) THEN 232 message_string = 'data_output_pr = ' // & 233 TRIM( data_output_pr(var_count) ) // & 234 ' is not implemented for ' // & 229 235 'large_scale_forcing = .FALSE.' 230 CALL message( 'lsf_nudging_check_data_output_pr', 'PA0393', & 231 1, 2, 0, 6, 0 ) 236 CALL message( 'lsf_nudging_check_data_output_pr', 'PA0393', 1, 2, 0, 6, 0 ) 232 237 ELSE 233 238 dopr_index(var_count) = 81 … … 238 243 239 244 CASE ( 'td_lsa_q' ) 240 IF ( .NOT.large_scale_forcing ) THEN241 message_string = 'data_output_pr = ' // &242 TRIM( data_output_pr(var_count) ) // &243 ' is not implemented for ' // &245 IF ( .NOT. large_scale_forcing ) THEN 246 message_string = 'data_output_pr = ' // & 247 TRIM( data_output_pr(var_count) ) // & 248 ' is not implemented for ' // & 244 249 'large_scale_forcing = .FALSE.' 245 CALL message( 'lsf_nudging_check_data_output_pr', 'PA0393', & 246 1, 2, 0, 6, 0 ) 250 CALL message( 'lsf_nudging_check_data_output_pr', 'PA0393', 1, 2, 0, 6, 0 ) 247 251 ELSE 248 252 dopr_index(var_count) = 82 … … 252 256 ENDIF 253 257 CASE ( 'td_sub_thetal' ) 254 IF ( .NOT.large_scale_forcing ) THEN255 message_string = 'data_output_pr = ' // &256 TRIM( data_output_pr(var_count) ) // &257 ' is not implemented for ' // &258 IF ( .NOT. large_scale_forcing ) THEN 259 message_string = 'data_output_pr = ' // & 260 TRIM( data_output_pr(var_count) ) // & 261 ' is not implemented for ' // & 258 262 'large_scale_forcing = .FALSE.' 259 CALL message( 'lsf_nudging_check_data_output_pr', 'PA0393', & 260 1, 2, 0, 6, 0 ) 263 CALL message( 'lsf_nudging_check_data_output_pr', 'PA0393', 1, 2, 0, 6, 0 ) 261 264 ELSE 262 265 dopr_index(var_count) = 83 … … 267 270 268 271 CASE ( 'td_sub_q' ) 269 IF ( .NOT.large_scale_forcing ) THEN270 message_string = 'data_output_pr = ' // &271 TRIM( data_output_pr(var_count) ) // &272 ' is not implemented for ' // &272 IF ( .NOT. large_scale_forcing ) THEN 273 message_string = 'data_output_pr = ' // & 274 TRIM( data_output_pr(var_count) ) // & 275 ' is not implemented for ' // & 273 276 'large_scale_forcing = .FALSE.' 274 CALL message( 'lsf_nudging_check_data_output_pr', 'PA0393', & 275 1, 2, 0, 6, 0 ) 277 CALL message( 'lsf_nudging_check_data_output_pr', 'PA0393', 1, 2, 0, 6, 0 ) 276 278 ELSE 277 279 dopr_index(var_count) = 84 … … 282 284 283 285 CASE ( 'td_nud_thetal' ) 284 IF ( .NOT.nudging ) THEN285 message_string = 'data_output_pr = ' // &286 TRIM( data_output_pr(var_count) ) // &287 ' is not implemented for ' // &286 IF ( .NOT. nudging ) THEN 287 message_string = 'data_output_pr = ' // & 288 TRIM( data_output_pr(var_count) ) // & 289 ' is not implemented for ' // & 288 290 'nudging = .FALSE.' 289 CALL message( 'lsf_nudging_check_data_output_pr', 'PA0394', & 290 1, 2, 0, 6, 0 ) 291 CALL message( 'lsf_nudging_check_data_output_pr', 'PA0394', 1, 2, 0, 6, 0 ) 291 292 ELSE 292 293 dopr_index(var_count) = 85 … … 297 298 298 299 CASE ( 'td_nud_q' ) 299 IF ( .NOT.nudging ) THEN300 message_string = 'data_output_pr = ' // &301 TRIM( data_output_pr(var_count) ) // &302 ' is not implemented for ' // &300 IF ( .NOT. nudging ) THEN 301 message_string = 'data_output_pr = ' // & 302 TRIM( data_output_pr(var_count) ) // & 303 ' is not implemented for ' // & 303 304 'nudging = .FALSE.' 304 CALL message( 'lsf_nudging_check_data_output_pr', 'PA0394', & 305 1, 2, 0, 6, 0 ) 305 CALL message( 'lsf_nudging_check_data_output_pr', 'PA0394', 1, 2, 0, 6, 0 ) 306 306 ELSE 307 307 dopr_index(var_count) = 86 … … 312 312 313 313 CASE ( 'td_nud_u' ) 314 IF ( .NOT.nudging ) THEN315 message_string = 'data_output_pr = ' // &316 TRIM( data_output_pr(var_count) ) // &317 ' is not implemented for ' // &314 IF ( .NOT. nudging ) THEN 315 message_string = 'data_output_pr = ' // & 316 TRIM( data_output_pr(var_count) ) // & 317 ' is not implemented for ' // & 318 318 'nudging = .FALSE.' 319 CALL message( 'lsf_nudging_check_data_output_pr', 'PA0394', & 320 1, 2, 0, 6, 0 ) 319 CALL message( 'lsf_nudging_check_data_output_pr', 'PA0394', 1, 2, 0, 6, 0 ) 321 320 ELSE 322 321 dopr_index(var_count) = 87 … … 327 326 328 327 CASE ( 'td_nud_v' ) 329 IF ( .NOT.nudging ) THEN330 message_string = 'data_output_pr = ' // &331 TRIM( data_output_pr(var_count) ) // &332 ' is not implemented for ' // &328 IF ( .NOT. nudging ) THEN 329 message_string = 'data_output_pr = ' // & 330 TRIM( data_output_pr(var_count) ) // & 331 ' is not implemented for ' // & 333 332 'nudging = .FALSE.' 334 CALL message( 'lsf_nudging_check_data_output_pr', 'PA0394', & 335 1, 2, 0, 6, 0 ) 333 CALL message( 'lsf_nudging_check_data_output_pr', 'PA0394', 1, 2, 0, 6, 0 ) 336 334 ELSE 337 335 dopr_index(var_count) = 88 … … 344 342 CASE DEFAULT 345 343 unit = 'illegal' 346 344 347 345 END SELECT 348 346 349 347 END SUBROUTINE lsf_nudging_check_data_output_pr 350 348 351 !------------------------------------------------------------------------------ !349 !--------------------------------------------------------------------------------------------------! 352 350 ! Description: 353 351 ! ------------ 354 352 !> @todo Missing subroutine description. 355 !------------------------------------------------------------------------------ !353 !--------------------------------------------------------------------------------------------------! 356 354 SUBROUTINE lsf_nudging_header ( io ) 357 355 … … 370 368 ELSE 371 369 WRITE ( io, 6 ) 372 ENDIF 370 ENDIF 373 371 ENDIF 374 372 … … 399 397 400 398 401 1 FORMAT (//' Large scale forcing and nudging:'/ &399 1 FORMAT (//' Large scale forcing and nudging:'/ & 402 400 ' -------------------------------'/) 403 401 2 FORMAT (' --> No large scale forcing from external is used (default) ') … … 407 405 6 FORMAT (' - large scale subsidence tendencies ') 408 406 7 FORMAT (' - and geostrophic wind components ug and vg') 409 8 FORMAT (' --> Large-scale vertical motion is used in the ', &407 8 FORMAT (' --> Large-scale vertical motion is used in the ', & 410 408 'prognostic equation(s) for') 411 409 9 FORMAT (' the scalar(s) only') … … 417 415 15 FORMAT (' - prescribed surface fluxes for humidity') 418 416 419 END SUBROUTINE lsf_nudging_header 420 421 !------------------------------------------------------------------------------ !417 END SUBROUTINE lsf_nudging_header 418 419 !--------------------------------------------------------------------------------------------------! 422 420 ! Description: 423 421 ! ------------ 424 422 !> @todo Missing subroutine description. 425 !------------------------------------------------------------------------------ !423 !--------------------------------------------------------------------------------------------------! 426 424 SUBROUTINE lsf_init 427 425 … … 437 435 438 436 REAL(wp) :: fac !< 437 REAL(wp) :: high_td_lsa_lpt !< 438 REAL(wp) :: high_td_lsa_q !< 439 REAL(wp) :: high_td_sub_lpt !< 440 REAL(wp) :: high_td_sub_q !< 439 441 REAL(wp) :: highheight !< 440 442 REAL(wp) :: highug_vert !< 441 443 REAL(wp) :: highvg_vert !< 442 444 REAL(wp) :: highwsubs_vert !< 445 REAL(wp) :: low_td_lsa_lpt !< 446 REAL(wp) :: low_td_lsa_q !< 447 REAL(wp) :: low_td_sub_lpt !< 448 REAL(wp) :: low_td_sub_q !< 443 449 REAL(wp) :: lowheight !< 444 450 REAL(wp) :: lowug_vert !< 445 451 REAL(wp) :: lowvg_vert !< 446 452 REAL(wp) :: lowwsubs_vert !< 447 REAL(wp) :: high_td_lsa_lpt !<448 REAL(wp) :: low_td_lsa_lpt !<449 REAL(wp) :: high_td_lsa_q !<450 REAL(wp) :: low_td_lsa_q !<451 REAL(wp) :: high_td_sub_lpt !<452 REAL(wp) :: low_td_sub_lpt !<453 REAL(wp) :: high_td_sub_q !<454 REAL(wp) :: low_td_sub_q !<455 453 REAL(wp) :: r_dummy !< 456 454 457 ALLOCATE( p_surf(0:nlsf), pt_surf(0:nlsf), q_surf(0:nlsf), &458 qsws_surf(0:nlsf), shf_surf(0:nlsf), &459 td_lsa_lpt(nzb:nzt+1,0:nlsf), td_lsa_q(nzb:nzt+1,0:nlsf), &460 td_sub_lpt(nzb:nzt+1,0:nlsf), td_sub_q(nzb:nzt+1,0:nlsf), &461 time_vert(0:nlsf), time_surf(0:nlsf), &462 ug_vert(nzb:nzt+1,0:nlsf), vg_vert(nzb:nzt+1,0:nlsf), &455 ALLOCATE( p_surf(0:nlsf), pt_surf(0:nlsf), q_surf(0:nlsf), & 456 qsws_surf(0:nlsf), shf_surf(0:nlsf), & 457 td_lsa_lpt(nzb:nzt+1,0:nlsf), td_lsa_q(nzb:nzt+1,0:nlsf), & 458 td_sub_lpt(nzb:nzt+1,0:nlsf), td_sub_q(nzb:nzt+1,0:nlsf), & 459 time_vert(0:nlsf), time_surf(0:nlsf), & 460 ug_vert(nzb:nzt+1,0:nlsf), vg_vert(nzb:nzt+1,0:nlsf), & 463 461 wsubs_vert(nzb:nzt+1,0:nlsf) ) 464 462 … … 470 468 471 469 ! 472 !-- Array for storing large scale forcing and nudging tendencies at each 473 !-- timestep for data output 470 !-- Array for storing large scale forcing and nudging tendencies at each timestep for data output 474 471 ALLOCATE( sums_ls_l(nzb:nzt+1,0:7) ) 475 472 sums_ls_l = 0.0_wp … … 477 474 ngp_sums_ls = (nz+2)*6 478 475 479 OPEN ( finput, FILE='LSF_DATA', STATUS='OLD', & 480 FORM='FORMATTED', IOSTAT=ierrn ) 476 OPEN ( finput, FILE='LSF_DATA', STATUS='OLD', FORM='FORMATTED', IOSTAT=ierrn ) 481 477 482 478 IF ( ierrn /= 0 ) THEN … … 504 500 DO WHILE ( time_surf(nt) < end_time ) 505 501 nt = nt + 1 506 READ ( finput, *, IOSTAT = ierrn ) time_surf(nt), shf_surf(nt), & 507 qsws_surf(nt), pt_surf(nt), & 508 q_surf(nt), p_surf(nt) 502 READ ( finput, *, IOSTAT = ierrn ) time_surf(nt), shf_surf(nt), qsws_surf(nt), & 503 pt_surf(nt), q_surf(nt), p_surf(nt) 509 504 510 505 IF ( ierrn /= 0 ) THEN 511 WRITE ( message_string, * ) 'No time dependent surface ' // &512 'variables in & LSF_DATA for end of run found'506 WRITE ( message_string, * ) 'No time dependent surface ' // & 507 'variables in & LSF_DATA for end of run found' 513 508 514 509 CALL message( 'ls_forcing', 'PA0363', 1, 2, 0, 6, 0 ) … … 517 512 518 513 IF ( time_surf(1) > end_time ) THEN 519 WRITE ( message_string, * ) 'Time dependent surface variables in ' // &520 '&LSF_DATA set in after end of ' , &514 WRITE ( message_string, * ) 'Time dependent surface variables in ' // & 515 '&LSF_DATA set in after end of ' , & 521 516 'simulation - lsf_surf is set to FALSE' 522 517 CALL message( 'ls_forcing', 'PA0371', 0, 0, 0, 6, 0 ) … … 539 534 ierrn = 1 ! not zero 540 535 ! 541 !-- Search for the next line consisting of "# time", 542 !-- from there onwards the profiles will be read543 DO WHILE ( .NOT. ( hash == "#" .AND. ierrn == 0 ) ) 536 !-- Search for the next line consisting of "# time", from there onwards the profiles will be 537 !-- read. 538 DO WHILE ( .NOT. ( hash == "#" .AND. ierrn == 0 ) ) 544 539 READ ( finput, *, IOSTAT=ierrn ) hash, time_vert(nt) 545 IF ( ierrn < 0 ) THEN 546 WRITE( message_string, * ) 'No time dependent vertical profiles', &547 ' in & LSF_DATA for end of run found'540 IF ( ierrn < 0 ) THEN 541 WRITE( message_string, * ) 'No time dependent vertical profiles', & 542 ' in & LSF_DATA for end of run found' 548 543 CALL message( 'ls_forcing', 'PA0372', 1, 2, 0, 6, 0 ) 549 544 ENDIF … … 552 547 IF ( nt == 1 .AND. time_vert(nt) > end_time ) EXIT 553 548 554 READ ( finput, *, IOSTAT=ierrn ) lowheight, lowug_vert, lowvg_vert, & 555 lowwsubs_vert, low_td_lsa_lpt, & 556 low_td_lsa_q, low_td_sub_lpt, & 549 READ ( finput, *, IOSTAT=ierrn ) lowheight, lowug_vert, lowvg_vert, lowwsubs_vert, & 550 low_td_lsa_lpt, low_td_lsa_q, low_td_sub_lpt, & 557 551 low_td_sub_q 558 552 IF ( ierrn /= 0 ) THEN … … 561 555 ENDIF 562 556 563 READ ( finput, *, IOSTAT=ierrn ) highheight, highug_vert, & 564 highvg_vert, highwsubs_vert, & 565 high_td_lsa_lpt, high_td_lsa_q, & 566 high_td_sub_lpt, high_td_sub_q 567 557 READ ( finput, *, IOSTAT=ierrn ) highheight, highug_vert, highvg_vert, highwsubs_vert, & 558 high_td_lsa_lpt, high_td_lsa_q, high_td_sub_lpt, & 559 high_td_sub_q 560 568 561 IF ( ierrn /= 0 ) THEN 569 562 message_string = 'errors in file LSF_DATA' … … 584 577 585 578 ierrn = 0 586 READ ( finput, *, IOSTAT=ierrn ) highheight, highug_vert, & 587 highvg_vert, highwsubs_vert, & 588 high_td_lsa_lpt, & 589 high_td_lsa_q, & 579 READ ( finput, *, IOSTAT=ierrn ) highheight, highug_vert, highvg_vert, & 580 highwsubs_vert, high_td_lsa_lpt, high_td_lsa_q, & 590 581 high_td_sub_lpt, high_td_sub_q 591 582 592 583 IF ( ierrn /= 0 ) THEN 593 WRITE( message_string, * ) 'zu(',k,') = ', zu(k), 'm ', &594 'is higher than the maximum height in LSF_DATA ',&595 'which is ', lowheight, 'm. Interpolation on PALM ',&596 'grid is not possible.'584 WRITE( message_string, * ) 'zu(',k,') = ', zu(k), 'm ', & 585 'is higher than the maximum height in LSF_DATA ', & 586 'which is ', lowheight, 'm. Interpolation on PALM ', & 587 'grid is not possible.' 597 588 CALL message( 'ls_forcing', 'PA0395', 1, 2, 0, 6, 0 ) 598 589 ENDIF … … 601 592 602 593 ! 603 !-- Interpolation of prescribed profiles in space 594 !-- Interpolation of prescribed profiles in space 604 595 fac = (highheight-zu(k))/(highheight - lowheight) 605 596 606 ug_vert(k,nt) = fac * lowug_vert & 607 + ( 1.0_wp - fac ) * highug_vert 608 vg_vert(k,nt) = fac * lowvg_vert & 609 + ( 1.0_wp - fac ) * highvg_vert 610 wsubs_vert(k,nt) = fac * lowwsubs_vert & 611 + ( 1.0_wp - fac ) * highwsubs_vert 612 613 td_lsa_lpt(k,nt) = fac * low_td_lsa_lpt & 614 + ( 1.0_wp - fac ) * high_td_lsa_lpt 615 td_lsa_q(k,nt) = fac * low_td_lsa_q & 616 + ( 1.0_wp - fac ) * high_td_lsa_q 617 td_sub_lpt(k,nt) = fac * low_td_sub_lpt & 618 + ( 1.0_wp - fac ) * high_td_sub_lpt 619 td_sub_q(k,nt) = fac * low_td_sub_q & 620 + ( 1.0_wp - fac ) * high_td_sub_q 597 ug_vert(k,nt) = fac * lowug_vert + ( 1.0_wp - fac ) * highug_vert 598 vg_vert(k,nt) = fac * lowvg_vert + ( 1.0_wp - fac ) * highvg_vert 599 wsubs_vert(k,nt) = fac * lowwsubs_vert + ( 1.0_wp - fac ) * highwsubs_vert 600 601 td_lsa_lpt(k,nt) = fac * low_td_lsa_lpt + ( 1.0_wp - fac ) * high_td_lsa_lpt 602 td_lsa_q(k,nt) = fac * low_td_lsa_q + ( 1.0_wp - fac ) * high_td_lsa_q 603 td_sub_lpt(k,nt) = fac * low_td_sub_lpt + ( 1.0_wp - fac ) * high_td_sub_lpt 604 td_sub_q(k,nt) = fac * low_td_sub_q + ( 1.0_wp - fac ) * high_td_sub_q 621 605 622 606 ENDDO 623 607 624 ENDDO 608 ENDDO 625 609 626 610 ! 627 611 !-- Large scale vertical velocity has to be zero at the surface 628 612 wsubs_vert(nzb,:) = 0.0_wp 629 613 630 614 IF ( time_vert(1) > end_time ) THEN 631 WRITE ( message_string, * ) 'Time dependent large scale profile ', &632 'forcing from&LSF_DATA sets in after end of ' ,&633 'simulation - lsf_vert is set to FALSE'615 WRITE ( message_string, * ) 'Time dependent large scale profile ', & 616 'forcing from&LSF_DATA sets in after end of ' , & 617 'simulation - lsf_vert is set to FALSE' 634 618 CALL message( 'ls_forcing', 'PA0373', 0, 0, 0, 6, 0 ) 635 619 lsf_vert = .FALSE. … … 640 624 END SUBROUTINE lsf_init 641 625 642 !------------------------------------------------------------------------------ !626 !--------------------------------------------------------------------------------------------------! 643 627 ! Description: 644 628 ! ------------ 645 629 !> @todo Missing subroutine description. 646 !------------------------------------------------------------------------------ !630 !--------------------------------------------------------------------------------------------------! 647 631 SUBROUTINE ls_forcing_surf ( time ) 648 632 … … 669 653 IF ( ibc_pt_b == 0 ) THEN 670 654 ! 671 !-- In case of Dirichlet boundary condition shf must not 672 !-- be set - it is calculated via MOST inprandtl_fluxes655 !-- In case of Dirichlet boundary condition shf must not be set - it is calculated via MOST in 656 !-- prandtl_fluxes 673 657 pt_surface = pt_surf(nt) + fac * ( pt_surf(nt+1) - pt_surf(nt) ) 674 658 675 659 ELSEIF ( ibc_pt_b == 1 ) THEN 676 660 ! 677 !-- In case of Neumann boundary condition pt_surface is needed for 661 !-- In case of Neumann boundary condition pt_surface is needed for 678 662 !-- calculation of reference density 679 dum_surf_flux = ( shf_surf(nt) + fac * & 680 ( shf_surf(nt+1) - shf_surf(nt) ) & 663 dum_surf_flux = ( shf_surf(nt) + fac * ( shf_surf(nt+1) - shf_surf(nt) ) & 681 664 ) * heatflux_input_conversion(nzb) 682 665 ! 683 !-- Save surface sensible heat flux on default, natural and urban surface 684 !-- type, if required 666 !-- Save surface sensible heat flux on default, natural and urban surface type, if required. 685 667 IF ( surf_def_h(0)%ns >= 1 ) surf_def_h(0)%shf(:) = dum_surf_flux 686 668 IF ( surf_lsm_h(0)%ns >= 1 ) surf_lsm_h(0)%shf(:) = dum_surf_flux … … 693 675 IF ( ibc_q_b == 0 ) THEN 694 676 ! 695 !-- In case of Dirichlet boundary condition qsws must not 696 !-- be set - it is calculated via MOSTin prandtl_fluxes677 !-- In case of Dirichlet boundary condition qsws must not be set - it is calculated via MOST 678 !-- in prandtl_fluxes 697 679 q_surface = q_surf(nt) + fac * ( q_surf(nt+1) - q_surf(nt) ) 698 680 699 681 ELSEIF ( ibc_q_b == 1 ) THEN 700 dum_surf_flux = ( qsws_surf(nt) + fac * & 701 ( qsws_surf(nt+1) - qsws_surf(nt) ) & 702 ) * waterflux_input_conversion(nzb) 703 ! 704 !-- Save surface sensible heat flux on default, natural and urban surface 705 !-- type, if required 682 dum_surf_flux = ( qsws_surf(nt) + fac * ( qsws_surf(nt+1) - qsws_surf(nt) ) & 683 ) * waterflux_input_conversion(nzb) 684 ! 685 !-- Save surface sensible heat flux on default, natural and urban surface type, if required 706 686 IF ( surf_def_h(0)%ns >= 1 ) surf_def_h(0)%qsws(:) = dum_surf_flux 707 687 IF ( surf_lsm_h(0)%ns >= 1 ) surf_lsm_h(0)%qsws(:) = dum_surf_flux … … 710 690 ENDIF 711 691 ! 712 !-- Surface heat- and waterflux will be written later onto surface elements 713 IF ( .NOT. neutral .AND. constant_heatflux .AND.&692 !-- Surface heat- and waterflux will be written later onto surface elements 693 IF ( .NOT. neutral .AND. constant_heatflux .AND. & 714 694 TRIM( initializing_actions ) /= 'read_restart_data' ) THEN 715 716 ENDIF 717 IF ( humidity .AND. constant_waterflux .AND. &695 surface_heatflux = shf_surf(1) 696 ENDIF 697 IF ( humidity .AND. constant_waterflux .AND. & 718 698 TRIM( initializing_actions ) /= 'read_restart_data' ) THEN 719 699 surface_waterflux = qsws_surf(1) 720 700 ENDIF 721 701 722 702 surface_pressure = p_surf(nt) + fac * ( p_surf(nt+1) - p_surf(nt) ) 723 703 724 END SUBROUTINE ls_forcing_surf 725 726 727 728 729 !------------------------------------------------------------------------------ !704 END SUBROUTINE ls_forcing_surf 705 706 707 708 709 !--------------------------------------------------------------------------------------------------! 730 710 ! Description: 731 711 ! ------------ 732 712 !> @todo Missing subroutine description. 733 !------------------------------------------------------------------------------ !713 !--------------------------------------------------------------------------------------------------! 734 714 SUBROUTINE ls_forcing_vert ( time ) 735 715 … … 758 738 759 739 IF ( large_scale_subsidence ) THEN 760 w_subs = wsubs_vert(:,nt) & 761 + fac * ( wsubs_vert(:,nt+1) - wsubs_vert(:,nt) ) 740 w_subs = wsubs_vert(:,nt) + fac * ( wsubs_vert(:,nt+1) - wsubs_vert(:,nt) ) 762 741 ENDIF 763 742 … … 765 744 766 745 767 !------------------------------------------------------------------------------ !746 !--------------------------------------------------------------------------------------------------! 768 747 ! Description: 769 748 ! ------------ 770 749 !> Call for all grid points 771 !------------------------------------------------------------------------------ !750 !--------------------------------------------------------------------------------------------------! 772 751 SUBROUTINE ls_advec ( time, prog_var ) 773 752 774 753 775 754 IMPLICIT NONE 776 755 777 CHARACTER (LEN=*) :: prog_var !< 778 779 REAL(wp), INTENT(in) :: time !<780 REAL(wp) :: fac !<781 782 INTEGER(iwp) :: i !<783 INTEGER(iwp) :: j !< 784 INTEGER(iwp) :: k !<785 INTEGER(iwp) :: nt !<786 787 ! 788 !-- Interpolation in time of LSF_DATA 756 CHARACTER (LEN=*) :: prog_var !< 757 758 INTEGER(iwp) :: i !< 759 INTEGER(iwp) :: j !< 760 INTEGER(iwp) :: k !< 761 INTEGER(iwp) :: nt !< 762 763 REAL(wp) :: fac !< 764 REAL(wp), INTENT(in) :: time !< 765 766 ! 767 !-- Interpolation in time of LSF_DATA 789 768 nt = 1 790 769 DO WHILE ( time > time_vert(nt) ) … … 798 777 799 778 ! 800 !-- Add horizontal large scale advection tendencies of pt and q 779 !-- Add horizontal large scale advection tendencies of pt and q 801 780 SELECT CASE ( prog_var ) 802 781 … … 806 785 DO j = nys, nyn 807 786 DO k = nzb+1, nzt 808 tend(k,j,i) = tend(k,j,i) + td_lsa_lpt(k,nt) + fac * &809 ( td_lsa_lpt(k,nt+1) - td_lsa_lpt(k,nt) ) *&810 MERGE( 1.0_wp, 0.0_wp,&811 BTEST( wall_flags_total_0(k,j,i), 0 ) )787 tend(k,j,i) = tend(k,j,i) + td_lsa_lpt(k,nt) + fac * & 788 ( td_lsa_lpt(k,nt+1) - td_lsa_lpt(k,nt) ) * & 789 MERGE( 1.0_wp, 0.0_wp, & 790 BTEST( wall_flags_total_0(k,j,i), 0 ) ) 812 791 ENDDO 813 792 ENDDO … … 819 798 DO j = nys, nyn 820 799 DO k = nzb+1, nzt 821 tend(k,j,i) = tend(k,j,i) + td_lsa_q(k,nt) + fac * &822 ( td_lsa_q(k,nt+1) - td_lsa_q(k,nt) ) *&823 MERGE( 1.0_wp, 0.0_wp,&824 BTEST( wall_flags_total_0(k,j,i), 0 ) )800 tend(k,j,i) = tend(k,j,i) + td_lsa_q(k,nt) + fac * & 801 ( td_lsa_q(k,nt+1) - td_lsa_q(k,nt) ) * & 802 MERGE( 1.0_wp, 0.0_wp, & 803 BTEST( wall_flags_total_0(k,j,i), 0 ) ) 825 804 ENDDO 826 805 ENDDO … … 840 819 DO j = nys, nyn 841 820 DO k = nzb+1, nzt 842 tend(k,j,i) = tend(k,j,i) + td_sub_lpt(k,nt) + fac * &843 ( td_sub_lpt(k,nt+1) - td_sub_lpt(k,nt) )*&844 MERGE( 1.0_wp, 0.0_wp,&845 BTEST( wall_flags_total_0(k,j,i), 0 ) )821 tend(k,j,i) = tend(k,j,i) + td_sub_lpt(k,nt) + fac * & 822 ( td_sub_lpt(k,nt+1) - td_sub_lpt(k,nt) ) * & 823 MERGE( 1.0_wp, 0.0_wp, & 824 BTEST( wall_flags_total_0(k,j,i), 0 ) ) 846 825 ENDDO 847 826 ENDDO 848 827 ENDDO 849 828 850 829 CASE ( 'q' ) 851 830 … … 853 832 DO j = nys, nyn 854 833 DO k = nzb+1, nzt 855 tend(k,j,i) = tend(k,j,i) + td_sub_q(k,nt) + fac * &856 ( td_sub_q(k,nt+1) - td_sub_q(k,nt) ) *&857 MERGE( 1.0_wp, 0.0_wp,&858 BTEST( wall_flags_total_0(k,j,i), 0 ) )834 tend(k,j,i) = tend(k,j,i) + td_sub_q(k,nt) + fac * & 835 ( td_sub_q(k,nt+1) - td_sub_q(k,nt) ) * & 836 MERGE( 1.0_wp, 0.0_wp, & 837 BTEST( wall_flags_total_0(k,j,i), 0 ) ) 859 838 ENDDO 860 839 ENDDO … … 868 847 869 848 870 !------------------------------------------------------------------------------ !849 !--------------------------------------------------------------------------------------------------! 871 850 ! Description: 872 851 ! ------------ 873 852 !> Call for grid point i,j 874 !------------------------------------------------------------------------------ !853 !--------------------------------------------------------------------------------------------------! 875 854 SUBROUTINE ls_advec_ij ( i, j, time, prog_var ) 876 855 877 856 IMPLICIT NONE 878 857 879 CHARACTER (LEN=*) :: prog_var !< 880 881 REAL(wp), INTENT(in) :: time !< 882 REAL(wp) :: fac !< 883 884 INTEGER(iwp) :: i !< 885 INTEGER(iwp) :: j !< 886 INTEGER(iwp) :: k !< 887 INTEGER(iwp) :: nt !< 888 889 ! 890 !-- Interpolation in time of LSF_DATA 858 CHARACTER (LEN=*) :: prog_var !< 859 860 INTEGER(iwp) :: i !< 861 INTEGER(iwp) :: j !< 862 INTEGER(iwp) :: k !< 863 INTEGER(iwp) :: nt !< 864 865 REAL(wp) :: fac !< 866 REAL(wp), INTENT(in) :: time !< 867 868 869 ! 870 !-- Interpolation in time of LSF_DATA 891 871 nt = 1 892 872 DO WHILE ( time > time_vert(nt) ) … … 900 880 901 881 ! 902 !-- Add horizontal large scale advection tendencies of pt and q 882 !-- Add horizontal large scale advection tendencies of pt and q 903 883 SELECT CASE ( prog_var ) 904 884 … … 906 886 907 887 DO k = nzb+1, nzt 908 tend(k,j,i) = tend(k,j,i) + td_lsa_lpt(k,nt) &909 + fac * ( td_lsa_lpt(k,nt+1) - td_lsa_lpt(k,nt) )* &910 MERGE( 1.0_wp, 0.0_wp, &911 BTEST( wall_flags_total_0(k,j,i), 0 ) )888 tend(k,j,i) = tend(k,j,i) + td_lsa_lpt(k,nt) & 889 + fac * ( td_lsa_lpt(k,nt+1) - td_lsa_lpt(k,nt) )* & 890 MERGE( 1.0_wp, 0.0_wp, & 891 BTEST( wall_flags_total_0(k,j,i), 0 ) ) 912 892 ENDDO 913 893 … … 915 895 916 896 DO k = nzb+1, nzt 917 tend(k,j,i) = tend(k,j,i) + td_lsa_q(k,nt) &918 + fac * ( td_lsa_q(k,nt+1) - td_lsa_q(k,nt) ) * &919 MERGE( 1.0_wp, 0.0_wp, &920 BTEST( wall_flags_total_0(k,j,i), 0 ) )897 tend(k,j,i) = tend(k,j,i) + td_lsa_q(k,nt) & 898 + fac * ( td_lsa_q(k,nt+1) - td_lsa_q(k,nt) ) * & 899 MERGE( 1.0_wp, 0.0_wp, & 900 BTEST( wall_flags_total_0(k,j,i), 0 ) ) 921 901 ENDDO 922 902 … … 932 912 933 913 DO k = nzb+1, nzt 934 tend(k,j,i) = tend(k,j,i) + td_sub_lpt(k,nt) + fac *&935 ( td_sub_lpt(k,nt+1) - td_sub_lpt(k,nt) ) *&936 MERGE( 1.0_wp, 0.0_wp,&937 BTEST( wall_flags_total_0(k,j,i), 0 ) )914 tend(k,j,i) = tend(k,j,i) + td_sub_lpt(k,nt) & 915 + fac * ( td_sub_lpt(k,nt+1) - td_sub_lpt(k,nt) ) * & 916 MERGE( 1.0_wp, 0.0_wp, & 917 BTEST( wall_flags_total_0(k,j,i), 0 ) ) 938 918 ENDDO 939 919 940 920 CASE ( 'q' ) 941 921 942 922 DO k = nzb+1, nzt 943 tend(k,j,i) = tend(k,j,i) + td_sub_q(k,nt) + fac *&944 ( td_sub_q(k,nt+1) - td_sub_q(k,nt) ) *&945 MERGE( 1.0_wp, 0.0_wp,&946 BTEST( wall_flags_total_0(k,j,i), 0 ) )923 tend(k,j,i) = tend(k,j,i) + td_sub_q(k,nt) & 924 + fac * ( td_sub_q(k,nt+1) - td_sub_q(k,nt) ) * & 925 MERGE( 1.0_wp, 0.0_wp, & 926 BTEST( wall_flags_total_0(k,j,i), 0 ) ) 947 927 ENDDO 948 928 … … 954 934 955 935 956 !------------------------------------------------------------------------------ !936 !--------------------------------------------------------------------------------------------------! 957 937 ! Description: 958 938 ! ------------ 959 939 !> @todo Missing subroutine description. 960 !------------------------------------------------------------------------------ !940 !--------------------------------------------------------------------------------------------------! 961 941 SUBROUTINE nudge_init 962 942 963 943 IMPLICIT NONE 964 944 945 CHARACTER(1) :: hash !< 965 946 966 947 INTEGER(iwp) :: finput = 90 !< … … 969 950 INTEGER(iwp) :: nt !< 970 951 971 CHARACTER(1) :: hash!<952 REAL(wp) :: fac !< 972 953 973 954 REAL(wp) :: highheight !< 974 955 REAL(wp) :: highqnudge !< 975 956 REAL(wp) :: highptnudge !< 957 REAL(wp) :: hightnudge !< 976 958 REAL(wp) :: highunudge !< 977 959 REAL(wp) :: highvnudge !< 978 960 REAL(wp) :: highwnudge !< 979 REAL(wp) :: hightnudge !<980 961 981 962 REAL(wp) :: lowheight !< 982 963 REAL(wp) :: lowqnudge !< 983 964 REAL(wp) :: lowptnudge !< 965 REAL(wp) :: lowtnudge !< 984 966 REAL(wp) :: lowunudge !< 985 967 REAL(wp) :: lowvnudge !< 986 968 REAL(wp) :: lowwnudge !< 987 REAL(wp) :: lowtnudge !< 988 989 REAL(wp) :: fac !< 990 991 ALLOCATE( ptnudge(nzb:nzt+1,1:ntnudge), qnudge(nzb:nzt+1,1:ntnudge), & 992 tnudge(nzb:nzt+1,1:ntnudge), unudge(nzb:nzt+1,1:ntnudge), & 969 970 971 ALLOCATE( ptnudge(nzb:nzt+1,1:ntnudge), qnudge(nzb:nzt+1,1:ntnudge), & 972 tnudge(nzb:nzt+1,1:ntnudge), unudge(nzb:nzt+1,1:ntnudge), & 993 973 vnudge(nzb:nzt+1,1:ntnudge), wnudge(nzb:nzt+1,1:ntnudge) ) 994 974 … … 1019 999 ierrn = 1 ! not zero 1020 1000 ! 1021 !-- Search for the next line consisting of "# time", 1022 !-- from there onwards the profiles will be read1023 DO WHILE ( .NOT. ( hash == "#" .AND. ierrn == 0 ) ) 1024 1001 !-- Search for the next line consisting of "# time", from there onwards the profiles will be 1002 !-- read. 1003 DO WHILE ( .NOT. ( hash == "#" .AND. ierrn == 0 ) ) 1004 1025 1005 READ ( finput, *, IOSTAT=ierrn ) hash, timenudge(nt) 1026 1006 IF ( ierrn < 0 ) EXIT rloop … … 1029 1009 1030 1010 ierrn = 0 1031 READ ( finput, *, IOSTAT=ierrn ) lowheight, lowtnudge, lowunudge, & 1032 lowvnudge, lowwnudge , lowptnudge, & 1033 lowqnudge 1011 READ ( finput, *, IOSTAT=ierrn ) lowheight, lowtnudge, lowunudge, lowvnudge, lowwnudge , & 1012 lowptnudge, lowqnudge 1034 1013 1035 1014 IF ( ierrn /= 0 ) THEN … … 1039 1018 1040 1019 ierrn = 0 1041 READ ( finput, *, IOSTAT=ierrn ) highheight, hightnudge, highunudge, & 1042 highvnudge, highwnudge , highptnudge, & 1043 highqnudge 1020 READ ( finput, *, IOSTAT=ierrn ) highheight, hightnudge, highunudge, highvnudge, & 1021 highwnudge , highptnudge, highqnudge 1044 1022 1045 1023 IF ( ierrn /= 0 ) THEN … … 1057 1035 lowptnudge = highptnudge 1058 1036 lowqnudge = highqnudge 1059 1037 1060 1038 ierrn = 0 1061 READ ( finput, *, IOSTAT=ierrn ) highheight , hightnudge , & 1062 highunudge , highvnudge , & 1063 highwnudge , highptnudge, & 1064 highqnudge 1039 READ ( finput, *, IOSTAT=ierrn ) highheight , hightnudge , highunudge , & 1040 highvnudge , highwnudge , highptnudge, highqnudge 1065 1041 IF (ierrn /= 0 ) THEN 1066 WRITE( message_string, * ) 'zu(',k,') = ', zu(k), 'm is ', &1067 'higher than the maximum height in NUDING_DATA which ',&1068 'is ', lowheight, 'm. Interpolation on PALM ', &1069 'grid is not possible.'1042 WRITE( message_string, * ) 'zu(',k,') = ', zu(k), 'm is ', & 1043 'higher than the maximum height in NUDING_DATA which ',& 1044 'is ', lowheight, 'm. Interpolation on PALM ', & 1045 'grid is not possible.' 1070 1046 CALL message( 'nudging', 'PA0364', 1, 2, 0, 6, 0 ) 1071 1047 ENDIF … … 1073 1049 1074 1050 ! 1075 !-- Interpolation of prescribed profiles in space 1051 !-- Interpolation of prescribed profiles in space 1076 1052 1077 1053 fac = ( highheight - zu(k) ) / ( highheight - lowheight ) … … 1099 1075 ENDIF 1100 1076 1101 WRITE( message_string, * ) 'Initial profiles of u, v, pt and q ', &1077 WRITE( message_string, * ) 'Initial profiles of u, v, pt and q ', & 1102 1078 'from NUDGING_DATA are used.' 1103 1079 CALL message( 'large_scale_forcing_nudging', 'PA0370', 0, 0, 0, 6, 0 ) … … 1107 1083 END SUBROUTINE nudge_init 1108 1084 1109 !------------------------------------------------------------------------------ !1085 !--------------------------------------------------------------------------------------------------! 1110 1086 ! Description: 1111 1087 ! ------------ 1112 1088 !> @todo Missing subroutine description. 1113 !------------------------------------------------------------------------------ !1089 !--------------------------------------------------------------------------------------------------! 1114 1090 SUBROUTINE calc_tnudge ( time ) 1115 1091 1116 1092 IMPLICIT NONE 1117 1093 1094 INTEGER(iwp) :: k !< 1095 INTEGER(iwp) :: nt !< 1118 1096 1119 1097 REAL(wp) :: dtm !< 1120 1098 REAL(wp) :: dtp !< 1121 1099 REAL(wp) :: time !< 1122 1123 INTEGER(iwp) :: k !<1124 INTEGER(iwp) :: nt !<1125 1100 1126 1101 nt = 1 … … 1141 1116 END SUBROUTINE calc_tnudge 1142 1117 1143 !------------------------------------------------------------------------------ !1118 !--------------------------------------------------------------------------------------------------! 1144 1119 ! Description: 1145 1120 ! ------------ 1146 1121 !> Call for all grid points 1147 !------------------------------------------------------------------------------ !1122 !--------------------------------------------------------------------------------------------------! 1148 1123 SUBROUTINE nudge ( time, prog_var ) 1149 1124 … … 1151 1126 1152 1127 CHARACTER (LEN=*) :: prog_var !< 1153 1154 REAL(wp) :: tmp_tend !<1155 REAL(wp) :: dtm !<1156 REAL(wp) :: dtp !<1157 REAL(wp) :: time !<1158 1128 1159 1129 INTEGER(iwp) :: i !< … … 1161 1131 INTEGER(iwp) :: k !< 1162 1132 INTEGER(iwp) :: nt !< 1133 1134 REAL(wp) :: dtm !< 1135 REAL(wp) :: dtp !< 1136 REAL(wp) :: time !< 1137 REAL(wp) :: tmp_tend !< 1163 1138 1164 1139 … … 1183 1158 DO k = nzb+1, nzt 1184 1159 1185 tmp_tend = - ( hom(k,1,1,0) - ( unudge(k,nt) * dtp + 1186 unudge(k,nt+1) * dtm ) )/ tmp_tnudge(k)1187 1188 tend(k,j,i) = tend(k,j,i) + tmp_tend *&1189 MERGE( 1.0_wp, 0.0_wp,&1190 BTEST( wall_flags_total_0(k,j,i), 1 ) )1191 1192 sums_ls_l(k,6) = sums_ls_l(k,6) + tmp_tend *&1193 weight_substep(intermediate_timestep_count)1160 tmp_tend = - ( hom(k,1,1,0) - ( unudge(k,nt) * dtp + unudge(k,nt+1) * dtm ) )& 1161 / tmp_tnudge(k) 1162 1163 tend(k,j,i) = tend(k,j,i) & 1164 + tmp_tend * MERGE( 1.0_wp, 0.0_wp, & 1165 BTEST( wall_flags_total_0(k,j,i), 1 ) ) 1166 1167 sums_ls_l(k,6) = sums_ls_l(k,6) & 1168 + tmp_tend * weight_substep(intermediate_timestep_count) 1194 1169 ENDDO 1195 1170 1196 1171 sums_ls_l(nzt+1,6) = sums_ls_l(nzt,6) 1197 1172 1198 1173 ENDDO 1199 1174 ENDDO … … 1206 1181 DO k = nzb+1, nzt 1207 1182 1208 tmp_tend = - ( hom(k,1,2,0) - ( vnudge(k,nt) * dtp + 1209 vnudge(k,nt+1) * dtm ) )/ tmp_tnudge(k)1210 1211 tend(k,j,i) = tend(k,j,i) + tmp_tend *&1212 MERGE( 1.0_wp, 0.0_wp,&1213 BTEST( wall_flags_total_0(k,j,i), 2 ) )1214 1215 sums_ls_l(k,7) = sums_ls_l(k,7) + tmp_tend *&1216 weight_substep(intermediate_timestep_count)1183 tmp_tend = - ( hom(k,1,2,0) - ( vnudge(k,nt) * dtp + vnudge(k,nt+1) * dtm ) )& 1184 / tmp_tnudge(k) 1185 1186 tend(k,j,i) = tend(k,j,i) & 1187 + tmp_tend * MERGE( 1.0_wp, 0.0_wp, & 1188 BTEST( wall_flags_total_0(k,j,i), 2 ) ) 1189 1190 sums_ls_l(k,7) = sums_ls_l(k,7) & 1191 + tmp_tend * weight_substep(intermediate_timestep_count) 1217 1192 ENDDO 1218 1193 1219 1194 sums_ls_l(nzt+1,7) = sums_ls_l(nzt,7) 1220 1195 … … 1229 1204 DO k = nzb+1, nzt 1230 1205 1231 tmp_tend = - ( hom(k,1,4,0) - ( ptnudge(k,nt) * dtp + &1206 tmp_tend = - ( hom(k,1,4,0) - ( ptnudge(k,nt) * dtp + & 1232 1207 ptnudge(k,nt+1) * dtm ) ) / tmp_tnudge(k) 1233 1208 1234 tend(k,j,i) = tend(k,j,i) + tmp_tend *&1235 MERGE( 1.0_wp, 0.0_wp,&1236 BTEST( wall_flags_total_0(k,j,i), 0 ) )1237 1238 sums_ls_l(k,4) = sums_ls_l(k,4) + tmp_tend *&1239 weight_substep(intermediate_timestep_count)1209 tend(k,j,i) = tend(k,j,i) & 1210 + tmp_tend * MERGE( 1.0_wp, 0.0_wp, & 1211 BTEST( wall_flags_total_0(k,j,i), 0 ) ) 1212 1213 sums_ls_l(k,4) = sums_ls_l(k,4) & 1214 + tmp_tend * weight_substep(intermediate_timestep_count) 1240 1215 ENDDO 1241 1216 … … 1252 1227 DO k = nzb+1, nzt 1253 1228 1254 tmp_tend = - ( hom(k,1,41,0) - ( qnudge(k,nt) * dtp + &1229 tmp_tend = - ( hom(k,1,41,0) - ( qnudge(k,nt) * dtp + & 1255 1230 qnudge(k,nt+1) * dtm ) ) / tmp_tnudge(k) 1256 1231 1257 tend(k,j,i) = tend(k,j,i) + tmp_tend *&1258 MERGE( 1.0_wp, 0.0_wp,&1259 BTEST( wall_flags_total_0(k,j,i), 0 ) )1260 1261 sums_ls_l(k,5) = sums_ls_l(k,5) + tmp_tend *&1262 weight_substep(intermediate_timestep_count)1232 tend(k,j,i) = tend(k,j,i) & 1233 + tmp_tend * MERGE( 1.0_wp, 0.0_wp, & 1234 BTEST( wall_flags_total_0(k,j,i), 0 ) ) 1235 1236 sums_ls_l(k,5) = sums_ls_l(k,5) & 1237 + tmp_tend * weight_substep(intermediate_timestep_count) 1263 1238 ENDDO 1264 1239 1265 1240 sums_ls_l(nzt+1,5) = sums_ls_l(nzt,5) 1266 1241 … … 1277 1252 1278 1253 1279 !------------------------------------------------------------------------------ !1254 !--------------------------------------------------------------------------------------------------! 1280 1255 ! Description: 1281 1256 ! ------------ 1282 1257 !> Call for grid point i,j 1283 !------------------------------------------------------------------------------ !1258 !--------------------------------------------------------------------------------------------------! 1284 1259 1285 1260 SUBROUTINE nudge_ij( i, j, time, prog_var ) … … 1289 1264 1290 1265 CHARACTER (LEN=*) :: prog_var !< 1291 1292 REAL(wp) :: tmp_tend !<1293 REAL(wp) :: dtm !<1294 REAL(wp) :: dtp !<1295 REAL(wp) :: time !<1296 1266 1297 1267 INTEGER(iwp) :: i !< … … 1299 1269 INTEGER(iwp) :: k !< 1300 1270 INTEGER(iwp) :: nt !< 1271 1272 REAL(wp) :: dtm !< 1273 REAL(wp) :: dtp !< 1274 REAL(wp) :: time !< 1275 REAL(wp) :: tmp_tend !< 1301 1276 1302 1277 … … 1318 1293 DO k = nzb+1, nzt 1319 1294 1320 tmp_tend = - ( hom(k,1,1,0) - ( unudge(k,nt) * dtp + 1321 unudge(k,nt+1) * dtm ) )/ tmp_tnudge(k)1322 1323 tend(k,j,i) = tend(k,j,i) + tmp_tend *&1324 MERGE( 1.0_wp, 0.0_wp,&1325 BTEST( wall_flags_total_0(k,j,i), 1 ) )1326 1327 sums_ls_l(k,6) = sums_ls_l(k,6) + tmp_tend&1328 * weight_substep(intermediate_timestep_count)1295 tmp_tend = - ( hom(k,1,1,0) - ( unudge(k,nt) * dtp + unudge(k,nt+1) * dtm ) ) & 1296 / tmp_tnudge(k) 1297 1298 tend(k,j,i) = tend(k,j,i) & 1299 + tmp_tend * MERGE( 1.0_wp, 0.0_wp, & 1300 BTEST( wall_flags_total_0(k,j,i), 1 ) ) 1301 1302 sums_ls_l(k,6) = sums_ls_l(k,6) & 1303 + tmp_tend * weight_substep(intermediate_timestep_count) 1329 1304 ENDDO 1330 1305 … … 1335 1310 DO k = nzb+1, nzt 1336 1311 1337 tmp_tend = - ( hom(k,1,2,0) - ( vnudge(k,nt) * dtp + 1338 vnudge(k,nt+1) * dtm ) )/ tmp_tnudge(k)1339 1340 tend(k,j,i) = tend(k,j,i) + tmp_tend *&1341 MERGE( 1.0_wp, 0.0_wp,&1342 BTEST( wall_flags_total_0(k,j,i), 2 ) )1343 1344 sums_ls_l(k,7) = sums_ls_l(k,7) + tmp_tend&1345 * weight_substep(intermediate_timestep_count)1312 tmp_tend = - ( hom(k,1,2,0) - ( vnudge(k,nt) * dtp + vnudge(k,nt+1) * dtm ) ) & 1313 / tmp_tnudge(k) 1314 1315 tend(k,j,i) = tend(k,j,i) & 1316 + tmp_tend * MERGE( 1.0_wp, 0.0_wp, & 1317 BTEST( wall_flags_total_0(k,j,i), 2 ) ) 1318 1319 sums_ls_l(k,7) = sums_ls_l(k,7) & 1320 + tmp_tend * weight_substep(intermediate_timestep_count) 1346 1321 ENDDO 1347 1322 … … 1352 1327 DO k = nzb+1, nzt 1353 1328 1354 tmp_tend = - ( hom(k,1,4,0) - ( ptnudge(k,nt) * dtp + 1355 ptnudge(k,nt+1) * dtm ) )/ tmp_tnudge(k)1356 1357 tend(k,j,i) = tend(k,j,i) + tmp_tend *&1358 MERGE( 1.0_wp, 0.0_wp,&1359 BTEST( wall_flags_total_0(k,j,i), 0 ) )1360 1361 sums_ls_l(k,4) = sums_ls_l(k,4) + tmp_tend&1362 * weight_substep(intermediate_timestep_count)1329 tmp_tend = - ( hom(k,1,4,0) - ( ptnudge(k,nt) * dtp + ptnudge(k,nt+1) * dtm ) ) & 1330 / tmp_tnudge(k) 1331 1332 tend(k,j,i) = tend(k,j,i) & 1333 + tmp_tend * MERGE( 1.0_wp, 0.0_wp, & 1334 BTEST( wall_flags_total_0(k,j,i), 0 ) ) 1335 1336 sums_ls_l(k,4) = sums_ls_l(k,4) & 1337 + tmp_tend * weight_substep(intermediate_timestep_count) 1363 1338 ENDDO 1364 1339 … … 1370 1345 DO k = nzb+1, nzt 1371 1346 1372 tmp_tend = - ( hom(k,1,41,0) - ( qnudge(k,nt) * dtp + 1373 qnudge(k,nt+1) * dtm ) )/ tmp_tnudge(k)1374 1375 tend(k,j,i) = tend(k,j,i) + tmp_tend *&1376 MERGE( 1.0_wp, 0.0_wp,&1377 BTEST( wall_flags_total_0(k,j,i), 0 ) )1378 1379 sums_ls_l(k,5) = sums_ls_l(k,5) + tmp_tend&1380 * weight_substep(intermediate_timestep_count)1347 tmp_tend = - ( hom(k,1,41,0) - ( qnudge(k,nt) * dtp + qnudge(k,nt+1) * dtm ) ) & 1348 / tmp_tnudge(k) 1349 1350 tend(k,j,i) = tend(k,j,i) & 1351 + tmp_tend * MERGE( 1.0_wp, 0.0_wp, & 1352 BTEST( wall_flags_total_0(k,j,i), 0 ) ) 1353 1354 sums_ls_l(k,5) = sums_ls_l(k,5) & 1355 + tmp_tend * weight_substep(intermediate_timestep_count) 1381 1356 ENDDO 1382 1357 … … 1393 1368 1394 1369 1395 !------------------------------------------------------------------------------ !1370 !--------------------------------------------------------------------------------------------------! 1396 1371 ! Description: 1397 1372 ! ------------ 1398 1373 !> @todo Missing subroutine description. 1399 !------------------------------------------------------------------------------ !1374 !--------------------------------------------------------------------------------------------------! 1400 1375 SUBROUTINE nudge_ref ( time ) 1401 1376 … … 1408 1383 1409 1384 ! 1410 !-- Interpolation in time of NUDGING_DATA for pt_init and q_init. This is 1411 !-- needed for correct upper boundary conditions for pt and q and in case that1412 ! large scale subsidence as well as scalar Rayleigh-damping are used1385 !-- Interpolation in time of NUDGING_DATA for pt_init and q_init. This is needed for correct 1386 !-- upper boundary conditions for pt and q and in case that large scale subsidence as well as 1387 !-- scalar Rayleigh-damping are used. 1413 1388 nt = 1 1414 1389 DO WHILE ( time > time_vert(nt) )
Note: See TracChangeset
for help on using the changeset viewer.
