Changeset 4583 for palm/trunk/SOURCE/diffusion_s.f90
- Timestamp:
- Jun 29, 2020 12:36:47 PM (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/diffusion_s.f90
r4360 r4583 1 1 !> @file diffusion_s.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$ 27 ! Introduction of wall_flags_total_0, which currently sets bits based on static 28 ! topography information used in wall_flags_static_0 29 ! 26 ! file re-formatted to follow the PALM coding standard 27 ! 28 ! 4360 2020-01-07 11:25:50Z suehring 29 ! Introduction of wall_flags_total_0, which currently sets bits based on static topography 30 ! information used in wall_flags_static_0 31 ! 30 32 ! 4329 2019-12-10 15:46:36Z motisi 31 33 ! Renamed wall_flags_0 to wall_flags_total_0 32 ! 34 ! 33 35 ! 4182 2019-08-22 15:20:23Z scharf 34 36 ! Corrected "Former revisions" section 35 ! 37 ! 36 38 ! 3927 2019-04-23 13:24:29Z raasch 37 39 ! pointer attribute removed from scalar 3d-array for performance reasons 38 ! 40 ! 39 41 ! 3761 2019-02-25 15:31:42Z raasch 40 42 ! unused variables removed 41 ! 43 ! 42 44 ! 3655 2019-01-07 16:51:22Z knoop 43 45 ! nopointer option removed … … 50 52 ! ------------ 51 53 !> Diffusion term of scalar quantities (temperature and water content) 52 !------------------------------------------------------------------------------ !54 !--------------------------------------------------------------------------------------------------! 53 55 MODULE diffusion_s_mod 54 56 55 57 56 58 PRIVATE … … 65 67 66 68 67 !------------------------------------------------------------------------------ !69 !--------------------------------------------------------------------------------------------------! 68 70 ! Description: 69 71 ! ------------ 70 72 !> Call for all grid points 71 !------------------------------------------------------------------------------ !72 SUBROUTINE diffusion_s( s, s_flux_def_h_up, s_flux_def_h_down, &73 s_flux_t, &74 s_flux_lsm_h_up, s_flux_usm_h_up, &75 s_flux_def_v_north, s_flux_def_v_south, &76 s_flux_def_v_east, s_flux_def_v_west, &77 s_flux_lsm_v_north, s_flux_lsm_v_south, &78 s_flux_lsm_v_east, s_flux_lsm_v_west, &79 s_flux_usm_v_north, s_flux_usm_v_south, &73 !--------------------------------------------------------------------------------------------------! 74 SUBROUTINE diffusion_s( s, s_flux_def_h_up, s_flux_def_h_down, & 75 s_flux_t, & 76 s_flux_lsm_h_up, s_flux_usm_h_up, & 77 s_flux_def_v_north, s_flux_def_v_south, & 78 s_flux_def_v_east, s_flux_def_v_west, & 79 s_flux_lsm_v_north, s_flux_lsm_v_south, & 80 s_flux_lsm_v_east, s_flux_lsm_v_west, & 81 s_flux_usm_v_north, s_flux_usm_v_south, & 80 82 s_flux_usm_v_east, s_flux_usm_v_west ) 81 83 82 USE arrays_3d, &84 USE arrays_3d, & 83 85 ONLY: ddzu, ddzw, kh, tend, drho_air, rho_air_zw 84 85 USE control_parameters, &86 87 USE control_parameters, & 86 88 ONLY: use_surface_fluxes, use_top_fluxes 87 88 USE grid_variables, &89 90 USE grid_variables, & 89 91 ONLY: ddx, ddx2, ddy, ddy2 90 91 USE indices, & 92 ONLY: nxl, nxlg, nxr, nxrg, nyn, nyng, nys, nysg, nzb, nzt, & 93 wall_flags_total_0 94 92 93 USE indices, & 94 ONLY: nxl, nxlg, nxr, nxrg, nyn, nyng, nys, nysg, nzb, nzt, wall_flags_total_0 95 95 96 USE kinds 96 97 97 USE surface_mod, & 98 ONLY : surf_def_h, surf_def_v, surf_lsm_h, surf_lsm_v, surf_usm_h, & 99 surf_usm_v 98 USE surface_mod, & 99 ONLY : surf_def_h, surf_def_v, surf_lsm_h, surf_lsm_v, surf_usm_h, surf_usm_v 100 100 101 101 IMPLICIT NONE … … 109 109 110 110 REAL(wp) :: flag !< flag to mask topography grid points 111 REAL(wp) :: mask_bottom !< flag to mask vertical upward-facing surface 112 REAL(wp) :: mask_east !< flag to mask vertical surface east of the grid point 111 REAL(wp) :: mask_bottom !< flag to mask vertical upward-facing surface 112 REAL(wp) :: mask_east !< flag to mask vertical surface east of the grid point 113 113 REAL(wp) :: mask_north !< flag to mask vertical surface north of the grid point 114 REAL(wp) :: mask_south !< flag to mask vertical surface south of the grid point 114 REAL(wp) :: mask_south !< flag to mask vertical surface south of the grid point 115 REAL(wp) :: mask_top !< flag to mask vertical downward-facing surface 115 116 REAL(wp) :: mask_west !< flag to mask vertical surface west of the grid point 116 REAL(wp) :: mask_top !< flag to mask vertical downward-facing surface 117 117 118 REAL(wp), DIMENSION(1:surf_def_h(1)%ns) :: s_flux_def_h_down !< flux at horizontal donwward-facing default-type surfaces 119 REAL(wp), DIMENSION(1:surf_def_h(0)%ns) :: s_flux_def_h_up !< flux at horizontal upward-facing default-type surfaces 120 REAL(wp), DIMENSION(1:surf_def_v(2)%ns) :: s_flux_def_v_east !< flux at east-facing vertical default-type surfaces 118 121 REAL(wp), DIMENSION(1:surf_def_v(0)%ns) :: s_flux_def_v_north !< flux at north-facing vertical default-type surfaces 119 122 REAL(wp), DIMENSION(1:surf_def_v(1)%ns) :: s_flux_def_v_south !< flux at south-facing vertical default-type surfaces 120 REAL(wp), DIMENSION(1:surf_def_v(2)%ns) :: s_flux_def_v_east !< flux at east-facing vertical default-type surfaces121 123 REAL(wp), DIMENSION(1:surf_def_v(3)%ns) :: s_flux_def_v_west !< flux at west-facing vertical default-type surfaces 122 REAL(wp), DIMENSION(1:surf_def_h(0)%ns) :: s_flux_def_h_up !< flux at horizontal upward-facing default-type surfaces123 REAL(wp), DIMENSION(1:surf_def_h(1)%ns) :: s_flux_def_h_down !< flux at horizontal donwward-facing default-type surfaces124 124 REAL(wp), DIMENSION(1:surf_lsm_h%ns) :: s_flux_lsm_h_up !< flux at horizontal upward-facing natural-type surfaces 125 REAL(wp), DIMENSION(1:surf_lsm_v(2)%ns) :: s_flux_lsm_v_east !< flux at east-facing vertical natural-type surfaces 125 126 REAL(wp), DIMENSION(1:surf_lsm_v(0)%ns) :: s_flux_lsm_v_north !< flux at north-facing vertical natural-type surfaces 126 127 REAL(wp), DIMENSION(1:surf_lsm_v(1)%ns) :: s_flux_lsm_v_south !< flux at south-facing vertical natural-type surfaces 127 REAL(wp), DIMENSION(1:surf_lsm_v(2)%ns) :: s_flux_lsm_v_east !< flux at east-facing vertical natural-type surfaces128 128 REAL(wp), DIMENSION(1:surf_lsm_v(3)%ns) :: s_flux_lsm_v_west !< flux at west-facing vertical natural-type surfaces 129 129 REAL(wp), DIMENSION(1:surf_usm_h%ns) :: s_flux_usm_h_up !< flux at horizontal upward-facing urban-type surfaces 130 REAL(wp), DIMENSION(1:surf_usm_v(2)%ns) :: s_flux_usm_v_east !< flux at east-facing vertical urban-type surfaces 130 131 REAL(wp), DIMENSION(1:surf_usm_v(0)%ns) :: s_flux_usm_v_north !< flux at north-facing vertical urban-type surfaces 131 132 REAL(wp), DIMENSION(1:surf_usm_v(1)%ns) :: s_flux_usm_v_south !< flux at south-facing vertical urban-type surfaces 132 REAL(wp), DIMENSION(1:surf_usm_v(2)%ns) :: s_flux_usm_v_east !< flux at east-facing vertical urban-type surfaces133 133 REAL(wp), DIMENSION(1:surf_usm_v(3)%ns) :: s_flux_usm_v_west !< flux at west-facing vertical urban-type surfaces 134 134 REAL(wp), DIMENSION(1:surf_def_h(2)%ns) :: s_flux_t !< flux at model top … … 164 164 ! 165 165 !-- Predetermine flag to mask topography and wall-bounded grid points 166 flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) ) 167 ! 168 !-- Predetermine flag to mask wall-bounded grid points, equivalent to 169 !-- former s_outerarray166 flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) ) 167 ! 168 !-- Predetermine flag to mask wall-bounded grid points, equivalent to former s_outer 169 !-- array 170 170 mask_west = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i-1), 0 ) ) 171 171 mask_east = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i+1), 0 ) ) … … 173 173 mask_north = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j+1,i), 0 ) ) 174 174 175 tend(k,j,i) = tend(k,j,i) & 176 + 0.5_wp * ( & 177 mask_east * ( kh(k,j,i) + kh(k,j,i+1) ) & 178 * ( s(k,j,i+1) - s(k,j,i) ) & 179 - mask_west * ( kh(k,j,i) + kh(k,j,i-1) ) & 180 * ( s(k,j,i) - s(k,j,i-1) ) & 181 ) * ddx2 * flag & 182 + 0.5_wp * ( & 183 mask_north * ( kh(k,j,i) + kh(k,j+1,i) ) & 184 * ( s(k,j+1,i) - s(k,j,i) ) & 185 - mask_south * ( kh(k,j,i) + kh(k,j-1,i) ) & 186 * ( s(k,j,i) - s(k,j-1,i) ) & 187 ) * ddy2 * flag 188 ENDDO 189 190 ! 191 !-- Apply prescribed horizontal wall heatflux where necessary. First, 192 !-- determine start and end index for respective (j,i)-index. Please 193 !-- note, in the flat case following loop will not be entered, as 194 !-- surf_s=1 and surf_e=0. Furtermore, note, no vertical natural surfaces 195 !-- so far. 196 !-- First, for default-type surfaces 175 tend(k,j,i) = tend(k,j,i) & 176 + 0.5_wp * ( & 177 mask_east * ( kh(k,j,i) + kh(k,j,i+1) ) & 178 * ( s(k,j,i+1) - s(k,j,i) ) & 179 - mask_west * ( kh(k,j,i) + kh(k,j,i-1) ) & 180 * ( s(k,j,i) - s(k,j,i-1) ) & 181 ) * ddx2 * flag & 182 + 0.5_wp * ( & 183 mask_north * ( kh(k,j,i) + kh(k,j+1,i) ) & 184 * ( s(k,j+1,i) - s(k,j,i) ) & 185 - mask_south * ( kh(k,j,i) + kh(k,j-1,i) ) & 186 * ( s(k,j,i) - s(k,j-1,i) ) & 187 ) * ddy2 * flag 188 ENDDO 189 190 ! 191 !-- Apply prescribed horizontal wall heatflux where necessary. First, determine start and 192 !-- end index for respective (j,i)-index. Please note, in the flat case following loop will 193 !-- not be entered, as surf_s=1 and surf_e=0. Furtermore, note, no vertical natural 194 !-- surfaces so far. 195 !-- First, for default-type surfaces. 197 196 !-- North-facing vertical default-type surfaces 198 197 surf_s = surf_def_v(0)%start_index(j,i) … … 294 293 295 294 ! 296 !-- Compute vertical diffusion. In case that surface fluxes have been 297 !-- prescribed or computed at bottom and/or top, index k starts/ends at 298 !-- nzb+2 or nzt-1, respectively. Model top is also mask if top flux 299 !-- is given. 295 !-- Compute vertical diffusion. In case that surface fluxes have been prescribed or 296 !-- computed at bottom and/or top, index k starts/ends at nzb+2 or nzt-1, respectively. 297 !-- Model top is also mask if top flux is given. 300 298 DO k = nzb+1, nzt 301 299 ! 302 !-- Determine flags to mask topography below and above. Flag 0 is 303 !-- used to mask topography in general, and flag 8 implies 304 !-- information about use_surface_fluxes. Flag 9 is used to control 305 !-- flux at model top. 306 mask_bottom = MERGE( 1.0_wp, 0.0_wp, & 307 BTEST( wall_flags_total_0(k-1,j,i), 8 ) ) 308 mask_top = MERGE( 1.0_wp, 0.0_wp, & 309 BTEST( wall_flags_total_0(k+1,j,i), 8 ) ) * & 310 MERGE( 1.0_wp, 0.0_wp, & 311 BTEST( wall_flags_total_0(k+1,j,i), 9 ) ) 312 flag = MERGE( 1.0_wp, 0.0_wp, & 313 BTEST( wall_flags_total_0(k,j,i), 0 ) ) 314 315 tend(k,j,i) = tend(k,j,i) & 316 + 0.5_wp * ( & 317 ( kh(k,j,i) + kh(k+1,j,i) ) * & 318 ( s(k+1,j,i)-s(k,j,i) ) * ddzu(k+1) & 319 * rho_air_zw(k) & 320 * mask_top & 321 - ( kh(k,j,i) + kh(k-1,j,i) ) * & 322 ( s(k,j,i)-s(k-1,j,i) ) * ddzu(k) & 323 * rho_air_zw(k-1) & 324 * mask_bottom & 325 ) * ddzw(k) * drho_air(k) & 300 !-- Determine flags to mask topography below and above. Flag 0 is used to mask 301 !-- topography in general, and flag 8 implies information about use_surface_fluxes. 302 !-- Flag 9 is used to control flux at model top. 303 mask_bottom = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k-1,j,i), 8 ) ) 304 mask_top = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k+1,j,i), 8 ) ) * & 305 MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k+1,j,i), 9 ) ) 306 flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) ) 307 308 tend(k,j,i) = tend(k,j,i) & 309 + 0.5_wp * ( & 310 ( kh(k,j,i) + kh(k+1,j,i) ) * & 311 ( s(k+1,j,i)-s(k,j,i) ) * ddzu(k+1) & 312 * rho_air_zw(k) & 313 * mask_top & 314 - ( kh(k,j,i) + kh(k-1,j,i) ) * & 315 ( s(k,j,i)-s(k-1,j,i) ) * ddzu(k) & 316 * rho_air_zw(k-1) & 317 * mask_bottom & 318 ) * ddzw(k) * drho_air(k) & 326 319 * flag 327 320 ENDDO … … 331 324 IF ( use_surface_fluxes ) THEN 332 325 ! 333 !-- Default-type surfaces, upward-facing 326 !-- Default-type surfaces, upward-facing 334 327 surf_s = surf_def_h(0)%start_index(j,i) 335 328 surf_e = surf_def_h(0)%end_index(j,i) … … 337 330 338 331 k = surf_def_h(0)%k(m) 339 tend(k,j,i) = tend(k,j,i) + s_flux_def_h_up(m) & 340 * ddzw(k) * drho_air(k) 332 tend(k,j,i) = tend(k,j,i) + s_flux_def_h_up(m) * ddzw(k) * drho_air(k) 341 333 342 334 ENDDO 343 335 ! 344 !-- Default-type surfaces, downward-facing 336 !-- Default-type surfaces, downward-facing 345 337 surf_s = surf_def_h(1)%start_index(j,i) 346 338 surf_e = surf_def_h(1)%end_index(j,i) … … 348 340 349 341 k = surf_def_h(1)%k(m) 350 tend(k,j,i) = tend(k,j,i) + s_flux_def_h_down(m) & 351 * ddzw(k) * drho_air(k) 342 tend(k,j,i) = tend(k,j,i) + s_flux_def_h_down(m) * ddzw(k) * drho_air(k) 352 343 353 344 ENDDO 354 345 ! 355 !-- Natural-type surfaces, upward-facing 346 !-- Natural-type surfaces, upward-facing 356 347 surf_s = surf_lsm_h%start_index(j,i) 357 348 surf_e = surf_lsm_h%end_index(j,i) … … 359 350 360 351 k = surf_lsm_h%k(m) 361 tend(k,j,i) = tend(k,j,i) + s_flux_lsm_h_up(m) & 362 * ddzw(k) * drho_air(k) 352 tend(k,j,i) = tend(k,j,i) + s_flux_lsm_h_up(m) * ddzw(k) * drho_air(k) 363 353 364 354 ENDDO 365 355 ! 366 !-- Urban-type surfaces, upward-facing 356 !-- Urban-type surfaces, upward-facing 367 357 surf_s = surf_usm_h%start_index(j,i) 368 358 surf_e = surf_usm_h%end_index(j,i) … … 370 360 371 361 k = surf_usm_h%k(m) 372 tend(k,j,i) = tend(k,j,i) + s_flux_usm_h_up(m) & 373 * ddzw(k) * drho_air(k) 362 tend(k,j,i) = tend(k,j,i) + s_flux_usm_h_up(m) * ddzw(k) * drho_air(k) 374 363 375 364 ENDDO … … 384 373 385 374 k = surf_def_h(2)%k(m) 386 tend(k,j,i) = tend(k,j,i) & 387 + ( - s_flux_t(m) ) * ddzw(k) * drho_air(k) 375 tend(k,j,i) = tend(k,j,i) + ( - s_flux_t(m) ) * ddzw(k) * drho_air(k) 388 376 ENDDO 389 377 ENDIF … … 394 382 END SUBROUTINE diffusion_s 395 383 396 !------------------------------------------------------------------------------ !384 !--------------------------------------------------------------------------------------------------! 397 385 ! Description: 398 386 ! ------------ 399 387 !> Call for grid point i,j 400 !------------------------------------------------------------------------------ !401 SUBROUTINE diffusion_s_ij( i, j, s, &402 s_flux_def_h_up, s_flux_def_h_down, &403 s_flux_t, &404 s_flux_lsm_h_up, s_flux_usm_h_up, &405 s_flux_def_v_north, s_flux_def_v_south, &406 s_flux_def_v_east, s_flux_def_v_west, &407 s_flux_lsm_v_north, s_flux_lsm_v_south, &408 s_flux_lsm_v_east, s_flux_lsm_v_west, &409 s_flux_usm_v_north, s_flux_usm_v_south, &410 s_flux_usm_v_east, s_flux_usm_v_west ) 411 412 USE arrays_3d, &388 !--------------------------------------------------------------------------------------------------! 389 SUBROUTINE diffusion_s_ij( i, j, s, & 390 s_flux_def_h_up, s_flux_def_h_down, & 391 s_flux_t, & 392 s_flux_lsm_h_up, s_flux_usm_h_up, & 393 s_flux_def_v_north, s_flux_def_v_south, & 394 s_flux_def_v_east, s_flux_def_v_west, & 395 s_flux_lsm_v_north, s_flux_lsm_v_south, & 396 s_flux_lsm_v_east, s_flux_lsm_v_west, & 397 s_flux_usm_v_north, s_flux_usm_v_south, & 398 s_flux_usm_v_east, s_flux_usm_v_west ) 399 400 USE arrays_3d, & 413 401 ONLY: ddzu, ddzw, kh, tend, drho_air, rho_air_zw 414 415 USE control_parameters, &402 403 USE control_parameters, & 416 404 ONLY: use_surface_fluxes, use_top_fluxes 417 418 USE grid_variables, &405 406 USE grid_variables, & 419 407 ONLY: ddx, ddx2, ddy, ddy2 420 421 USE indices, &408 409 USE indices, & 422 410 ONLY: nxlg, nxrg, nyng, nysg, nzb, nzt, wall_flags_total_0 423 411 424 412 USE kinds 425 413 426 USE surface_mod, & 427 ONLY : surf_def_h, surf_def_v, surf_lsm_h, surf_lsm_v, surf_usm_h, & 428 surf_usm_v 414 USE surface_mod, & 415 ONLY : surf_def_h, surf_def_v, surf_lsm_h, surf_lsm_v, surf_usm_h, surf_usm_v 429 416 430 417 IMPLICIT NONE … … 438 425 439 426 REAL(wp) :: flag !< flag to mask topography grid points 440 REAL(wp) :: mask_bottom !< flag to mask vertical upward-facing surface 441 REAL(wp) :: mask_east !< flag to mask vertical surface east of the grid point 427 REAL(wp) :: mask_bottom !< flag to mask vertical upward-facing surface 428 REAL(wp) :: mask_east !< flag to mask vertical surface east of the grid point 442 429 REAL(wp) :: mask_north !< flag to mask vertical surface north of the grid point 443 REAL(wp) :: mask_south !< flag to mask vertical surface south of the grid point 430 REAL(wp) :: mask_south !< flag to mask vertical surface south of the grid point 431 REAL(wp) :: mask_top !< flag to mask vertical downward-facing surface 444 432 REAL(wp) :: mask_west !< flag to mask vertical surface west of the grid point 445 REAL(wp) :: mask_top !< flag to mask vertical downward-facing surface 446 433 434 REAL(wp), DIMENSION(1:surf_def_h(1)%ns) :: s_flux_def_h_down !< flux at horizontal donwward-facing default-type surfaces 435 REAL(wp), DIMENSION(1:surf_def_h(0)%ns) :: s_flux_def_h_up !< flux at horizontal upward-facing default-type surfaces 436 REAL(wp), DIMENSION(1:surf_def_v(2)%ns) :: s_flux_def_v_east !< flux at east-facing vertical default-type surfaces 447 437 REAL(wp), DIMENSION(1:surf_def_v(0)%ns) :: s_flux_def_v_north !< flux at north-facing vertical default-type surfaces 448 438 REAL(wp), DIMENSION(1:surf_def_v(1)%ns) :: s_flux_def_v_south !< flux at south-facing vertical default-type surfaces 449 REAL(wp), DIMENSION(1:surf_def_v(2)%ns) :: s_flux_def_v_east !< flux at east-facing vertical default-type surfaces450 439 REAL(wp), DIMENSION(1:surf_def_v(3)%ns) :: s_flux_def_v_west !< flux at west-facing vertical default-type surfaces 451 REAL(wp), DIMENSION(1:surf_def_h(0)%ns) :: s_flux_def_h_up !< flux at horizontal upward-facing default-type surfaces452 REAL(wp), DIMENSION(1:surf_def_h(1)%ns) :: s_flux_def_h_down !< flux at horizontal donwward-facing default-type surfaces453 440 REAL(wp), DIMENSION(1:surf_lsm_h%ns) :: s_flux_lsm_h_up !< flux at horizontal upward-facing natural-type surfaces 441 REAL(wp), DIMENSION(1:surf_lsm_v(2)%ns) :: s_flux_lsm_v_east !< flux at east-facing vertical urban-type surfaces 454 442 REAL(wp), DIMENSION(1:surf_lsm_v(0)%ns) :: s_flux_lsm_v_north !< flux at north-facing vertical urban-type surfaces 455 443 REAL(wp), DIMENSION(1:surf_lsm_v(1)%ns) :: s_flux_lsm_v_south !< flux at south-facing vertical urban-type surfaces 456 REAL(wp), DIMENSION(1:surf_lsm_v(2)%ns) :: s_flux_lsm_v_east !< flux at east-facing vertical urban-type surfaces457 444 REAL(wp), DIMENSION(1:surf_lsm_v(3)%ns) :: s_flux_lsm_v_west !< flux at west-facing vertical urban-type surfaces 458 445 REAL(wp), DIMENSION(1:surf_usm_h%ns) :: s_flux_usm_h_up !< flux at horizontal upward-facing urban-type surfaces 446 REAL(wp), DIMENSION(1:surf_usm_v(2)%ns) :: s_flux_usm_v_east !< flux at east-facing vertical urban-type surfaces 459 447 REAL(wp), DIMENSION(1:surf_usm_v(0)%ns) :: s_flux_usm_v_north !< flux at north-facing vertical urban-type surfaces 460 448 REAL(wp), DIMENSION(1:surf_usm_v(1)%ns) :: s_flux_usm_v_south !< flux at south-facing vertical urban-type surfaces 461 REAL(wp), DIMENSION(1:surf_usm_v(2)%ns) :: s_flux_usm_v_east !< flux at east-facing vertical urban-type surfaces462 449 REAL(wp), DIMENSION(1:surf_usm_v(3)%ns) :: s_flux_usm_v_west !< flux at west-facing vertical urban-type surfaces 463 450 REAL(wp), DIMENSION(1:surf_def_h(2)%ns) :: s_flux_t !< flux at model top … … 465 452 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) :: s !< treated scalar 466 453 454 467 455 ! 468 456 !-- Compute horizontal diffusion … … 470 458 ! 471 459 !-- Predetermine flag to mask topography and wall-bounded grid points 472 flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) ) 473 ! 474 !-- Predetermine flag to mask wall-bounded grid points, equivalent to 475 !-- former s_outer array 460 flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) ) 461 ! 462 !-- Predetermine flag to mask wall-bounded grid points, equivalent to former s_outer array 476 463 mask_west = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i-1), 0 ) ) 477 464 mask_east = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i+1), 0 ) ) … … 479 466 mask_north = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j+1,i), 0 ) ) 480 467 ! 481 !-- Finally, determine flag to mask both topography itself as well 482 !-- as wall-bounded gridpoints, which will be treated further below483 484 tend(k,j,i) = tend(k,j,i) &485 + 0.5_wp * ( &486 mask_east * ( kh(k,j,i) + kh(k,j,i+1) ) &487 * ( s(k,j,i+1) - s(k,j,i) ) &488 - mask_west * ( kh(k,j,i) + kh(k,j,i-1) ) &489 * ( s(k,j,i) - s(k,j,i-1) ) &490 ) * ddx2 * flag &491 + 0.5_wp * ( &492 mask_north * ( kh(k,j,i) + kh(k,j+1,i) ) &493 * ( s(k,j+1,i) - s(k,j,i) ) &494 - mask_south * ( kh(k,j,i) + kh(k,j-1,i) ) &495 * ( s(k,j,i) - s(k,j-1,i) ) &468 !-- Finally, determine flag to mask both topography itself as well as wall-bounded grid 469 !-- points, which will be treated further below 470 471 tend(k,j,i) = tend(k,j,i) & 472 + 0.5_wp * ( & 473 mask_east * ( kh(k,j,i) + kh(k,j,i+1) ) & 474 * ( s(k,j,i+1) - s(k,j,i) ) & 475 - mask_west * ( kh(k,j,i) + kh(k,j,i-1) ) & 476 * ( s(k,j,i) - s(k,j,i-1) ) & 477 ) * ddx2 * flag & 478 + 0.5_wp * ( & 479 mask_north * ( kh(k,j,i) + kh(k,j+1,i) ) & 480 * ( s(k,j+1,i) - s(k,j,i) ) & 481 - mask_south * ( kh(k,j,i) + kh(k,j-1,i) ) & 482 * ( s(k,j,i) - s(k,j-1,i) ) & 496 483 ) * ddy2 * flag 497 484 ENDDO 498 485 499 486 ! 500 !-- Apply prescribed horizontal wall heatflux where necessary. First, 501 !-- determine start and end index for respective (j,i)-index. Please 502 !-- note, in the flat case following loops will not be entered, as 503 !-- surf_s=1 and surf_e=0. Furtermore, note, no vertical natural surfaces 504 !-- so far. 505 !-- First, for default-type surfaces 487 !-- Apply prescribed horizontal wall heatflux where necessary. First, determine start and end 488 !-- index for respective (j,i)-index. Please note, in the flat case following loops will not be 489 !-- entered, as surf_s=1 and surf_e=0. Furtermore, note, no vertical natural surfaces so far. 490 !-- First, for default-type surfaces. 506 491 !-- North-facing vertical default-type surfaces 507 492 surf_s = surf_def_v(0)%start_index(j,i) … … 604 589 605 590 ! 606 !-- Compute vertical diffusion. In case that surface fluxes have been 607 !-- prescribed or computed at bottom and/or top, index k starts/ends at 608 !-- nzb+2 or nzt-1, respectively. Model top is also mask if top flux 609 !-- is given. 591 !-- Compute vertical diffusion. In case that surface fluxes have been prescribed or computed at 592 !-- bottom and/or top, index k starts/ends at nzb+2 or nzt-1, respectively. Model top is also 593 !-- mask if top flux is given. 610 594 DO k = nzb+1, nzt 611 595 ! 612 !-- Determine flags to mask topography below and above. Flag 0 is 613 !-- used to mask topography in general, and flag 8 implies 614 !-- information about use_surface_fluxes. Flag 9 is used to control 615 !-- flux at model top. 616 mask_bottom = MERGE( 1.0_wp, 0.0_wp, & 617 BTEST( wall_flags_total_0(k-1,j,i), 8 ) ) 618 mask_top = MERGE( 1.0_wp, 0.0_wp, & 619 BTEST( wall_flags_total_0(k+1,j,i), 8 ) ) * & 620 MERGE( 1.0_wp, 0.0_wp, & 621 BTEST( wall_flags_total_0(k+1,j,i), 9 ) ) 622 flag = MERGE( 1.0_wp, 0.0_wp, & 623 BTEST( wall_flags_total_0(k,j,i), 0 ) ) 624 625 tend(k,j,i) = tend(k,j,i) & 626 + 0.5_wp * ( & 627 ( kh(k,j,i) + kh(k+1,j,i) ) * & 628 ( s(k+1,j,i)-s(k,j,i) ) * ddzu(k+1) & 629 * rho_air_zw(k) & 630 * mask_top & 631 - ( kh(k,j,i) + kh(k-1,j,i) ) * & 632 ( s(k,j,i)-s(k-1,j,i) ) * ddzu(k) & 633 * rho_air_zw(k-1) & 634 * mask_bottom & 635 ) * ddzw(k) * drho_air(k) & 596 !-- Determine flags to mask topography below and above. Flag 0 is used to mask topography in 597 !-- general, and flag 8 implies information about use_surface_fluxes. Flag 9 is used to 598 !-- control flux at model top. 599 mask_bottom = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k-1,j,i), 8 ) ) 600 mask_top = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k+1,j,i), 8 ) ) * & 601 MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k+1,j,i), 9 ) ) 602 flag = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) ) 603 604 tend(k,j,i) = tend(k,j,i) & 605 + 0.5_wp * ( & 606 ( kh(k,j,i) + kh(k+1,j,i) ) * & 607 ( s(k+1,j,i)-s(k,j,i) ) * ddzu(k+1) & 608 * rho_air_zw(k) & 609 * mask_top & 610 - ( kh(k,j,i) + kh(k-1,j,i) ) * & 611 ( s(k,j,i)-s(k-1,j,i) ) * ddzu(k) & 612 * rho_air_zw(k-1) & 613 * mask_bottom & 614 ) * ddzw(k) * drho_air(k) & 636 615 * flag 637 616 ENDDO … … 649 628 k = surf_def_h(0)%k(m) 650 629 651 tend(k,j,i) = tend(k,j,i) + s_flux_def_h_up(m) & 652 * ddzw(k) * drho_air(k) 630 tend(k,j,i) = tend(k,j,i) + s_flux_def_h_up(m) * ddzw(k) * drho_air(k) 653 631 ENDDO 654 632 ! … … 660 638 k = surf_def_h(1)%k(m) 661 639 662 tend(k,j,i) = tend(k,j,i) + s_flux_def_h_down(m) & 663 * ddzw(k) * drho_air(k) 640 tend(k,j,i) = tend(k,j,i) + s_flux_def_h_down(m) * ddzw(k) * drho_air(k) 664 641 ENDDO 665 642 ! … … 670 647 k = surf_lsm_h%k(m) 671 648 672 tend(k,j,i) = tend(k,j,i) + s_flux_lsm_h_up(m) & 673 * ddzw(k) * drho_air(k) 649 tend(k,j,i) = tend(k,j,i) + s_flux_lsm_h_up(m) * ddzw(k) * drho_air(k) 674 650 ENDDO 675 651 ! … … 680 656 k = surf_usm_h%k(m) 681 657 682 tend(k,j,i) = tend(k,j,i) + s_flux_usm_h_up(m) & 683 * ddzw(k) * drho_air(k) 658 tend(k,j,i) = tend(k,j,i) + s_flux_usm_h_up(m) * ddzw(k) * drho_air(k) 684 659 ENDDO 685 660 ENDIF … … 692 667 693 668 k = surf_def_h(2)%k(m) 694 tend(k,j,i) = tend(k,j,i) & 695 + ( - s_flux_t(m) ) * ddzw(k) * drho_air(k) 669 tend(k,j,i) = tend(k,j,i) + ( - s_flux_t(m) ) * ddzw(k) * drho_air(k) 696 670 ENDDO 697 671 ENDIF
Note: See TracChangeset
for help on using the changeset viewer.