Changeset 4649 for palm/trunk/SOURCE/prognostic_equations.f90
- Timestamp:
- Aug 25, 2020 12:11:17 PM (4 years ago)
- File:
-
- 1 edited
-
palm/trunk/SOURCE/prognostic_equations.f90 (modified) (59 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/prognostic_equations.f90
r4370 r4649 1 1 !> @file prognostic_equations.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 !--------------------------------------------------------------------------------------------------! 18 ! 19 19 ! 20 20 ! Current revisions: 21 ! ----------------- -22 ! 23 ! 21 ! ----------------- 22 ! 23 ! 24 24 ! Former revisions: 25 25 ! ----------------- 26 26 ! $Id$ 27 ! vector directives added to force vectorization on Intel19 compiler 28 ! 27 ! File re-formatted to follow the PALM coding standard 28 ! 29 ! 30 ! 4370 2020-01-10 14:00:44Z raasch 31 ! Vector directives added to force vectorization on Intel19 compiler 32 ! 29 33 ! 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 ! 34 ! Introduction of wall_flags_total_0, which currently sets bits based on static topography 35 ! information used in wall_flags_static_0 36 ! 33 37 ! 4329 2019-12-10 15:46:36Z motisi 34 38 ! Renamed wall_flags_0 to wall_flags_static_0 35 ! 39 ! 36 40 ! 4182 2019-08-22 15:20:23Z scharf 37 41 ! Corrected "Former revisions" section 38 ! 42 ! 39 43 ! 4110 2019-07-22 17:05:21Z suehring 40 ! pass integer flag array to WS scalar advection routine which is now necessary41 ! as the flags may differ for scalars, e.g. pt can be cyclic while chemical42 ! species may be non-cyclic. Further, pass boundary flags.43 ! 44 ! Pass integer flag array to WS scalar advection routine which is now necessary as the flags may 45 ! differ for scalars, e.g. pt can be cyclic while chemical species may be non-cyclic. Further, 46 ! pass boundary flags. 47 ! 44 48 ! 4109 2019-07-22 17:00:34Z suehring 45 ! Application of monotonic flux limiter for the vertical scalar advection 46 ! up to the topography top (only for the cache-optimized version at the 47 ! moment). Please note, at the moment the limiter is only applied for passive 48 ! scalars. 49 ! 49 ! Application of monotonic flux limiter for the vertical scalar advection up to the topography top 50 ! (only for the cache-optimized version at the moment). Please note, at the moment the limiter is 51 ! only applied for passive scalars. 52 ! 50 53 ! 4048 2019-06-21 21:00:21Z knoop 51 54 ! Moved tcm_prognostic_equations to module_interface 52 ! 55 ! 53 56 ! 3987 2019-05-22 09:52:13Z kanani 54 57 ! Introduce alternative switch for debug output during timestepping 55 ! 58 ! 56 59 ! 3956 2019-05-07 12:32:52Z monakurppa 57 60 ! Removed salsa calls. 58 ! 61 ! 59 62 ! 3931 2019-04-24 16:34:28Z schwenkel 60 63 ! Correct/complete module_interface introduction for chemistry model 61 ! 64 ! 62 65 ! 3899 2019-04-16 14:05:27Z monakurppa 63 66 ! Corrections in the OpenMP version of salsa 64 67 ! 65 68 ! 3887 2019 -04-12 08:47:41Z schwenkel 66 ! Implicit Bugfix for chemistry model, loop for non_transport_physics over 67 ! ghost points is avoided.Instead introducing module_interface_exchange_horiz.68 ! 69 ! Implicit Bugfix for chemistry model, loop for non_transport_physics over ghost points is avoided. 70 ! Instead introducing module_interface_exchange_horiz. 71 ! 69 72 ! 3885 2019-04-11 11:29:34Z kanani 70 ! Changes related to global restructuring of location messages and introduction 71 ! of additional debugmessages72 ! 73 ! Changes related to global restructuring of location messages and introduction of additional debug 74 ! messages 75 ! 73 76 ! 3881 2019-04-10 09:31:22Z suehring 74 77 ! Bugfix in OpenMP directive 75 ! 78 ! 76 79 ! 3880 2019-04-08 21:43:02Z knoop 77 80 ! Moved wtm_tendencies to module_interface_actions 78 ! 81 ! 79 82 ! 3874 2019-04-08 16:53:48Z knoop 80 83 ! Added non_transport_physics module interfaces and moved bcm code into it 81 ! 84 ! 82 85 ! 3872 2019-04-08 15:03:06Z knoop 83 86 ! Moving prognostic equations of bcm into bulk_cloud_model_mod 84 ! 87 ! 85 88 ! 3864 2019-04-05 09:01:56Z monakurppa 86 89 ! Modifications made for salsa: 87 ! - salsa_prognostic_equations moved to salsa_mod (and the call to 88 ! module_interface_mod) 89 ! - Renamed nbins --> nbins_aerosol, ncc_tot --> ncomponents_mass and 90 ! ngast --> ngases_salsa and loop indices b, c and sg to ib, ic and ig 91 ! 90 ! - salsa_prognostic_equations moved to salsa_mod (and the call to module_interface_mod) 91 ! - Renamed nbins --> nbins_aerosol, ncc_tot --> ncomponents_mass and ngast --> ngases_salsa and 92 ! loop indices b, c and sg to ib, ic and ig 93 ! 92 94 ! 3840 2019-03-29 10:35:52Z knoop 93 ! added USE chem_gasphase_mod for nspec, nspec and spc_names94 ! 95 ! Added USE chem_gasphase_mod for nspec, nspec and spc_names 96 ! 95 97 ! 3820 2019-03-27 11:53:41Z forkel 96 ! renamed do_depo to deposition_dry (ecc)97 ! 98 ! Renamed do_depo to deposition_dry (ecc) 99 ! 98 100 ! 3797 2019-03-15 11:15:38Z forkel 99 ! Call chem_integegrate in OpenMP loop 100 ! 101 ! 101 ! Call chem_integegrate in OpenMP loop (ketelsen) 102 ! 103 ! 102 104 ! 3771 2019-02-28 12:19:33Z raasch 103 ! preprocessor directivs fro rrtmg added104 ! 105 ! Preprocessor directivs fro rrtmg added 106 ! 105 107 ! 3761 2019-02-25 15:31:42Z raasch 106 ! unused variable removed107 ! 108 ! Unused variable removed 109 ! 108 110 ! 3719 2019-02-06 13:10:18Z kanani 109 111 ! Cleaned up chemistry cpu measurements 110 ! 112 ! 111 113 ! 3684 2019-01-20 20:20:58Z knoop 112 114 ! OpenACC port for SPEC … … 115 117 ! Initial revision 116 118 ! 117 ! 119 !--------------------------------------------------------------------------------------------------! 118 120 ! Description: 119 121 ! ------------ 120 122 !> Solving the prognostic equations. 121 !------------------------------------------------------------------------------ !123 !--------------------------------------------------------------------------------------------------! 122 124 MODULE prognostic_equations_mod 123 125 124 USE advec_s_bc_mod, &126 USE advec_s_bc_mod, & 125 127 ONLY: advec_s_bc 126 128 127 USE advec_s_pw_mod, &129 USE advec_s_pw_mod, & 128 130 ONLY: advec_s_pw 129 131 130 USE advec_s_up_mod, &132 USE advec_s_up_mod, & 131 133 ONLY: advec_s_up 132 134 133 USE advec_u_pw_mod, &135 USE advec_u_pw_mod, & 134 136 ONLY: advec_u_pw 135 137 136 USE advec_u_up_mod, &138 USE advec_u_up_mod, & 137 139 ONLY: advec_u_up 138 140 139 USE advec_v_pw_mod, &141 USE advec_v_pw_mod, & 140 142 ONLY: advec_v_pw 141 143 142 USE advec_v_up_mod, &144 USE advec_v_up_mod, & 143 145 ONLY: advec_v_up 144 146 145 USE advec_w_pw_mod, &147 USE advec_w_pw_mod, & 146 148 ONLY: advec_w_pw 147 149 148 USE advec_w_up_mod, &150 USE advec_w_up_mod, & 149 151 ONLY: advec_w_up 150 152 151 USE advec_ws, & 152 ONLY: advec_s_ws, advec_u_ws, advec_v_ws, advec_w_ws 153 154 USE arrays_3d, & 155 ONLY: diss_l_e, diss_l_pt, diss_l_q, & 156 diss_l_s, diss_l_sa, diss_s_e, & 157 diss_s_pt, diss_s_q, diss_s_s, diss_s_sa, & 158 e, e_p, flux_s_e, flux_s_pt, flux_s_q, & 159 flux_s_s, flux_s_sa, flux_l_e, & 160 flux_l_pt, flux_l_q, flux_l_s, & 161 flux_l_sa, pt, ptdf_x, ptdf_y, pt_init, & 162 pt_p, prho, q, q_init, q_p, rdf, rdf_sc, & 163 ref_state, rho_ocean, s, s_init, s_p, tend, te_m, & 164 tpt_m, tq_m, ts_m, tu_m, tv_m, tw_m, u, & 165 ug, u_init, u_p, v, vg, vpt, v_init, v_p, w, w_p 166 167 USE buoyancy_mod, & 153 USE advec_ws, & 154 ONLY: advec_s_ws, & 155 advec_u_ws, & 156 advec_v_ws, & 157 advec_w_ws 158 159 USE arrays_3d, & 160 ONLY: diss_l_e, & 161 diss_l_pt, & 162 diss_l_q, & 163 diss_l_s, & 164 diss_l_sa, & 165 diss_s_e, & 166 diss_s_pt, & 167 diss_s_q, & 168 diss_s_s, & 169 diss_s_sa, & 170 e, & 171 e_p, & 172 flux_s_e, & 173 flux_s_pt, & 174 flux_s_q, & 175 flux_s_s, & 176 flux_s_sa, & 177 flux_l_e, & 178 flux_l_pt, & 179 flux_l_q, & 180 flux_l_s, & 181 flux_l_sa, & 182 pt, & 183 ptdf_x, & 184 ptdf_y, & 185 pt_init, & 186 pt_p, & 187 prho, & 188 q, & 189 q_init, & 190 q_p, & 191 rdf, & 192 rdf_sc, & 193 ref_state, & 194 rho_ocean, & 195 s, & 196 s_init, & 197 s_p, & 198 tend, & 199 te_m, & 200 tpt_m, & 201 tq_m, & 202 ts_m, & 203 tu_m, & 204 tv_m, & 205 tw_m, & 206 u, & 207 ug, & 208 u_init, & 209 u_p, & 210 v, & 211 vg, & 212 vpt, & 213 v_init, & 214 v_p, & 215 w, & 216 w_p 217 218 USE buoyancy_mod, & 168 219 ONLY: buoyancy 169 220 170 USE control_parameters, & 171 ONLY: bc_dirichlet_l, & 172 bc_dirichlet_n, & 173 bc_dirichlet_r, & 174 bc_dirichlet_s, & 175 bc_radiation_l, & 176 bc_radiation_n, & 177 bc_radiation_r, & 178 bc_radiation_s, & 179 constant_diffusion, & 180 debug_output_timestep, & 181 dp_external, dp_level_ind_b, dp_smooth_factor, dpdxy, dt_3d, & 182 humidity, intermediate_timestep_count, & 183 intermediate_timestep_count_max, large_scale_forcing, & 184 large_scale_subsidence, & 185 monotonic_limiter_z, & 186 neutral, nudging, & 187 ocean_mode, passive_scalar, plant_canopy, pt_reference, & 188 scalar_advec, scalar_advec, simulated_time, sloping_surface, & 189 timestep_scheme, tsc, use_subsidence_tendencies, & 190 use_upstream_for_tke, wind_turbine, ws_scheme_mom, & 191 ws_scheme_sca, urban_surface, land_surface, & 192 time_since_reference_point, salsa 193 194 USE coriolis_mod, & 221 USE control_parameters, & 222 ONLY: bc_dirichlet_l, & 223 bc_dirichlet_n, & 224 bc_dirichlet_r, & 225 bc_dirichlet_s, & 226 bc_radiation_l, & 227 bc_radiation_n, & 228 bc_radiation_r, & 229 bc_radiation_s, & 230 constant_diffusion, & 231 debug_output_timestep, & 232 dp_external, & 233 dp_level_ind_b, & 234 dp_smooth_factor, & 235 dpdxy, & 236 dt_3d, & 237 humidity, & 238 intermediate_timestep_count, & 239 intermediate_timestep_count_max, & 240 land_surface, & 241 large_scale_forcing, & 242 large_scale_subsidence, & 243 monotonic_limiter_z, & 244 neutral, & 245 nudging, & 246 ocean_mode, & 247 passive_scalar, & 248 plant_canopy, & 249 pt_reference, & 250 salsa, & 251 scalar_advec, & 252 scalar_advec, & 253 simulated_time, & 254 sloping_surface, & 255 time_since_reference_point, & 256 timestep_scheme, & 257 tsc, & 258 urban_surface, & 259 use_subsidence_tendencies, & 260 use_upstream_for_tke, & 261 wind_turbine, & 262 ws_scheme_mom, & 263 ws_scheme_sca 264 265 266 267 268 269 USE coriolis_mod, & 195 270 ONLY: coriolis 196 271 197 USE cpulog, & 198 ONLY: cpu_log, log_point, log_point_s 199 200 USE diffusion_s_mod, & 272 USE cpulog, & 273 ONLY: cpu_log, & 274 log_point, & 275 log_point_s 276 277 USE diffusion_s_mod, & 201 278 ONLY: diffusion_s 202 279 203 USE diffusion_u_mod, &280 USE diffusion_u_mod, & 204 281 ONLY: diffusion_u 205 282 206 USE diffusion_v_mod, &283 USE diffusion_v_mod, & 207 284 ONLY: diffusion_v 208 285 209 USE diffusion_w_mod, &286 USE diffusion_w_mod, & 210 287 ONLY: diffusion_w 211 288 212 USE indices, & 213 ONLY: advc_flags_s, & 214 nbgp, nxl, nxlg, nxlu, nxr, nxrg, nyn, nyng, nys, nysg, nysv, & 215 nzb, nzt, wall_flags_total_0 289 USE indices, & 290 ONLY: advc_flags_s, & 291 nbgp, & 292 nxl, & 293 nxlg, & 294 nxlu, & 295 nxr, & 296 nxrg, & 297 nyn, & 298 nyng, & 299 nys, & 300 nysg, & 301 nysv, & 302 nzb, & 303 nzt, & 304 wall_flags_total_0 216 305 217 306 USE kinds 218 307 219 USE lsf_nudging_mod, & 220 ONLY: ls_advec, nudge 221 222 USE module_interface, & 223 ONLY: module_interface_actions, & 224 module_interface_non_advective_processes, & 225 module_interface_exchange_horiz, & 308 USE lsf_nudging_mod, & 309 ONLY: ls_advec, & 310 nudge 311 312 USE module_interface, & 313 ONLY: module_interface_actions, & 314 module_interface_exchange_horiz, & 315 module_interface_non_advective_processes, & 226 316 module_interface_prognostic_equations 227 317 228 USE ocean_mod, & 229 ONLY: stokes_drift_terms, stokes_force, & 230 wave_breaking, wave_breaking_term 231 232 USE plant_canopy_model_mod, & 233 ONLY: cthf, pcm_tendency 318 USE ocean_mod, & 319 ONLY: stokes_drift_terms, & 320 stokes_force, & 321 wave_breaking, & 322 wave_breaking_term 323 324 USE plant_canopy_model_mod, & 325 ONLY: cthf, & 326 pcm_tendency 234 327 235 328 #if defined( __rrtmg ) 236 USE radiation_model_mod, & 237 ONLY: radiation, radiation_tendency, & 329 USE radiation_model_mod, & 330 ONLY: radiation, & 331 radiation_tendency, & 238 332 skip_time_do_radiation 239 333 #endif 240 334 241 USE statistics, &335 USE statistics, & 242 336 ONLY: hom 243 337 244 USE subsidence_mod, &338 USE subsidence_mod, & 245 339 ONLY: subsidence 246 340 247 USE surface_mod, & 248 ONLY : surf_def_h, surf_def_v, surf_lsm_h, surf_lsm_v, surf_usm_h, & 341 USE surface_mod, & 342 ONLY : surf_def_h, & 343 surf_def_v, & 344 surf_lsm_h, & 345 surf_lsm_v, & 346 surf_usm_h, & 249 347 surf_usm_v 250 348 … … 266 364 267 365 268 !------------------------------------------------------------------------------ !366 !--------------------------------------------------------------------------------------------------! 269 367 ! Description: 270 368 ! ------------ 271 !> Version with one optimized loop over all equations. It is only allowed to 272 !> be called for theWicker and Skamarock or Piascek-Williams advection scheme.369 !> Version with one optimized loop over all equations. It is only allowed to be called for the 370 !> Wicker and Skamarock or Piascek-Williams advection scheme. 273 371 !> 274 !> Here the calls of most subroutines are embedded in two DO loops over i and j, 275 !> so communication between CPUs is not allowed (does not make sense) within 276 !> these loops. 372 !> Here the calls of most subroutines are embedded in two DO loops over i and j, so communication 373 !> between CPUs is not allowed (does not make sense) within these loops. 277 374 !> 278 375 !> (Optimized to avoid cache missings, i.e. for Power4/5-architectures.) 279 !------------------------------------------------------------------------------ !376 !--------------------------------------------------------------------------------------------------! 280 377 281 378 SUBROUTINE prognostic_equations_cache … … 307 404 ENDDO 308 405 ! 309 !-- Module Inferface for exchange horiz after non_advective_processes but before 310 !-- advection.Therefore, non_advective_processes must not run for ghost points.406 !-- Module Inferface for exchange horiz after non_advective_processes but before advection. 407 !-- Therefore, non_advective_processes must not run for ghost points. 311 408 !$OMP END PARALLEL 312 409 CALL module_interface_exchange_horiz() … … 322 419 323 420 ! 324 !-- Store the first loop index. It differs for each thread and is required 325 !-- later in advec_ws 421 !-- Store the first loop index. It differs for each thread and is required later in advec_ws 326 422 IF ( loop_start ) THEN 327 423 loop_start = .FALSE. … … 331 427 DO j = nys, nyn 332 428 ! 333 !-- Tendency terms for u-velocity component. Please note, in case of 334 !-- non-cyclic boundary conditions the grid point i=0 is excluded from335 !-- the prognostic equations for theu-component.429 !-- Tendency terms for u-velocity component. Please note, in case of non-cyclic boundary 430 !-- conditions the grid point i=0 is excluded from the prognostic equations for the 431 !-- u-component. 336 432 IF ( i >= nxlu ) THEN 337 433 … … 377 473 DO k = nzb+1, nzt 378 474 379 u_p(k,j,i) = u(k,j,i) + ( dt_3d * & 380 ( tsc(2) * tend(k,j,i) + & 381 tsc(3) * tu_m(k,j,i) ) & 382 - tsc(5) * rdf(k) & 383 * ( u(k,j,i) - u_init(k) ) & 384 ) * MERGE( 1.0_wp, 0.0_wp, & 385 BTEST( wall_flags_total_0(k,j,i), 1 )& 386 ) 475 u_p(k,j,i) = u(k,j,i) + ( dt_3d * ( tsc(2) * tend(k,j,i) + tsc(3) * tu_m(k,j,i) ) & 476 - tsc(5) * rdf(k) * ( u(k,j,i) - u_init(k) ) ) & 477 * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 1 ) ) 387 478 ENDDO 388 479 389 480 ! 390 481 !-- Add turbulence generated by wave breaking (in ocean mode only) 391 IF ( wave_breaking .AND. & 392 intermediate_timestep_count == intermediate_timestep_count_max )& 393 THEN 482 IF ( wave_breaking .AND. & 483 intermediate_timestep_count == intermediate_timestep_count_max ) THEN 394 484 CALL wave_breaking_term( i, j, 1 ) 395 485 ENDIF … … 402 492 tu_m(k,j,i) = tend(k,j,i) 403 493 ENDDO 404 ELSEIF ( intermediate_timestep_count < & 405 intermediate_timestep_count_max ) THEN 494 ELSEIF ( intermediate_timestep_count < intermediate_timestep_count_max ) THEN 406 495 DO k = nzb+1, nzt 407 tu_m(k,j,i) = -9.5625_wp * tend(k,j,i) & 408 + 5.3125_wp * tu_m(k,j,i) 496 tu_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * tu_m(k,j,i) 409 497 ENDDO 410 498 ENDIF … … 413 501 ENDIF 414 502 ! 415 !-- Tendency terms for v-velocity component. Please note, in case of 416 !-- non-cyclic boundary conditions the grid point j=0 is excluded from417 !-- the prognostic equations for the v-component. !--503 !-- Tendency terms for v-velocity component. Please note, in case of non-cyclic boundary 504 !-- conditions the grid point j=0 is excluded from the prognostic equations for the 505 !-- v-component. 418 506 IF ( j >= nysv ) THEN 419 507 … … 455 543 !-- Prognostic equation for v-velocity component 456 544 DO k = nzb+1, nzt 457 v_p(k,j,i) = v(k,j,i) + ( dt_3d * & 458 ( tsc(2) * tend(k,j,i) + & 459 tsc(3) * tv_m(k,j,i) ) & 460 - tsc(5) * rdf(k) & 461 * ( v(k,j,i) - v_init(k) )& 462 ) * MERGE( 1.0_wp, 0.0_wp, & 463 BTEST( wall_flags_total_0(k,j,i), 2 )& 464 ) 545 v_p(k,j,i) = v(k,j,i) + ( dt_3d * ( tsc(2) * tend(k,j,i) + tsc(3) * tv_m(k,j,i) ) & 546 - tsc(5) * rdf(k) * ( v(k,j,i) - v_init(k) ) ) & 547 * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 2 ) ) 465 548 ENDDO 466 549 467 550 ! 468 551 !-- Add turbulence generated by wave breaking (in ocean mode only) 469 IF ( wave_breaking .AND. & 470 intermediate_timestep_count == intermediate_timestep_count_max )& 471 THEN 552 IF ( wave_breaking .AND. & 553 intermediate_timestep_count == intermediate_timestep_count_max ) THEN 472 554 CALL wave_breaking_term( i, j, 2 ) 473 555 ENDIF … … 480 562 tv_m(k,j,i) = tend(k,j,i) 481 563 ENDDO 482 ELSEIF ( intermediate_timestep_count < & 483 intermediate_timestep_count_max ) THEN 564 ELSEIF ( intermediate_timestep_count < intermediate_timestep_count_max ) THEN 484 565 DO k = nzb+1, nzt 485 tv_m(k,j,i) = -9.5625_wp * tend(k,j,i) & 486 + 5.3125_wp * tv_m(k,j,i) 566 tv_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * tv_m(k,j,i) 487 567 ENDDO 488 568 ENDIF … … 530 610 !-- Prognostic equation for w-velocity component 531 611 DO k = nzb+1, nzt-1 532 w_p(k,j,i) = w(k,j,i) + ( dt_3d * & 533 ( tsc(2) * tend(k,j,i) + & 534 tsc(3) * tw_m(k,j,i) ) & 535 - tsc(5) * rdf(k) * w(k,j,i) & 536 ) * MERGE( 1.0_wp, 0.0_wp, & 537 BTEST( wall_flags_total_0(k,j,i), 3 )& 538 ) 612 w_p(k,j,i) = w(k,j,i) + ( dt_3d * ( tsc(2) * tend(k,j,i) + tsc(3) * tw_m(k,j,i) ) & 613 - tsc(5) * rdf(k) * w(k,j,i) ) & 614 * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 3 ) ) 539 615 ENDDO 540 616 … … 546 622 tw_m(k,j,i) = tend(k,j,i) 547 623 ENDDO 548 ELSEIF ( intermediate_timestep_count < & 549 intermediate_timestep_count_max ) THEN 624 ELSEIF ( intermediate_timestep_count < intermediate_timestep_count_max ) THEN 550 625 DO k = nzb+1, nzt-1 551 tw_m(k,j,i) = -9.5625_wp * tend(k,j,i) & 552 + 5.3125_wp * tw_m(k,j,i) 626 tw_m(k,j,i) = -9.5625_wp * tend(k,j,i) + 5.3125_wp * tw_m(k,j,i) 553 627 ENDDO 554 628 ENDIF … … 563 637 IF ( timestep_scheme(1:5) == 'runge' ) THEN 564 638 IF ( ws_scheme_sca ) THEN 565 CALL advec_s_ws( advc_flags_s, & 566 i, j, pt, 'pt', flux_s_pt, diss_s_pt, & 567 flux_l_pt, diss_l_pt, i_omp_start, tn, & 568 bc_dirichlet_l .OR. bc_radiation_l, & 569 bc_dirichlet_n .OR. bc_radiation_n, & 570 bc_dirichlet_r .OR. bc_radiation_r, & 639 CALL advec_s_ws( advc_flags_s, i, j, pt, 'pt', flux_s_pt, diss_s_pt, & 640 flux_l_pt, diss_l_pt, i_omp_start, tn, & 641 bc_dirichlet_l .OR. bc_radiation_l, & 642 bc_dirichlet_n .OR. bc_radiation_n, & 643 bc_dirichlet_r .OR. bc_radiation_r, & 571 644 bc_dirichlet_s .OR. bc_radiation_s ) 572 645 ELSE … … 576 649 CALL advec_s_up( i, j, pt ) 577 650 ENDIF 578 CALL diffusion_s( i, j, pt, & 579 surf_def_h(0)%shf, surf_def_h(1)%shf, & 580 surf_def_h(2)%shf, & 581 surf_lsm_h%shf, surf_usm_h%shf, & 582 surf_def_v(0)%shf, surf_def_v(1)%shf, & 583 surf_def_v(2)%shf, surf_def_v(3)%shf, & 584 surf_lsm_v(0)%shf, surf_lsm_v(1)%shf, & 585 surf_lsm_v(2)%shf, surf_lsm_v(3)%shf, & 586 surf_usm_v(0)%shf, surf_usm_v(1)%shf, & 651 CALL diffusion_s( i, j, pt, surf_def_h(0)%shf, surf_def_h(1)%shf, surf_def_h(2)%shf, & 652 surf_lsm_h%shf, surf_usm_h%shf, surf_def_v(0)%shf, & 653 surf_def_v(1)%shf, surf_def_v(2)%shf, surf_def_v(3)%shf, & 654 surf_lsm_v(0)%shf, surf_lsm_v(1)%shf, surf_lsm_v(2)%shf, & 655 surf_lsm_v(3)%shf, surf_usm_v(0)%shf, surf_usm_v(1)%shf, & 587 656 surf_usm_v(2)%shf, surf_usm_v(3)%shf ) 588 657 589 658 ! 590 659 !-- Consideration of heat sources within the plant canopy 591 IF ( plant_canopy .AND. 592 (cthf /= 0.0_wp .OR. urban_surface .OR. land_surface) )THEN660 IF ( plant_canopy .AND. (cthf /= 0.0_wp .OR. urban_surface .OR. land_surface) ) & 661 THEN 593 662 CALL pcm_tendency( i, j, 4 ) 594 663 ENDIF … … 606 675 ! 607 676 !-- If required, compute effect of large-scale subsidence/ascent 608 IF ( large_scale_subsidence .AND. & 609 .NOT. use_subsidence_tendencies ) THEN 677 IF ( large_scale_subsidence .AND. .NOT. use_subsidence_tendencies ) THEN 610 678 CALL subsidence( i, j, tend, pt, pt_init, 2 ) 611 679 ENDIF … … 614 682 ! 615 683 !-- If required, add tendency due to radiative heating/cooling 616 IF ( radiation .AND. & 617 simulated_time > skip_time_do_radiation ) THEN 684 IF ( radiation .AND. simulated_time > skip_time_do_radiation ) THEN 618 685 CALL radiation_tendency ( i, j, tend ) 619 686 ENDIF … … 624 691 !-- Prognostic equation for potential temperature 625 692 DO k = nzb+1, nzt 626 pt_p(k,j,i) = pt(k,j,i) + ( dt_3d * & 627 ( tsc(2) * tend(k,j,i) + & 628 tsc(3) * tpt_m(k,j,i) ) & 629 - tsc(5) & 630 * ( pt(k,j,i) - pt_init(k) ) & 631 * ( rdf_sc(k) + ptdf_x(i) & 632 + ptdf_y(j) ) & 633 ) & 634 * MERGE( 1.0_wp, 0.0_wp, & 635 BTEST( wall_flags_total_0(k,j,i), 0 )& 636 ) 693 pt_p(k,j,i) = pt(k,j,i) + & 694 ( dt_3d * ( tsc(2) * tend(k,j,i) + tsc(3) * tpt_m(k,j,i) ) - tsc(5) & 695 * ( pt(k,j,i) - pt_init(k) ) * ( rdf_sc(k) + ptdf_x(i) & 696 + ptdf_y(j) ) ) & 697 * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) ) 637 698 ENDDO 638 699 … … 644 705 tpt_m(k,j,i) = tend(k,j,i) 645 706 ENDDO 646 ELSEIF ( intermediate_timestep_count < & 647 intermediate_timestep_count_max ) THEN 707 ELSEIF ( intermediate_timestep_count < intermediate_timestep_count_max ) THEN 648 708 DO k = nzb+1, nzt 649 tpt_m(k,j,i) = -9.5625_wp * tend(k,j,i) + & 650 5.3125_wp * tpt_m(k,j,i) 709 tpt_m(k,j,i) = - 9.5625_wp * tend(k,j,i) + 5.3125_wp * tpt_m(k,j,i) 651 710 ENDDO 652 711 ENDIF … … 662 721 !-- Tendency-terms for total water content / scalar 663 722 tend(:,j,i) = 0.0_wp 664 IF ( timestep_scheme(1:5) == 'runge' ) & 665 THEN 723 IF ( timestep_scheme(1:5) == 'runge' ) THEN 666 724 IF ( ws_scheme_sca ) THEN 667 CALL advec_s_ws( advc_flags_s, & 668 i, j, q, 'q', flux_s_q, & 669 diss_s_q, flux_l_q, diss_l_q, & 670 i_omp_start, tn, & 671 bc_dirichlet_l .OR. bc_radiation_l, & 672 bc_dirichlet_n .OR. bc_radiation_n, & 673 bc_dirichlet_r .OR. bc_radiation_r, & 725 CALL advec_s_ws( advc_flags_s, i, j, q, 'q', flux_s_q, diss_s_q, flux_l_q, & 726 diss_l_q, i_omp_start, tn, & 727 bc_dirichlet_l .OR. bc_radiation_l, & 728 bc_dirichlet_n .OR. bc_radiation_n, & 729 bc_dirichlet_r .OR. bc_radiation_r, & 674 730 bc_dirichlet_s .OR. bc_radiation_s ) 675 731 ELSE … … 679 735 CALL advec_s_up( i, j, q ) 680 736 ENDIF 681 CALL diffusion_s( i, j, q, & 682 surf_def_h(0)%qsws, surf_def_h(1)%qsws, & 683 surf_def_h(2)%qsws, & 684 surf_lsm_h%qsws, surf_usm_h%qsws, & 685 surf_def_v(0)%qsws, surf_def_v(1)%qsws, & 686 surf_def_v(2)%qsws, surf_def_v(3)%qsws, & 687 surf_lsm_v(0)%qsws, surf_lsm_v(1)%qsws, & 688 surf_lsm_v(2)%qsws, surf_lsm_v(3)%qsws, & 689 surf_usm_v(0)%qsws, surf_usm_v(1)%qsws, & 690 surf_usm_v(2)%qsws, surf_usm_v(3)%qsws ) 737 CALL diffusion_s( i, j, q, surf_def_h(0)%qsws, surf_def_h(1)%qsws, & 738 surf_def_h(2)%qsws, surf_lsm_h%qsws, surf_usm_h%qsws, & 739 surf_def_v(0)%qsws, surf_def_v(1)%qsws, surf_def_v(2)%qsws, & 740 surf_def_v(3)%qsws, surf_lsm_v(0)%qsws, surf_lsm_v(1)%qsws, & 741 surf_lsm_v(2)%qsws, surf_lsm_v(3)%qsws, surf_usm_v(0)%qsws, & 742 surf_usm_v(1)%qsws, surf_usm_v(2)%qsws, surf_usm_v(3)%qsws ) 691 743 692 744 ! … … 706 758 ! 707 759 !-- If required compute influence of large-scale subsidence/ascent 708 IF ( large_scale_subsidence .AND. & 709 .NOT. use_subsidence_tendencies ) THEN 760 IF ( large_scale_subsidence .AND. .NOT. use_subsidence_tendencies ) THEN 710 761 CALL subsidence( i, j, tend, q, q_init, 3 ) 711 762 ENDIF … … 716 767 !-- Prognostic equation for total water content / scalar 717 768 DO k = nzb+1, nzt 718 q_p(k,j,i) = q(k,j,i) + ( dt_3d * & 719 ( tsc(2) * tend(k,j,i) + & 720 tsc(3) * tq_m(k,j,i) ) & 721 - tsc(5) * rdf_sc(k) * & 722 ( q(k,j,i) - q_init(k) ) & 723 ) & 724 * MERGE( 1.0_wp, 0.0_wp, & 725 BTEST( wall_flags_total_0(k,j,i), 0 )& 726 ) 769 q_p(k,j,i) = q(k,j,i) & 770 + ( dt_3d * ( tsc(2) * tend(k,j,i) + tsc(3) * tq_m(k,j,i) ) - tsc(5) & 771 * rdf_sc(k) * ( q(k,j,i) - q_init(k) ) ) & 772 * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) ) 727 773 IF ( q_p(k,j,i) < 0.0_wp ) q_p(k,j,i) = 0.1_wp * q(k,j,i) 728 774 ENDDO … … 738 784 intermediate_timestep_count_max ) THEN 739 785 DO k = nzb+1, nzt 740 tq_m(k,j,i) = -9.5625_wp * tend(k,j,i) + & 741 5.3125_wp * tq_m(k,j,i) 786 tq_m(k,j,i) = - 9.5625_wp * tend(k,j,i) + 5.3125_wp * tq_m(k,j,i) 742 787 ENDDO 743 788 ENDIF … … 752 797 !-- Tendency-terms for total water content / scalar 753 798 tend(:,j,i) = 0.0_wp 754 IF ( timestep_scheme(1:5) == 'runge' ) & 755 THEN 799 IF ( timestep_scheme(1:5) == 'runge' ) THEN 756 800 IF ( ws_scheme_sca ) THEN 757 801 ! 758 !-- For scalar advection apply monotonic flux limiter near 759 !-- topography. 760 CALL advec_s_ws( advc_flags_s, & 761 i, j, s, 's', flux_s_s, & 762 diss_s_s, flux_l_s, diss_l_s, i_omp_start, & 763 tn, & 764 bc_dirichlet_l .OR. bc_radiation_l, & 765 bc_dirichlet_n .OR. bc_radiation_n, & 766 bc_dirichlet_r .OR. bc_radiation_r, & 767 bc_dirichlet_s .OR. bc_radiation_s, & 802 !-- For scalar advection apply monotonic flux limiter near topography. 803 CALL advec_s_ws( advc_flags_s, i, j, s, 's', flux_s_s, diss_s_s, flux_l_s, & 804 diss_l_s, i_omp_start, tn, & 805 bc_dirichlet_l .OR. bc_radiation_l, & 806 bc_dirichlet_n .OR. bc_radiation_n, & 807 bc_dirichlet_r .OR. bc_radiation_r, & 808 bc_dirichlet_s .OR. bc_radiation_s, & 768 809 monotonic_limiter_z ) 769 810 ELSE … … 773 814 CALL advec_s_up( i, j, s ) 774 815 ENDIF 775 CALL diffusion_s( i, j, s, & 776 surf_def_h(0)%ssws, surf_def_h(1)%ssws, & 777 surf_def_h(2)%ssws, & 778 surf_lsm_h%ssws, surf_usm_h%ssws, & 779 surf_def_v(0)%ssws, surf_def_v(1)%ssws, & 780 surf_def_v(2)%ssws, surf_def_v(3)%ssws, & 781 surf_lsm_v(0)%ssws, surf_lsm_v(1)%ssws, & 782 surf_lsm_v(2)%ssws, surf_lsm_v(3)%ssws, & 783 surf_usm_v(0)%ssws, surf_usm_v(1)%ssws, & 784 surf_usm_v(2)%ssws, surf_usm_v(3)%ssws ) 816 CALL diffusion_s( i, j, s, surf_def_h(0)%ssws, surf_def_h(1)%ssws, & 817 surf_def_h(2)%ssws, surf_lsm_h%ssws, surf_usm_h%ssws, & 818 surf_def_v(0)%ssws, surf_def_v(1)%ssws, surf_def_v(2)%ssws, & 819 surf_def_v(3)%ssws, surf_lsm_v(0)%ssws, surf_lsm_v(1)%ssws, & 820 surf_lsm_v(2)%ssws, surf_lsm_v(3)%ssws, surf_usm_v(0)%ssws, & 821 surf_usm_v(1)%ssws, surf_usm_v(2)%ssws, surf_usm_v(3)%ssws ) 785 822 786 823 ! … … 799 836 800 837 ! 801 !-- If required compute influence of large-scale subsidence/ascent. 802 !-- Note, the last argument is of no meaning in this case, as it is 803 !-- only used in conjunction with large_scale_forcing, which is to 804 !-- date not implemented for scalars. 805 IF ( large_scale_subsidence .AND. & 806 .NOT. use_subsidence_tendencies .AND. & 838 !-- If required compute influence of large-scale subsidence/ascent. Note, the last argument 839 !-- is of no meaning in this case, as it is only used in conjunction with 840 !-- large_scale_forcing, which is to date not implemented for scalars. 841 IF ( large_scale_subsidence .AND. .NOT. use_subsidence_tendencies .AND. & 807 842 .NOT. large_scale_forcing ) THEN 808 843 CALL subsidence( i, j, tend, s, s_init, 3 ) … … 814 849 !-- Prognostic equation for scalar 815 850 DO k = nzb+1, nzt 816 s_p(k,j,i) = s(k,j,i) + ( dt_3d * & 817 ( tsc(2) * tend(k,j,i) + & 818 tsc(3) * ts_m(k,j,i) ) & 819 - tsc(5) * rdf_sc(k) & 820 * ( s(k,j,i) - s_init(k) )& 821 ) & 822 * MERGE( 1.0_wp, 0.0_wp, & 823 BTEST( wall_flags_total_0(k,j,i), 0 )& 824 ) 851 s_p(k,j,i) = s(k,j,i) & 852 + ( dt_3d * ( tsc(2) * tend(k,j,i) + tsc(3) * ts_m(k,j,i) ) & 853 - tsc(5) * rdf_sc(k) * ( s(k,j,i) - s_init(k) ) ) & 854 * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) ) 825 855 IF ( s_p(k,j,i) < 0.0_wp ) s_p(k,j,i) = 0.1_wp * s(k,j,i) 826 856 ENDDO … … 833 863 ts_m(k,j,i) = tend(k,j,i) 834 864 ENDDO 835 ELSEIF ( intermediate_timestep_count < & 836 intermediate_timestep_count_max ) THEN 865 ELSEIF ( intermediate_timestep_count < intermediate_timestep_count_max ) THEN 837 866 DO k = nzb+1, nzt 838 ts_m(k,j,i) = -9.5625_wp * tend(k,j,i) + & 839 5.3125_wp * ts_m(k,j,i) 867 ts_m(k,j,i) = - 9.5625_wp * tend(k,j,i) + 5.3125_wp * ts_m(k,j,i) 840 868 ENDDO 841 869 ENDIF … … 859 887 860 888 861 !------------------------------------------------------------------------------ !889 !--------------------------------------------------------------------------------------------------! 862 890 ! Description: 863 891 ! ------------ 864 892 !> Version for vector machines 865 !------------------------------------------------------------------------------ !893 !--------------------------------------------------------------------------------------------------! 866 894 867 895 SUBROUTINE prognostic_equations_vector 868 896 869 897 870 INTEGER(iwp) :: i 871 INTEGER(iwp) :: j 872 INTEGER(iwp) :: k 873 874 REAL(wp) 898 INTEGER(iwp) :: i !< 899 INTEGER(iwp) :: j !< 900 INTEGER(iwp) :: k !< 901 902 REAL(wp) :: sbt !< 875 903 876 904 … … 880 908 CALL module_interface_non_advective_processes 881 909 ! 882 !-- Module Inferface for exchange horiz after non_advective_processes but before 883 !-- advection.Therefore, non_advective_processes must not run for ghost points.910 !-- Module Inferface for exchange horiz after non_advective_processes but before advection. 911 !-- Therefore, non_advective_processes must not run for ghost points. 884 912 CALL module_interface_exchange_horiz() 885 913 ! … … 939 967 DO i = nxlu, nxr 940 968 DO j = nys, nyn 941 !following directive is required to vectorize on Intel19 969 ! 970 !-- Following directive is required to vectorize on Intel19 942 971 !DIR$ IVDEP 943 972 DO k = nzb+1, nzt 944 u_p(k,j,i) = u(k,j,i) + ( dt_3d * ( tsc(2) * tend(k,j,i) + & 945 tsc(3) * tu_m(k,j,i) ) & 946 - tsc(5) * rdf(k) * & 947 ( u(k,j,i) - u_init(k) ) & 948 ) * MERGE( 1.0_wp, 0.0_wp, & 949 BTEST( wall_flags_total_0(k,j,i), 1 ) & 950 ) 973 u_p(k,j,i) = u(k,j,i) & 974 + ( dt_3d * ( tsc(2) * tend(k,j,i) + tsc(3) * tu_m(k,j,i) ) & 975 - tsc(5) * rdf(k) * ( u(k,j,i) - u_init(k) ) ) & 976 * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 1 ) ) 951 977 ENDDO 952 978 ENDDO … … 955 981 ! 956 982 !-- Add turbulence generated by wave breaking (in ocean mode only) 957 IF ( wave_breaking .AND. & 958 intermediate_timestep_count == intermediate_timestep_count_max ) & 959 THEN 983 IF ( wave_breaking .AND. intermediate_timestep_count == intermediate_timestep_count_max ) THEN 960 984 CALL wave_breaking_term( 1 ) 961 985 ENDIF … … 974 998 ENDDO 975 999 ENDDO 976 ELSEIF ( intermediate_timestep_count < & 977 intermediate_timestep_count_max ) THEN 1000 ELSEIF ( intermediate_timestep_count < intermediate_timestep_count_max ) THEN 978 1001 !$ACC PARALLEL LOOP COLLAPSE(3) PRIVATE(i, j, k) & 979 1002 !$ACC PRESENT(tend, tu_m) … … 981 1004 DO j = nys, nyn 982 1005 DO k = nzb+1, nzt 983 tu_m(k,j,i) = -9.5625_wp * tend(k,j,i) & 984 + 5.3125_wp * tu_m(k,j,i) 1006 tu_m(k,j,i) = - 9.5625_wp * tend(k,j,i) + 5.3125_wp * tu_m(k,j,i) 985 1007 ENDDO 986 1008 ENDDO … … 1044 1066 DO i = nxl, nxr 1045 1067 DO j = nysv, nyn 1046 !following directive is required to vectorize on Intel19 1068 ! 1069 !-- Following directive is required to vectorize on Intel19 1047 1070 !DIR$ IVDEP 1048 1071 DO k = nzb+1, nzt 1049 v_p(k,j,i) = v(k,j,i) + ( dt_3d * ( tsc(2) * tend(k,j,i) + & 1050 tsc(3) * tv_m(k,j,i) ) & 1051 - tsc(5) * rdf(k) * & 1052 ( v(k,j,i) - v_init(k) ) & 1053 ) * MERGE( 1.0_wp, 0.0_wp, & 1054 BTEST( wall_flags_total_0(k,j,i), 2 )& 1055 ) 1072 v_p(k,j,i) = v(k,j,i) & 1073 + ( dt_3d * ( tsc(2) * tend(k,j,i) + tsc(3) * tv_m(k,j,i) ) - tsc(5) & 1074 * rdf(k) * ( v(k,j,i) - v_init(k) ) ) & 1075 * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 2 ) ) 1056 1076 ENDDO 1057 1077 ENDDO … … 1060 1080 ! 1061 1081 !-- Add turbulence generated by wave breaking (in ocean mode only) 1062 IF ( wave_breaking .AND. & 1063 intermediate_timestep_count == intermediate_timestep_count_max ) & 1064 THEN 1082 IF ( wave_breaking .AND. intermediate_timestep_count == intermediate_timestep_count_max ) THEN 1065 1083 CALL wave_breaking_term( 2 ) 1066 1084 ENDIF … … 1079 1097 ENDDO 1080 1098 ENDDO 1081 ELSEIF ( intermediate_timestep_count < & 1082 intermediate_timestep_count_max ) THEN 1099 ELSEIF ( intermediate_timestep_count < intermediate_timestep_count_max ) THEN 1083 1100 !$ACC PARALLEL LOOP COLLAPSE(3) PRIVATE(i, j, k) & 1084 1101 !$ACC PRESENT(tend, tv_m) … … 1086 1103 DO j = nysv, nyn 1087 1104 DO k = nzb+1, nzt 1088 tv_m(k,j,i) = -9.5625_wp * tend(k,j,i) & 1089 + 5.3125_wp * tv_m(k,j,i) 1105 tv_m(k,j,i) = - 9.5625_wp * tend(k,j,i) + 5.3125_wp * tv_m(k,j,i) 1090 1106 ENDDO 1091 1107 ENDDO … … 1145 1161 DO i = nxl, nxr 1146 1162 DO j = nys, nyn 1147 !following directive is required to vectorize on Intel19 1163 ! 1164 !-- Following directive is required to vectorize on Intel19 1148 1165 !DIR$ IVDEP 1149 1166 DO k = nzb+1, nzt-1 1150 w_p(k,j,i) = w(k,j,i) + ( dt_3d * ( tsc(2) * tend(k,j,i) + & 1151 tsc(3) * tw_m(k,j,i) ) & 1152 - tsc(5) * rdf(k) * w(k,j,i) & 1153 ) * MERGE( 1.0_wp, 0.0_wp, & 1154 BTEST( wall_flags_total_0(k,j,i), 3 )& 1155 ) 1167 w_p(k,j,i) = w(k,j,i) & 1168 + ( dt_3d * ( tsc(2) * tend(k,j,i) + tsc(3) * tw_m(k,j,i) ) & 1169 - tsc(5) * rdf(k) * w(k,j,i) ) & 1170 * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 3 ) ) 1156 1171 ENDDO 1157 1172 ENDDO … … 1171 1186 ENDDO 1172 1187 ENDDO 1173 ELSEIF ( intermediate_timestep_count < & 1174 intermediate_timestep_count_max ) THEN 1188 ELSEIF ( intermediate_timestep_count < intermediate_timestep_count_max ) THEN 1175 1189 !$ACC PARALLEL LOOP COLLAPSE(3) PRIVATE(i, j, k) & 1176 1190 !$ACC PRESENT(tend, tw_m) … … 1178 1192 DO j = nys, nyn 1179 1193 DO k = nzb+1, nzt-1 1180 tw_m(k,j,i) = -9.5625_wp * tend(k,j,i) & 1181 + 5.3125_wp * tw_m(k,j,i) 1194 tw_m(k,j,i) = - 9.5625_wp * tend(k,j,i) + 5.3125_wp * tw_m(k,j,i) 1182 1195 ENDDO 1183 1196 ENDDO … … 1218 1231 IF ( timestep_scheme(1:5) == 'runge' ) THEN 1219 1232 IF ( ws_scheme_sca ) THEN 1220 CALL advec_s_ws( advc_flags_s, pt, 'pt', &1221 bc_dirichlet_l .OR. bc_radiation_l, &1222 bc_dirichlet_n .OR. bc_radiation_n, &1223 bc_dirichlet_r .OR. bc_radiation_r, &1233 CALL advec_s_ws( advc_flags_s, pt, 'pt', & 1234 bc_dirichlet_l .OR. bc_radiation_l, & 1235 bc_dirichlet_n .OR. bc_radiation_n, & 1236 bc_dirichlet_r .OR. bc_radiation_r, & 1224 1237 bc_dirichlet_s .OR. bc_radiation_s ) 1225 1238 ELSE … … 1231 1244 ENDIF 1232 1245 1233 CALL diffusion_s( pt, & 1234 surf_def_h(0)%shf, surf_def_h(1)%shf, & 1235 surf_def_h(2)%shf, & 1236 surf_lsm_h%shf, surf_usm_h%shf, & 1237 surf_def_v(0)%shf, surf_def_v(1)%shf, & 1238 surf_def_v(2)%shf, surf_def_v(3)%shf, & 1239 surf_lsm_v(0)%shf, surf_lsm_v(1)%shf, & 1240 surf_lsm_v(2)%shf, surf_lsm_v(3)%shf, & 1241 surf_usm_v(0)%shf, surf_usm_v(1)%shf, & 1242 surf_usm_v(2)%shf, surf_usm_v(3)%shf ) 1246 CALL diffusion_s( pt, surf_def_h(0)%shf, surf_def_h(1)%shf, surf_def_h(2)%shf, & 1247 surf_lsm_h%shf, surf_usm_h%shf, surf_def_v(0)%shf, surf_def_v(1)%shf, & 1248 surf_def_v(2)%shf, surf_def_v(3)%shf, surf_lsm_v(0)%shf, & 1249 surf_lsm_v(1)%shf, surf_lsm_v(2)%shf, surf_lsm_v(3)%shf, & 1250 surf_usm_v(0)%shf, surf_usm_v(1)%shf, surf_usm_v(2)%shf, & 1251 surf_usm_v(3)%shf ) 1243 1252 1244 1253 ! 1245 1254 !-- Consideration of heat sources within the plant canopy 1246 IF ( plant_canopy .AND. & 1247 (cthf /= 0.0_wp .OR. urban_surface .OR. land_surface) ) THEN 1255 IF ( plant_canopy .AND. (cthf /= 0.0_wp .OR. urban_surface .OR. land_surface) ) THEN 1248 1256 CALL pcm_tendency( 4 ) 1249 1257 ENDIF … … 1261 1269 ! 1262 1270 !-- If required compute influence of large-scale subsidence/ascent 1263 IF ( large_scale_subsidence .AND. & 1264 .NOT. use_subsidence_tendencies ) THEN 1271 IF ( large_scale_subsidence .AND. .NOT. use_subsidence_tendencies ) THEN 1265 1272 CALL subsidence( tend, pt, pt_init, 2 ) 1266 1273 ENDIF … … 1269 1276 ! 1270 1277 !-- If required, add tendency due to radiative heating/cooling 1271 IF ( radiation .AND. & 1272 simulated_time > skip_time_do_radiation ) THEN 1278 IF ( radiation .AND. simulated_time > skip_time_do_radiation ) THEN 1273 1279 CALL radiation_tendency ( tend ) 1274 1280 ENDIF … … 1286 1292 DO i = nxl, nxr 1287 1293 DO j = nys, nyn 1288 !following directive is required to vectorize on Intel19 1294 ! 1295 !-- Following directive is required to vectorize on Intel19 1289 1296 !DIR$ IVDEP 1290 1297 DO k = nzb+1, nzt 1291 pt_p(k,j,i) = pt(k,j,i) + ( dt_3d * ( sbt * tend(k,j,i) + & 1292 tsc(3) * tpt_m(k,j,i) ) & 1293 - tsc(5) * & 1294 ( pt(k,j,i) - pt_init(k) ) *& 1295 ( rdf_sc(k) + ptdf_x(i) + ptdf_y(j) )& 1296 ) & 1297 * MERGE( 1.0_wp, 0.0_wp, & 1298 BTEST( wall_flags_total_0(k,j,i), 0 ) & 1299 ) 1298 pt_p(k,j,i) = pt(k,j,i) & 1299 + ( dt_3d * ( sbt * tend(k,j,i) + tsc(3) * tpt_m(k,j,i) ) & 1300 - tsc(5) * ( pt(k,j,i) - pt_init(k) ) & 1301 * ( rdf_sc(k) + ptdf_x(i) + ptdf_y(j) ) ) & 1302 * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) ) 1300 1303 ENDDO 1301 1304 ENDDO … … 1314 1317 ENDDO 1315 1318 ENDDO 1316 ELSEIF ( intermediate_timestep_count < & 1317 intermediate_timestep_count_max ) THEN 1319 ELSEIF ( intermediate_timestep_count < intermediate_timestep_count_max ) THEN 1318 1320 !$ACC PARALLEL LOOP COLLAPSE(3) PRIVATE(i, j, k) & 1319 1321 !$ACC PRESENT(tend, tpt_m) … … 1321 1323 DO j = nys, nyn 1322 1324 DO k = nzb+1, nzt 1323 tpt_m(k,j,i) = -9.5625_wp * tend(k,j,i) + & 1324 5.3125_wp * tpt_m(k,j,i) 1325 tpt_m(k,j,i) = - 9.5625_wp * tend(k,j,i) + 5.3125_wp * tpt_m(k,j,i) 1325 1326 ENDDO 1326 1327 ENDDO … … 1360 1361 IF ( timestep_scheme(1:5) == 'runge' ) THEN 1361 1362 IF ( ws_scheme_sca ) THEN 1362 CALL advec_s_ws( advc_flags_s, q, 'q', &1363 bc_dirichlet_l .OR. bc_radiation_l, &1364 bc_dirichlet_n .OR. bc_radiation_n, &1365 bc_dirichlet_r .OR. bc_radiation_r, &1363 CALL advec_s_ws( advc_flags_s, q, 'q', & 1364 bc_dirichlet_l .OR. bc_radiation_l, & 1365 bc_dirichlet_n .OR. bc_radiation_n, & 1366 bc_dirichlet_r .OR. bc_radiation_r, & 1366 1367 bc_dirichlet_s .OR. bc_radiation_s ) 1367 1368 ELSE … … 1373 1374 ENDIF 1374 1375 1375 CALL diffusion_s( q, & 1376 surf_def_h(0)%qsws, surf_def_h(1)%qsws, & 1377 surf_def_h(2)%qsws, & 1378 surf_lsm_h%qsws, surf_usm_h%qsws, & 1379 surf_def_v(0)%qsws, surf_def_v(1)%qsws, & 1380 surf_def_v(2)%qsws, surf_def_v(3)%qsws, & 1381 surf_lsm_v(0)%qsws, surf_lsm_v(1)%qsws, & 1382 surf_lsm_v(2)%qsws, surf_lsm_v(3)%qsws, & 1383 surf_usm_v(0)%qsws, surf_usm_v(1)%qsws, & 1384 surf_usm_v(2)%qsws, surf_usm_v(3)%qsws ) 1376 CALL diffusion_s( q, surf_def_h(0)%qsws, surf_def_h(1)%qsws, surf_def_h(2)%qsws, & 1377 surf_lsm_h%qsws, surf_usm_h%qsws, surf_def_v(0)%qsws, surf_def_v(1)%qsws, & 1378 surf_def_v(2)%qsws, surf_def_v(3)%qsws, surf_lsm_v(0)%qsws, & 1379 surf_lsm_v(1)%qsws, surf_lsm_v(2)%qsws, surf_lsm_v(3)%qsws, & 1380 surf_usm_v(0)%qsws, surf_usm_v(1)%qsws, surf_usm_v(2)%qsws, & 1381 surf_usm_v(3)%qsws ) 1385 1382 1386 1383 ! … … 1400 1397 ! 1401 1398 !-- If required compute influence of large-scale subsidence/ascent 1402 IF ( large_scale_subsidence .AND. & 1403 .NOT. use_subsidence_tendencies ) THEN 1399 IF ( large_scale_subsidence .AND. .NOT. use_subsidence_tendencies ) THEN 1404 1400 CALL subsidence( tend, q, q_init, 3 ) 1405 1401 ENDIF … … 1411 1407 DO i = nxl, nxr 1412 1408 DO j = nys, nyn 1413 !following directive is required to vectorize on Intel19 1409 ! 1410 !-- Following directive is required to vectorize on Intel19 1414 1411 !DIR$ IVDEP 1415 1412 DO k = nzb+1, nzt 1416 q_p(k,j,i) = q(k,j,i) + ( dt_3d * ( sbt * tend(k,j,i) + & 1417 tsc(3) * tq_m(k,j,i) ) & 1418 - tsc(5) * rdf_sc(k) * & 1419 ( q(k,j,i) - q_init(k) ) & 1420 ) * MERGE( 1.0_wp, 0.0_wp, & 1421 BTEST( wall_flags_total_0(k,j,i), 0 ) & 1422 ) 1413 q_p(k,j,i) = q(k,j,i) & 1414 + ( dt_3d * ( sbt * tend(k,j,i) + tsc(3) * tq_m(k,j,i) ) & 1415 - tsc(5) * rdf_sc(k) * ( q(k,j,i) - q_init(k) ) ) & 1416 * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) ) 1423 1417 IF ( q_p(k,j,i) < 0.0_wp ) q_p(k,j,i) = 0.1_wp * q(k,j,i) 1424 1418 ENDDO … … 1437 1431 ENDDO 1438 1432 ENDDO 1439 ELSEIF ( intermediate_timestep_count < & 1440 intermediate_timestep_count_max ) THEN 1433 ELSEIF ( intermediate_timestep_count < intermediate_timestep_count_max ) THEN 1441 1434 DO i = nxl, nxr 1442 1435 DO j = nys, nyn 1443 1436 DO k = nzb+1, nzt 1444 tq_m(k,j,i) = -9.5625_wp * tend(k,j,i) & 1445 + 5.3125_wp * tq_m(k,j,i) 1437 tq_m(k,j,i) = - 9.5625_wp * tend(k,j,i) + 5.3125_wp * tq_m(k,j,i) 1446 1438 ENDDO 1447 1439 ENDDO … … 1480 1472 IF ( timestep_scheme(1:5) == 'runge' ) THEN 1481 1473 IF ( ws_scheme_sca ) THEN 1482 CALL advec_s_ws( advc_flags_s, s, 's', &1483 bc_dirichlet_l .OR. bc_radiation_l, &1484 bc_dirichlet_n .OR. bc_radiation_n, &1485 bc_dirichlet_r .OR. bc_radiation_r, &1474 CALL advec_s_ws( advc_flags_s, s, 's', & 1475 bc_dirichlet_l .OR. bc_radiation_l, & 1476 bc_dirichlet_n .OR. bc_radiation_n, & 1477 bc_dirichlet_r .OR. bc_radiation_r, & 1486 1478 bc_dirichlet_s .OR. bc_radiation_s ) 1487 1479 ELSE … … 1493 1485 ENDIF 1494 1486 1495 CALL diffusion_s( s, & 1496 surf_def_h(0)%ssws, surf_def_h(1)%ssws, & 1497 surf_def_h(2)%ssws, & 1498 surf_lsm_h%ssws, surf_usm_h%ssws, & 1499 surf_def_v(0)%ssws, surf_def_v(1)%ssws, & 1500 surf_def_v(2)%ssws, surf_def_v(3)%ssws, & 1501 surf_lsm_v(0)%ssws, surf_lsm_v(1)%ssws, & 1502 surf_lsm_v(2)%ssws, surf_lsm_v(3)%ssws, & 1503 surf_usm_v(0)%ssws, surf_usm_v(1)%ssws, & 1504 surf_usm_v(2)%ssws, surf_usm_v(3)%ssws ) 1487 CALL diffusion_s( s, surf_def_h(0)%ssws, surf_def_h(1)%ssws, surf_def_h(2)%ssws, & 1488 surf_lsm_h%ssws, surf_usm_h%ssws, surf_def_v(0)%ssws, surf_def_v(1)%ssws, & 1489 surf_def_v(2)%ssws, surf_def_v(3)%ssws, surf_lsm_v(0)%ssws, & 1490 surf_lsm_v(1)%ssws, surf_lsm_v(2)%ssws, surf_lsm_v(3)%ssws, & 1491 surf_usm_v(0)%ssws, surf_usm_v(1)%ssws, surf_usm_v(2)%ssws, & 1492 surf_usm_v(3)%ssws ) 1505 1493 1506 1494 ! … … 1521 1509 !-- If required compute influence of large-scale subsidence/ascent. 1522 1510 !-- Not implemented for scalars so far. 1523 IF ( large_scale_subsidence .AND. & 1524 .NOT. use_subsidence_tendencies .AND. & 1511 IF ( large_scale_subsidence .AND. .NOT. use_subsidence_tendencies .AND. & 1525 1512 .NOT. large_scale_forcing ) THEN 1526 1513 CALL subsidence( tend, s, s_init, 3 ) … … 1533 1520 DO i = nxl, nxr 1534 1521 DO j = nys, nyn 1535 !following directive is required to vectorize on Intel19 1522 ! 1523 !-- Following directive is required to vectorize on Intel19 1536 1524 !DIR$ IVDEP 1537 1525 DO k = nzb+1, nzt 1538 s_p(k,j,i) = s(k,j,i) + ( dt_3d * ( sbt * tend(k,j,i) + & 1539 tsc(3) * ts_m(k,j,i) ) & 1540 - tsc(5) * rdf_sc(k) * & 1541 ( s(k,j,i) - s_init(k) ) & 1542 ) & 1543 * MERGE( 1.0_wp, 0.0_wp, & 1544 BTEST( wall_flags_total_0(k,j,i), 0 ) & 1545 ) 1526 s_p(k,j,i) = s(k,j,i) & 1527 + ( dt_3d * ( sbt * tend(k,j,i) + tsc(3) * ts_m(k,j,i) ) & 1528 - tsc(5) * rdf_sc(k) * ( s(k,j,i) - s_init(k) ) ) & 1529 * MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_total_0(k,j,i), 0 ) ) 1546 1530 IF ( s_p(k,j,i) < 0.0_wp ) s_p(k,j,i) = 0.1_wp * s(k,j,i) 1547 1531 ENDDO … … 1560 1544 ENDDO 1561 1545 ENDDO 1562 ELSEIF ( intermediate_timestep_count < & 1563 intermediate_timestep_count_max ) THEN 1546 ELSEIF ( intermediate_timestep_count < intermediate_timestep_count_max ) THEN 1564 1547 DO i = nxl, nxr 1565 1548 DO j = nys, nyn 1566 1549 DO k = nzb+1, nzt 1567 ts_m(k,j,i) = -9.5625_wp * tend(k,j,i) & 1568 + 5.3125_wp * ts_m(k,j,i) 1550 ts_m(k,j,i) = - 9.5625_wp * tend(k,j,i) + 5.3125_wp * ts_m(k,j,i) 1569 1551 ENDDO 1570 1552 ENDDO
Note: See TracChangeset
for help on using the changeset viewer.
