Changeset 3182 for palm/trunk/SOURCE/large_scale_forcing_nudging_mod.f90
- Timestamp:
- Jul 27, 2018 1:36:03 PM (6 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/large_scale_forcing_nudging_mod.f90
r3049 r3182 20 20 ! Current revisions: 21 21 ! ------------------ 22 ! 22 ! * Adjustment to new Inifor version: 23 ! - No vertical interpolation/extrapolation of lateral boundary data required 24 ! any more (Inifor can treat grid stretching now 25 ! - Revise initialization in case of COSMO forcing 26 ! * Rename variables and subroutines for offline nesting 23 27 ! 24 28 ! Former revisions: … … 69 73 70 74 USE arrays_3d, & 71 ONLY: dzw, e, heatflux_input_conversion, pt, pt_init, q, q_init, s,&72 tend, u, u_init, ug, v, v_init, vg, w, w_subs,&75 ONLY: dzw, e, diss, heatflux_input_conversion, pt, pt_init, q, & 76 q_init, s, tend, u, u_init, ug, v, v_init, vg, w, w_subs, & 73 77 waterflux_input_conversion, zu, zw 74 78 75 79 USE control_parameters, & 76 ONLY: bc_lr, bc_ns, bc_pt_b, bc_q_b, constant_diffusion, & 80 ONLY: bc_dirichlet_l, bc_dirichlet_n, bc_dirichlet_r, bc_dirichlet_s, & 81 bc_lr, bc_ns, bc_pt_b, bc_q_b, constant_diffusion, & 77 82 constant_heatflux, constant_waterflux, & 78 data_output_pr, dt_3d, end_time, forcing, & 79 force_bound_l, force_bound_n, force_bound_r, force_bound_s, & 83 data_output_pr, dt_3d, end_time, & 80 84 humidity, initializing_actions, intermediate_timestep_count, & 81 85 ibc_pt_b, ibc_q_b, & 82 86 large_scale_forcing, large_scale_subsidence, lsf_surf, lsf_vert,& 83 lsf_exception, message_string, ne utral, nudging, passive_scalar,&84 pt_surface, ocean, q_surface, surface_heatflux,&85 surface_ pressure, surface_waterflux, topography,&86 use_subsidence_tendencies87 lsf_exception, message_string, nesting_offline, neutral, & 88 nudging, passive_scalar, pt_surface, ocean, q_surface, & 89 surface_heatflux, surface_pressure, surface_waterflux, & 90 topography, use_subsidence_tendencies 87 91 88 92 USE grid_variables 89 90 USE pegrid91 93 92 94 USE indices, & … … 96 98 USE kinds 97 99 100 USE netcdf_data_input_mod, & 101 ONLY: nest_offl 102 103 USE pegrid 104 98 105 USE surface_mod, & 99 106 ONLY: surf_def_h, surf_lsm_h, surf_usm_h … … 101 108 USE statistics, & 102 109 ONLY: hom, statistic_regions, sums_ls_l, weight_substep 103 104 USE netcdf_data_input_mod, &105 ONLY: force, netcdf_data_input_interpolate106 110 107 111 INTEGER(iwp) :: nlsf = 1000 !< maximum number of profiles in LSF_DATA (large scale forcing) … … 140 144 ! 141 145 !-- Public subroutines 142 PUBLIC ls_forcing_surf, ls_forcing_vert, ls_advec, lsf_init,&146 PUBLIC calc_tnudge, ls_forcing_surf, ls_forcing_vert, ls_advec, lsf_init, & 143 147 lsf_nudging_check_parameters, nudge_init, & 144 148 lsf_nudging_check_data_output_pr, lsf_nudging_header, & 145 calc_tnudge, nudge, nudge_ref, forcing_bc_mass_conservation, & 146 forcing_bc 149 lsf_nesting_offline, lsf_nesting_offline_mass_conservation, & 150 nudge, nudge_ref 151 147 152 ! 148 153 !-- Public variables 149 154 PUBLIC qsws_surf, shf_surf, td_lsa_lpt, td_lsa_q, td_sub_lpt, & 150 td_sub_q, time_vert , force155 td_sub_q, time_vert 151 156 152 157 … … 167 172 ! Description: 168 173 ! ------------ 169 !> @todo Missing subroutine description. 170 !------------------------------------------------------------------------------! 171 SUBROUTINE forcing_bc_mass_conservation 174 !> In this subroutine a constant mass within the model domain is guaranteed. 175 !> Larger-scale models may be based on a compressible equation system, which is 176 !> not consistent with PALMs incompressible equation system. In order to avoid 177 !> a decrease or increase of mass during the simulation, non-divergent flow 178 !> through the lateral and top boundaries is compensated by the vertical wind 179 !> component at the top boundary. 180 !------------------------------------------------------------------------------! 181 SUBROUTINE lsf_nesting_offline_mass_conservation 172 182 173 183 USE control_parameters, & … … 176 186 IMPLICIT NONE 177 187 178 INTEGER(iwp) :: i !< 179 INTEGER(iwp) :: j !< 180 INTEGER(iwp) :: k !< 181 182 REAL(wp) :: w_correct !<183 REAL(wp), DIMENSION(1:3) :: volume_flow_l !< 188 INTEGER(iwp) :: i !< grid index in x-direction 189 INTEGER(iwp) :: j !< grid index in y-direction 190 INTEGER(iwp) :: k !< grid index in z-direction 191 192 REAL(wp) :: w_correct !< vertical velocity increment required to compensate non-divergent flow through the boundaries 193 REAL(wp), DIMENSION(1:3) :: volume_flow_l !< local volume flow 184 194 185 195 volume_flow = 0.0_wp … … 188 198 d_area_t = 1.0_wp / ( ( nx + 1 ) * dx * ( ny + 1 ) * dy ) 189 199 190 IF ( force_bound_l ) THEN200 IF ( bc_dirichlet_l ) THEN 191 201 i = nxl 192 202 DO j = nys, nyn … … 198 208 ENDDO 199 209 ENDIF 200 IF ( force_bound_r ) THEN210 IF ( bc_dirichlet_r ) THEN 201 211 i = nxr+1 202 212 DO j = nys, nyn … … 208 218 ENDDO 209 219 ENDIF 210 IF ( force_bound_s ) THEN220 IF ( bc_dirichlet_s ) THEN 211 221 j = nys 212 222 DO i = nxl, nxr … … 218 228 ENDDO 219 229 ENDIF 220 IF ( force_bound_n ) THEN230 IF ( bc_dirichlet_n ) THEN 221 231 j = nyn+1 222 232 DO i = nxl, nxr … … 255 265 ENDDO 256 266 257 write(9,*) "w correction", w_correct 258 flush(9) 259 260 END SUBROUTINE forcing_bc_mass_conservation 267 END SUBROUTINE lsf_nesting_offline_mass_conservation 261 268 262 269 … … 264 271 ! Description: 265 272 ! ------------ 266 !> @todo Missing subroutine description. 267 !------------------------------------------------------------------------------! 268 SUBROUTINE forcing_bc 273 !> Set the lateral and top boundary conditions in case the PALM domain is 274 !> nested offline in a mesoscale model. 275 !------------------------------------------------------------------------------! 276 SUBROUTINE lsf_nesting_offline 269 277 270 278 USE control_parameters, & 271 ONLY: force_bound_l, force_bound_n, force_bound_r, force_bound_s, & 272 humidity, neutral, passive_scalar, simulated_time 273 274 USE netcdf_data_input_mod, & 275 ONLY: force 279 ONLY: bc_dirichlet_l, bc_dirichlet_n, bc_dirichlet_r, & 280 bc_dirichlet_s, humidity, neutral, passive_scalar, rans_mode,& 281 rans_tke_e, time_since_reference_point 276 282 277 283 IMPLICIT NONE … … 284 290 REAL(wp) :: ddt_lsf !< inverse value of time resolution of forcing data 285 291 REAL(wp) :: t_ref !< time past since last reference step 286 287 ! 288 !-- If required, interpolate and/or extrapolate data vertically. This is 289 !-- required as Inifor outputs only equidistant vertical data. 290 IF ( ANY( zu(1:nzt+1) /= force%zu_atmos(1:force%nzu) ) ) THEN 291 IF ( .NOT. force%interpolated ) THEN 292 293 DO t = 0, 1 294 IF ( force_bound_l ) THEN 295 CALL netcdf_data_input_interpolate( force%u_left(t,:,:), & 296 zu(nzb+1:nzt+1), & 297 force%zu_atmos ) 298 CALL netcdf_data_input_interpolate( force%v_left(t,:,:), & 299 zu(nzb+1:nzt+1), & 300 force%zu_atmos ) 301 CALL netcdf_data_input_interpolate( force%w_left(t,:,:), & 302 zw(nzb+1:nzt+1), & 303 force%zw_atmos ) 304 IF ( .NOT. neutral ) & 305 CALL netcdf_data_input_interpolate( force%pt_left(t,:,:),& 306 zu(nzb+1:nzt+1), & 307 force%zu_atmos ) 308 IF ( humidity ) & 309 CALL netcdf_data_input_interpolate( force%q_left(t,:,:), & 310 zu(nzb+1:nzt+1), & 311 force%zu_atmos ) 312 ENDIF 313 IF ( force_bound_r ) THEN 314 CALL netcdf_data_input_interpolate( force%u_right(t,:,:), & 315 zu(nzb+1:nzt+1), & 316 force%zu_atmos ) 317 CALL netcdf_data_input_interpolate( force%v_right(t,:,:), & 318 zu(nzb+1:nzt+1), & 319 force%zu_atmos ) 320 CALL netcdf_data_input_interpolate( force%w_right(t,:,:), & 321 zw(nzb+1:nzt+1), & 322 force%zw_atmos ) 323 IF ( .NOT. neutral ) & 324 CALL netcdf_data_input_interpolate( force%pt_right(t,:,:),& 325 zu(nzb+1:nzt+1), & 326 force%zu_atmos ) 327 IF ( humidity ) & 328 CALL netcdf_data_input_interpolate( force%q_right(t,:,:),& 329 zu(nzb+1:nzt+1), & 330 force%zu_atmos ) 331 ENDIF 332 IF ( force_bound_n ) THEN 333 CALL netcdf_data_input_interpolate( force%u_north(t,:,:), & 334 zu(nzb+1:nzt+1), & 335 force%zu_atmos ) 336 CALL netcdf_data_input_interpolate( force%v_north(t,:,:), & 337 zu(nzb+1:nzt+1), & 338 force%zu_atmos ) 339 CALL netcdf_data_input_interpolate( force%w_north(t,:,:), & 340 zw(nzb+1:nzt+1), & 341 force%zw_atmos ) 342 IF ( .NOT. neutral ) & 343 CALL netcdf_data_input_interpolate( force%pt_north(t,:,:),& 344 zu(nzb+1:nzt+1), & 345 force%zu_atmos ) 346 IF ( humidity ) & 347 CALL netcdf_data_input_interpolate( force%q_north(t,:,:),& 348 zu(nzb+1:nzt+1), & 349 force%zu_atmos ) 350 ENDIF 351 IF ( force_bound_s ) THEN 352 CALL netcdf_data_input_interpolate( force%u_south(t,:,:), & 353 zu(nzb+1:nzt+1), & 354 force%zu_atmos ) 355 CALL netcdf_data_input_interpolate( force%v_south(t,:,:), & 356 zu(nzb+1:nzt+1), & 357 force%zu_atmos ) 358 CALL netcdf_data_input_interpolate( force%w_south(t,:,:), & 359 zw(nzb+1:nzt+1), & 360 force%zw_atmos ) 361 IF ( .NOT. neutral ) & 362 CALL netcdf_data_input_interpolate( force%pt_south(t,:,:),& 363 zu(nzb+1:nzt+1), & 364 force%zu_atmos ) 365 IF ( humidity ) & 366 CALL netcdf_data_input_interpolate( force%q_south(t,:,:),& 367 zu(nzb+1:nzt+1), & 368 force%zu_atmos ) 369 ENDIF 370 ENDDO 371 ! 372 !-- Note, no interpolation of top boundary. Just use initial value. 373 !-- No physical meaningful extrapolation possible if only one layer is 374 !-- given. 375 376 force%interpolated = .TRUE. 377 ENDIF 378 ENDIF 379 292 380 293 ! 381 294 !-- Calculate time interval of forcing data 382 ddt_lsf = 1.0_wp / ( force%time(force%tind_p) - force%time(force%tind) ) 295 ddt_lsf = 1.0_wp / ( nest_offl%time(nest_offl%tind_p) - & 296 nest_offl%time(nest_offl%tind) ) 383 297 ! 384 298 !-- Calculate reziproke time past since last reference step. Please note, 385 !-- as simulated time is still not updated, the actual time here is 386 !-- simulated time + dt_3d 387 t_ref = simulated_time + dt_3d - force%time(force%tind) 388 389 IF ( force_bound_l ) THEN 390 391 DO j = nys, nyn 392 DO k = nzb+1, nzt+1 393 u(k,j,nxlg:nxl) = force%u_left(0,k,j) + ddt_lsf * t_ref * & 394 ( force%u_left(1,k,j) - force%u_left(0,k,j) ) * & 395 MERGE( 1.0_wp, 0.0_wp, & 396 BTEST( wall_flags_0(k,j,nxlg:nxl), 1 ) ) 397 ENDDO 398 ENDDO 299 !-- the time coordinate is still not updated, so that the actual time need 300 !-- to be incremented by dt_3d. Moreover, note that the simulation time 301 !-- passed since simulation start is time_since_reference_point, not 302 !-- simulated_time! 303 t_ref = time_since_reference_point + dt_3d - & 304 nest_offl%time(nest_offl%tind) 305 306 IF ( bc_dirichlet_l ) THEN 399 307 400 308 DO j = nys, nyn 401 309 DO k = nzb+1, nzt 402 w(k,j,nxlg:nxl-1) = force%w_left(0,k,j) + ddt_lsf * t_ref * & 403 ( force%w_left(1,k,j) - force%w_left(0,k,j) ) * & 404 MERGE( 1.0_wp, 0.0_wp, & 405 BTEST( wall_flags_0(k,j,nxlg:nxl-1), 3 ) ) 310 u(k,j,nxlg:nxl) = nest_offl%u_left(0,k,j) + ddt_lsf * t_ref * & 311 ( nest_offl%u_left(1,k,j) - nest_offl%u_left(0,k,j) ) * & 312 MERGE( 1.0_wp, 0.0_wp, & 313 BTEST( wall_flags_0(k,j,nxlg:nxl), 1 ) ) 314 ENDDO 315 ENDDO 316 317 DO j = nys, nyn 318 DO k = nzb+1, nzt-1 319 w(k,j,nxlg:nxl-1) = nest_offl%w_left(0,k,j) + ddt_lsf * t_ref *& 320 ( nest_offl%w_left(1,k,j) - nest_offl%w_left(0,k,j) ) *& 321 MERGE( 1.0_wp, 0.0_wp, & 322 BTEST( wall_flags_0(k,j,nxlg:nxl-1), 3 ) ) 406 323 ENDDO 407 324 ENDDO 408 325 409 326 DO j = nysv, nyn 410 DO k = nzb+1, nzt +1411 v(k,j,nxlg:nxl-1) = force%v_left(0,k,j) + ddt_lsf * t_ref *&412 ( force%v_left(1,k,j) - force%v_left(0,k,j) ) *&413 MERGE( 1.0_wp, 0.0_wp,&414 327 DO k = nzb+1, nzt 328 v(k,j,nxlg:nxl-1) = nest_offl%v_left(0,k,j) + ddt_lsf * t_ref *& 329 ( nest_offl%v_left(1,k,j) - nest_offl%v_left(0,k,j) ) *& 330 MERGE( 1.0_wp, 0.0_wp, & 331 BTEST( wall_flags_0(k,j,nxlg:nxl-1), 2 ) ) 415 332 ENDDO 416 333 ENDDO … … 418 335 IF ( .NOT. neutral ) THEN 419 336 DO j = nys, nyn 420 DO k = nzb+1, nzt +1421 pt(k,j,nxlg:nxl-1) = force%pt_left(0,k,j) + ddt_lsf *&422 t_ref *&423 ( force%pt_left(1,k,j) - force%pt_left(0,k,j) )337 DO k = nzb+1, nzt 338 pt(k,j,nxlg:nxl-1) = nest_offl%pt_left(0,k,j) + ddt_lsf * & 339 t_ref * & 340 ( nest_offl%pt_left(1,k,j) - nest_offl%pt_left(0,k,j) ) 424 341 425 342 ENDDO … … 429 346 IF ( humidity ) THEN 430 347 DO j = nys, nyn 431 DO k = nzb+1, nzt +1432 q(k,j,nxlg:nxl-1) = force%q_left(0,k,j) + ddt_lsf *&433 t_ref *&434 ( force%q_left(1,k,j) - force%q_left(0,k,j) )348 DO k = nzb+1, nzt 349 q(k,j,nxlg:nxl-1) = nest_offl%q_left(0,k,j) + ddt_lsf * & 350 t_ref * & 351 ( nest_offl%q_left(1,k,j) - nest_offl%q_left(0,k,j) ) 435 352 436 353 ENDDO … … 440 357 ENDIF 441 358 442 IF ( force_bound_r ) THEN 443 444 DO j = nys, nyn 445 DO k = nzb+1, nzt+1 446 u(k,j,nxr+1:nxrg) = force%u_right(0,k,j) + ddt_lsf * t_ref * & 447 ( force%u_right(1,k,j) - force%u_right(0,k,j) ) * & 448 MERGE( 1.0_wp, 0.0_wp, & 449 BTEST( wall_flags_0(k,j,nxr+1:nxrg), 1 ) ) 450 451 ENDDO 452 ENDDO 359 IF ( bc_dirichlet_r ) THEN 360 453 361 DO j = nys, nyn 454 362 DO k = nzb+1, nzt 455 w(k,j,nxr+1:nxrg) = force%w_right(0,k,j) + ddt_lsf * t_ref * & 456 ( force%w_right(1,k,j) - force%w_right(0,k,j) ) * & 457 MERGE( 1.0_wp, 0.0_wp, & 458 BTEST( wall_flags_0(k,j,nxr+1:nxrg), 3 ) ) 363 u(k,j,nxr+1:nxrg) = nest_offl%u_right(0,k,j) + ddt_lsf * t_ref *& 364 ( nest_offl%u_right(1,k,j) - nest_offl%u_right(0,k,j) ) *& 365 MERGE( 1.0_wp, 0.0_wp, & 366 BTEST( wall_flags_0(k,j,nxr+1:nxrg), 1 ) ) 367 368 ENDDO 369 ENDDO 370 DO j = nys, nyn 371 DO k = nzb+1, nzt-1 372 w(k,j,nxr+1:nxrg) = nest_offl%w_right(0,k,j) + ddt_lsf * t_ref *& 373 ( nest_offl%w_right(1,k,j) - nest_offl%w_right(0,k,j) ) *& 374 MERGE( 1.0_wp, 0.0_wp, & 375 BTEST( wall_flags_0(k,j,nxr+1:nxrg), 3 ) ) 459 376 ENDDO 460 377 ENDDO 461 378 462 379 DO j = nysv, nyn 463 DO k = nzb+1, nzt +1464 v(k,j,nxr+1:nxrg) = force%v_right(0,k,j) + ddt_lsf * t_ref *&465 ( force%v_right(1,k,j) - force%v_right(0,k,j) ) *&466 MERGE( 1.0_wp, 0.0_wp,&467 380 DO k = nzb+1, nzt 381 v(k,j,nxr+1:nxrg) = nest_offl%v_right(0,k,j) + ddt_lsf * t_ref *& 382 ( nest_offl%v_right(1,k,j) - nest_offl%v_right(0,k,j) ) *& 383 MERGE( 1.0_wp, 0.0_wp, & 384 BTEST( wall_flags_0(k,j,nxr+1:nxrg), 2 ) ) 468 385 ENDDO 469 386 ENDDO … … 471 388 IF ( .NOT. neutral ) THEN 472 389 DO j = nys, nyn 473 DO k = nzb+1, nzt +1474 pt(k,j,nxr+1:nxrg) = force%pt_right(0,k,j) + ddt_lsf *&475 t_ref *&476 ( force%pt_right(1,k,j) - force%pt_right(0,k,j) )390 DO k = nzb+1, nzt 391 pt(k,j,nxr+1:nxrg) = nest_offl%pt_right(0,k,j) + ddt_lsf * & 392 t_ref * & 393 ( nest_offl%pt_right(1,k,j) - nest_offl%pt_right(0,k,j) ) 477 394 478 395 ENDDO … … 482 399 IF ( humidity ) THEN 483 400 DO j = nys, nyn 484 DO k = nzb+1, nzt +1485 q(k,j,nxr+1:nxrg) = force%q_right(0,k,j) + ddt_lsf *&486 t_ref *&487 ( force%q_right(1,k,j) - force%q_right(0,k,j) )401 DO k = nzb+1, nzt 402 q(k,j,nxr+1:nxrg) = nest_offl%q_right(0,k,j) + ddt_lsf * & 403 t_ref * & 404 ( nest_offl%q_right(1,k,j) - nest_offl%q_right(0,k,j) ) 488 405 489 406 ENDDO … … 493 410 ENDIF 494 411 495 IF ( force_bound_s ) THEN 496 497 DO i = nxl, nxr 498 DO k = nzb+1, nzt+1 499 v(k,nysg:nys,i) = force%v_south(0,k,i) + ddt_lsf * t_ref * & 500 ( force%v_south(1,k,i) - force%v_south(0,k,i) ) * & 501 MERGE( 1.0_wp, 0.0_wp, & 502 BTEST( wall_flags_0(k,nysg:nys,i), 2 ) ) 503 ENDDO 504 ENDDO 412 IF ( bc_dirichlet_s ) THEN 505 413 506 414 DO i = nxl, nxr 507 415 DO k = nzb+1, nzt 508 w(k,nysg:nys-1,i) = force%w_south(0,k,i) + ddt_lsf * t_ref * & 509 ( force%w_south(1,k,i) - force%w_south(0,k,i) ) * & 510 MERGE( 1.0_wp, 0.0_wp, & 416 v(k,nysg:nys,i) = nest_offl%v_south(0,k,i) + ddt_lsf * t_ref *& 417 ( nest_offl%v_south(1,k,i) - nest_offl%v_south(0,k,i) ) *& 418 MERGE( 1.0_wp, 0.0_wp, & 419 BTEST( wall_flags_0(k,nysg:nys,i), 2 ) ) 420 ENDDO 421 ENDDO 422 423 DO i = nxl, nxr 424 DO k = nzb+1, nzt-1 425 w(k,nysg:nys-1,i) = nest_offl%w_south(0,k,i) + ddt_lsf * t_ref *& 426 ( nest_offl%w_south(1,k,i) - nest_offl%w_south(0,k,i) ) *& 427 MERGE( 1.0_wp, 0.0_wp, & 511 428 BTEST( wall_flags_0(k,nysg:nys-1,i), 3 ) ) 512 429 ENDDO … … 514 431 515 432 DO i = nxlu, nxr 516 DO k = nzb+1, nzt +1517 u(k,nysg:nys-1,i) = force%u_south(0,k,i) + ddt_lsf * t_ref *&518 ( force%u_south(1,k,i) - force%u_south(0,k,i) ) *&519 MERGE( 1.0_wp, 0.0_wp, &433 DO k = nzb+1, nzt 434 u(k,nysg:nys-1,i) = nest_offl%u_south(0,k,i) + ddt_lsf * t_ref *& 435 ( nest_offl%u_south(1,k,i) - nest_offl%u_south(0,k,i) ) *& 436 MERGE( 1.0_wp, 0.0_wp, & 520 437 BTEST( wall_flags_0(k,nysg:nys-1,i), 1 ) ) 521 438 ENDDO … … 524 441 IF ( .NOT. neutral ) THEN 525 442 DO i = nxl, nxr 526 DO k = nzb+1, nzt +1527 pt(k,nysg:nys-1,i) = force%pt_south(0,k,i) + ddt_lsf *&528 t_ref *&529 ( force%pt_south(1,k,i) - force%pt_south(0,k,i) )443 DO k = nzb+1, nzt 444 pt(k,nysg:nys-1,i) = nest_offl%pt_south(0,k,i) + ddt_lsf * & 445 t_ref * & 446 ( nest_offl%pt_south(1,k,i) - nest_offl%pt_south(0,k,i) ) 530 447 531 448 ENDDO … … 535 452 IF ( humidity ) THEN 536 453 DO i = nxl, nxr 537 DO k = nzb+1, nzt +1538 q(k,nysg:nys-1,i) = force%q_south(0,k,i) + ddt_lsf *&539 t_ref *&540 ( force%q_south(1,k,i) - force%q_south(0,k,i) )454 DO k = nzb+1, nzt 455 q(k,nysg:nys-1,i) = nest_offl%q_south(0,k,i) + ddt_lsf * & 456 t_ref * & 457 ( nest_offl%q_south(1,k,i) - nest_offl%q_south(0,k,i) ) 541 458 542 459 ENDDO … … 546 463 ENDIF 547 464 548 IF ( force_bound_n ) THEN 549 550 DO i = nxl, nxr 551 DO k = nzb+1, nzt+1 552 v(k,nyn+1:nyng,i) = force%v_north(0,k,i) + ddt_lsf * t_ref * & 553 ( force%v_north(1,k,i) - force%v_north(0,k,i) ) * & 554 MERGE( 1.0_wp, 0.0_wp, & 555 BTEST( wall_flags_0(k,nyn+1:nyng,i), 2 ) ) 556 ENDDO 557 ENDDO 465 IF ( bc_dirichlet_n ) THEN 466 558 467 DO i = nxl, nxr 559 468 DO k = nzb+1, nzt 560 w(k,nyn+1:nyng,i) = force%w_north(0,k,i) + ddt_lsf * t_ref * & 561 ( force%w_north(1,k,i) - force%w_north(0,k,i) ) * & 562 MERGE( 1.0_wp, 0.0_wp, & 469 v(k,nyn+1:nyng,i) = nest_offl%v_north(0,k,i) + ddt_lsf * t_ref *& 470 ( nest_offl%v_north(1,k,i) - nest_offl%v_north(0,k,i) ) *& 471 MERGE( 1.0_wp, 0.0_wp, & 472 BTEST( wall_flags_0(k,nyn+1:nyng,i), 2 ) ) 473 ENDDO 474 ENDDO 475 DO i = nxl, nxr 476 DO k = nzb+1, nzt-1 477 w(k,nyn+1:nyng,i) = nest_offl%w_north(0,k,i) + ddt_lsf * t_ref *& 478 ( nest_offl%w_north(1,k,i) - nest_offl%w_north(0,k,i) ) *& 479 MERGE( 1.0_wp, 0.0_wp, & 563 480 BTEST( wall_flags_0(k,nyn+1:nyng,i), 3 ) ) 564 481 ENDDO … … 566 483 567 484 DO i = nxlu, nxr 568 DO k = nzb+1, nzt +1569 u(k,nyn+1:nyng,i) = force%u_north(0,k,i) + ddt_lsf * t_ref *&570 ( force%u_north(1,k,i) - force%u_north(0,k,i) ) *&571 MERGE( 1.0_wp, 0.0_wp, &485 DO k = nzb+1, nzt 486 u(k,nyn+1:nyng,i) = nest_offl%u_north(0,k,i) + ddt_lsf * t_ref *& 487 ( nest_offl%u_north(1,k,i) - nest_offl%u_north(0,k,i) ) *& 488 MERGE( 1.0_wp, 0.0_wp, & 572 489 BTEST( wall_flags_0(k,nyn+1:nyng,i), 1 ) ) 573 490 … … 577 494 IF ( .NOT. neutral ) THEN 578 495 DO i = nxl, nxr 579 DO k = nzb+1, nzt +1580 pt(k,nyn+1:nyng,i) = force%pt_north(0,k,i) + ddt_lsf *&581 t_ref *&582 ( force%pt_north(1,k,i) - force%pt_north(0,k,i) )496 DO k = nzb+1, nzt 497 pt(k,nyn+1:nyng,i) = nest_offl%pt_north(0,k,i) + ddt_lsf * & 498 t_ref * & 499 ( nest_offl%pt_north(1,k,i) - nest_offl%pt_north(0,k,i) ) 583 500 584 501 ENDDO … … 588 505 IF ( humidity ) THEN 589 506 DO i = nxl, nxr 590 DO k = nzb+1, nzt +1591 q(k,nyn+1:nyng,i) = force%q_north(0,k,i) + ddt_lsf *&592 t_ref *&593 ( force%q_north(1,k,i) - force%q_north(0,k,i) )507 DO k = nzb+1, nzt 508 q(k,nyn+1:nyng,i) = nest_offl%q_north(0,k,i) + ddt_lsf * & 509 t_ref * & 510 ( nest_offl%q_north(1,k,i) - nest_offl%q_north(0,k,i) ) 594 511 595 512 ENDDO … … 600 517 ! 601 518 !-- Top boundary. 602 !-- Please note, only map Inifor data on model top in case the numeric is603 !-- identical to the Inifor grid. At the top boundary an extrapolation is604 !-- not possible.605 519 DO i = nxlu, nxr 606 520 DO j = nys, nyn 607 u(nzt+1,j,i) = force%u_top(0,j,i) + ddt_lsf * t_ref *&608 ( force%u_top(1,j,i) - force%u_top(0,j,i) ) *&609 MERGE( 1.0_wp, 0.0_wp, &521 u(nzt+1,j,i) = nest_offl%u_top(0,j,i) + ddt_lsf * t_ref * & 522 ( nest_offl%u_top(1,j,i) - nest_offl%u_top(0,j,i) ) * & 523 MERGE( 1.0_wp, 0.0_wp, & 610 524 BTEST( wall_flags_0(nzt+1,j,i), 1 ) ) 611 525 ENDDO … … 614 528 DO i = nxl, nxr 615 529 DO j = nysv, nyn 616 v(nzt+1,j,i) = force%v_top(0,j,i) + ddt_lsf * t_ref *&617 ( force%v_top(1,j,i) - force%v_top(0,j,i) ) *&530 v(nzt+1,j,i) = nest_offl%v_top(0,j,i) + ddt_lsf * t_ref * & 531 ( nest_offl%v_top(1,j,i) - nest_offl%v_top(0,j,i) ) * & 618 532 MERGE( 1.0_wp, 0.0_wp, & 619 533 BTEST( wall_flags_0(nzt+1,j,i), 2 ) ) … … 623 537 DO i = nxl, nxr 624 538 DO j = nys, nyn 625 w(nzt:nzt+1,j,i) = force%w_top(0,j,i) + ddt_lsf * t_ref *&626 ( force%w_top(1,j,i) - force%w_top(0,j,i) ) *&539 w(nzt:nzt+1,j,i) = nest_offl%w_top(0,j,i) + ddt_lsf * t_ref * & 540 ( nest_offl%w_top(1,j,i) - nest_offl%w_top(0,j,i) ) * & 627 541 MERGE( 1.0_wp, 0.0_wp, & 628 542 BTEST( wall_flags_0(nzt:nzt+1,j,i), 3 ) ) … … 634 548 DO i = nxl, nxr 635 549 DO j = nys, nyn 636 pt(nzt+1,j,i) = force%pt_top(0,j,i) + ddt_lsf * t_ref *&637 ( force%pt_top(1,j,i) - force%pt_top(0,j,i) )550 pt(nzt+1,j,i) = nest_offl%pt_top(0,j,i) + ddt_lsf * t_ref * & 551 ( nest_offl%pt_top(1,j,i) - nest_offl%pt_top(0,j,i) ) 638 552 ENDDO 639 553 ENDDO … … 643 557 DO i = nxl, nxr 644 558 DO j = nys, nyn 645 q(nzt+1,j,i) = force%q_top(0,j,i) + ddt_lsf * t_ref *&646 ( force%q_top(1,j,i) - force%q_top(0,j,i) )559 q(nzt+1,j,i) = nest_offl%q_top(0,j,i) + ddt_lsf * t_ref * & 560 ( nest_offl%q_top(1,j,i) - nest_offl%q_top(0,j,i) ) 647 561 ENDDO 648 562 ENDDO … … 651 565 !-- At the edges( left-south, left-north, right-south and right-north) set 652 566 !-- data on ghost points. 653 IF ( force_bound_l .AND. force_bound_s ) THEN567 IF ( bc_dirichlet_l .AND. bc_dirichlet_s ) THEN 654 568 DO i = 1, nbgp 655 569 u(:,nys-i,nxlg:nxl) = u(:,nys,nxlg:nxl) 656 570 w(:,nys-i,nxlg:nxl-1) = w(:,nys,nxlg:nxl-1) 657 571 IF ( .NOT. neutral ) pt(:,nys-i,nxlg:nxl-1) = pt(:,nys,nxlg:nxl-1) 658 IF ( humidity )q(:,nys-i,nxlg:nxl-1) = q(:,nys,nxlg:nxl-1)572 IF ( humidity ) q(:,nys-i,nxlg:nxl-1) = q(:,nys,nxlg:nxl-1) 659 573 ENDDO 660 574 DO i = 1, nbgp+1 … … 662 576 ENDDO 663 577 ENDIF 664 IF ( force_bound_l .AND. force_bound_n ) THEN578 IF ( bc_dirichlet_l .AND. bc_dirichlet_n ) THEN 665 579 DO i = 1, nbgp 666 580 u(:,nyn+i,nxlg:nxl) = u(:,nyn,nxlg:nxl) … … 668 582 w(:,nyn+i,nxlg:nxl-1) = w(:,nyn,nxlg:nxl-1) 669 583 IF ( .NOT. neutral ) pt(:,nyn+i,nxlg:nxl-1) = pt(:,nyn,nxlg:nxl-1) 670 IF ( humidity )q(:,nyn+i,nxlg:nxl-1) = q(:,nyn,nxlg:nxl-1)671 ENDDO 672 ENDIF 673 IF ( force_bound_r .AND. force_bound_s ) THEN584 IF ( humidity ) q(:,nyn+i,nxlg:nxl-1) = q(:,nyn,nxlg:nxl-1) 585 ENDDO 586 ENDIF 587 IF ( bc_dirichlet_r .AND. bc_dirichlet_s ) THEN 674 588 DO i = 1, nbgp 675 589 u(:,nys-i,nxr+1:nxrg) = u(:,nys,nxr+1:nxrg) 676 590 w(:,nys-i,nxr+1:nxrg) = w(:,nys,nxr+1:nxrg) 677 591 IF ( .NOT. neutral ) pt(:,nys-i,nxr+1:nxrg) = pt(:,nys,nxr+1:nxrg) 678 IF ( humidity )q(:,nys-i,nxr+1:nxrg) = q(:,nys,nxr+1:nxrg)592 IF ( humidity ) q(:,nys-i,nxr+1:nxrg) = q(:,nys,nxr+1:nxrg) 679 593 ENDDO 680 594 DO i = 1, nbgp+1 … … 682 596 ENDDO 683 597 ENDIF 684 IF ( force_bound_r .AND. force_bound_n ) THEN598 IF ( bc_dirichlet_r .AND. bc_dirichlet_n ) THEN 685 599 DO i = 1, nbgp 686 600 u(:,nyn+i,nxr+1:nxrg) = u(:,nyn,nxr+1:nxrg) … … 688 602 w(:,nyn+i,nxr+1:nxrg) = w(:,nyn,nxr+1:nxrg) 689 603 IF ( .NOT. neutral ) pt(:,nyn+i,nxr+1:nxrg) = pt(:,nyn,nxr+1:nxrg) 690 IF ( humidity ) q(:,nyn+i,nxr+1:nxrg) = q(:,nyn,nxr+1:nxrg) 691 ENDDO 692 ENDIF 693 ! 694 !-- Moreover, set Neumann boundary condition for subgrid-scale TKE and 695 !-- passive scalar 604 IF ( humidity ) q(:,nyn+i,nxr+1:nxrg) = q(:,nyn,nxr+1:nxrg) 605 ENDDO 606 ENDIF 607 ! 608 !-- Moreover, set Neumann boundary condition for subgrid-scale TKE, 609 !-- passive scalar, dissipation, and chemical species if required 610 IF ( rans_mode .AND. rans_tke_e ) THEN 611 IF ( bc_dirichlet_l ) diss(:,:,nxl-1) = diss(:,:,nxl) 612 IF ( bc_dirichlet_r ) diss(:,:,nxr+1) = diss(:,:,nxr) 613 IF ( bc_dirichlet_s ) diss(:,nys-1,:) = diss(:,nys,:) 614 IF ( bc_dirichlet_n ) diss(:,nyn+1,:) = diss(:,nyn,:) 615 ENDIF 696 616 IF ( .NOT. constant_diffusion ) THEN 697 IF ( force_bound_l ) e(:,:,nxl-1) = e(:,:,nxl)698 IF ( force_bound_r ) e(:,:,nxr+1) = e(:,:,nxr)699 IF ( force_bound_s ) e(:,nys-1,:) = e(:,nys,:)700 IF ( force_bound_n ) e(:,nyn+1,:) = e(:,nyn,:)617 IF ( bc_dirichlet_l ) e(:,:,nxl-1) = e(:,:,nxl) 618 IF ( bc_dirichlet_r ) e(:,:,nxr+1) = e(:,:,nxr) 619 IF ( bc_dirichlet_s ) e(:,nys-1,:) = e(:,nys,:) 620 IF ( bc_dirichlet_n ) e(:,nyn+1,:) = e(:,nyn,:) 701 621 e(nzt+1,:,:) = e(nzt,:,:) 702 622 ENDIF 703 623 IF ( passive_scalar ) THEN 704 IF ( force_bound_l ) s(:,:,nxl-1) = s(:,:,nxl) 705 IF ( force_bound_r ) s(:,:,nxr+1) = s(:,:,nxr) 706 IF ( force_bound_s ) s(:,nys-1,:) = s(:,nys,:) 707 IF ( force_bound_n ) s(:,nyn+1,:) = s(:,nyn,:) 708 ENDIF 624 IF ( bc_dirichlet_l ) s(:,:,nxl-1) = s(:,:,nxl) 625 IF ( bc_dirichlet_r ) s(:,:,nxr+1) = s(:,:,nxr) 626 IF ( bc_dirichlet_s ) s(:,nys-1,:) = s(:,nys,:) 627 IF ( bc_dirichlet_n ) s(:,nyn+1,:) = s(:,nyn,:) 628 ENDIF 629 709 630 710 631 … … 720 641 !-- treatment of fluxes. 721 642 !-- For the moment, comment this out! 722 ! surface_pressure = force%surface_pressure(force%tind) + &643 ! surface_pressure = nest_offl%surface_pressure(nest_offl%tind) + & 723 644 ! ddt_lsf * t_ref * & 724 ! ( force%surface_pressure(force%tind_p) &725 ! - force%surface_pressure(force%tind) )726 727 END SUBROUTINE forcing_bc645 ! ( nest_offl%surface_pressure(nest_offl%tind_p) & 646 ! - nest_offl%surface_pressure(nest_offl%tind) ) 647 648 END SUBROUTINE lsf_nesting_offline 728 649 729 650 !------------------------------------------------------------------------------! … … 1040 961 REAL(wp) :: r_dummy !< 1041 962 1042 IF ( forcing) THEN963 IF ( nesting_offline ) THEN 1043 964 ! 1044 965 !-- Allocate arrays for geostrophic wind components. Arrays will … … 1047 968 !-- case of cyclic boundary conditions. 1048 969 IF ( bc_lr_cyc .AND. bc_ns_cyc ) THEN 1049 ALLOCATE( force%ug(0:1,nzb:nzt+1) )1050 ALLOCATE( force%vg(0:1,nzb:nzt+1) )970 ALLOCATE( nest_offl%ug(0:1,nzb:nzt+1) ) 971 ALLOCATE( nest_offl%vg(0:1,nzb:nzt+1) ) 1051 972 ENDIF 1052 973 ! 1053 974 !-- Allocate arrays for reading boundary values. Arrays will incorporate 2 1054 975 !-- time levels in order to interpolate in between. 1055 IF ( force_bound_l ) THEN1056 ALLOCATE( force%u_left(0:1,nzb+1:nzt+1,nys:nyn) )1057 ALLOCATE( force%v_left(0:1,nzb+1:nzt+1,nysv:nyn) )1058 ALLOCATE( force%w_left(0:1,nzb+1:nzt,nys:nyn))1059 IF ( humidity ) ALLOCATE( force%q_left(0:1,nzb+1:nzt+1,nys:nyn) )1060 IF ( .NOT. neutral ) ALLOCATE( force%pt_left(0:1,nzb+1:nzt+1,nys:nyn) )1061 ENDIF 1062 IF ( force_bound_r ) THEN1063 ALLOCATE( force%u_right(0:1,nzb+1:nzt+1,nys:nyn) )1064 ALLOCATE( force%v_right(0:1,nzb+1:nzt+1,nysv:nyn) )1065 ALLOCATE( force%w_right(0:1,nzb+1:nzt,nys:nyn))1066 IF ( humidity ) ALLOCATE( force%q_right(0:1,nzb+1:nzt+1,nys:nyn) )1067 IF ( .NOT. neutral ) ALLOCATE( force%pt_right(0:1,nzb+1:nzt+1,nys:nyn) )1068 ENDIF 1069 IF ( force_bound_n ) THEN1070 ALLOCATE( force%u_north(0:1,nzb+1:nzt+1,nxlu:nxr) )1071 ALLOCATE( force%v_north(0:1,nzb+1:nzt+1,nxl:nxr) )1072 ALLOCATE( force%w_north(0:1,nzb+1:nzt,nxl:nxr))1073 IF ( humidity ) ALLOCATE( force%q_north(0:1,nzb+1:nzt+1,nxl:nxr) )1074 IF ( .NOT. neutral ) ALLOCATE( force%pt_north(0:1,nzb+1:nzt+1,nxl:nxr) )1075 ENDIF 1076 IF ( force_bound_s ) THEN1077 ALLOCATE( force%u_south(0:1,nzb+1:nzt+1,nxlu:nxr) )1078 ALLOCATE( force%v_south(0:1,nzb+1:nzt+1,nxl:nxr) )1079 ALLOCATE( force%w_south(0:1,nzb+1:nzt,nxl:nxr) )1080 IF ( humidity ) ALLOCATE( force%q_south(0:1,nzb+1:nzt+1,nxl:nxr) )1081 IF ( .NOT. neutral ) ALLOCATE( force%pt_south(0:1,nzb+1:nzt+1,nxl:nxr) )976 IF ( bc_dirichlet_l ) THEN 977 ALLOCATE( nest_offl%u_left(0:1,nzb+1:nzt,nys:nyn) ) 978 ALLOCATE( nest_offl%v_left(0:1,nzb+1:nzt,nysv:nyn) ) 979 ALLOCATE( nest_offl%w_left(0:1,nzb+1:nzt-1,nys:nyn) ) 980 IF ( humidity ) ALLOCATE( nest_offl%q_left(0:1,nzb+1:nzt,nys:nyn) ) 981 IF ( .NOT. neutral ) ALLOCATE( nest_offl%pt_left(0:1,nzb+1:nzt,nys:nyn) ) 982 ENDIF 983 IF ( bc_dirichlet_r ) THEN 984 ALLOCATE( nest_offl%u_right(0:1,nzb+1:nzt,nys:nyn) ) 985 ALLOCATE( nest_offl%v_right(0:1,nzb+1:nzt,nysv:nyn) ) 986 ALLOCATE( nest_offl%w_right(0:1,nzb+1:nzt-1,nys:nyn) ) 987 IF ( humidity ) ALLOCATE( nest_offl%q_right(0:1,nzb+1:nzt,nys:nyn) ) 988 IF ( .NOT. neutral ) ALLOCATE( nest_offl%pt_right(0:1,nzb+1:nzt,nys:nyn) ) 989 ENDIF 990 IF ( bc_dirichlet_n ) THEN 991 ALLOCATE( nest_offl%u_north(0:1,nzb+1:nzt,nxlu:nxr) ) 992 ALLOCATE( nest_offl%v_north(0:1,nzb+1:nzt,nxl:nxr) ) 993 ALLOCATE( nest_offl%w_north(0:1,nzb+1:nzt-1,nxl:nxr) ) 994 IF ( humidity ) ALLOCATE( nest_offl%q_north(0:1,nzb+1:nzt,nxl:nxr) ) 995 IF ( .NOT. neutral ) ALLOCATE( nest_offl%pt_north(0:1,nzb+1:nzt,nxl:nxr) ) 996 ENDIF 997 IF ( bc_dirichlet_s ) THEN 998 ALLOCATE( nest_offl%u_south(0:1,nzb+1:nzt,nxlu:nxr) ) 999 ALLOCATE( nest_offl%v_south(0:1,nzb+1:nzt,nxl:nxr) ) 1000 ALLOCATE( nest_offl%w_south(0:1,nzb+1:nzt-1,nxl:nxr) ) 1001 IF ( humidity ) ALLOCATE( nest_offl%q_south(0:1,nzb+1:nzt,nxl:nxr) ) 1002 IF ( .NOT. neutral ) ALLOCATE( nest_offl%pt_south(0:1,nzb+1:nzt,nxl:nxr) ) 1082 1003 ENDIF 1083 1004 1084 ALLOCATE( force%u_top(0:1,nys:nyn,nxlu:nxr) ) 1085 ALLOCATE( force%v_top(0:1,nysv:nyn,nxl:nxr) ) 1086 ALLOCATE( force%w_top(0:1,nys:nyn,nxl:nxr) ) 1087 IF ( humidity ) ALLOCATE( force%q_top(0:1,nys:nyn,nxl:nxr) ) 1088 IF ( .NOT. neutral ) ALLOCATE( force%pt_top(0:1,nys:nyn,nxl:nxr) ) 1089 1090 ! 1091 !-- Initial call of input. Time array, initial 3D data of u, v, w, 1092 !-- potential temperature, as well as mixing ratio, will be read. 1093 !-- Moreover, data at lateral and top boundary will be read. 1005 ALLOCATE( nest_offl%u_top(0:1,nys:nyn,nxlu:nxr) ) 1006 ALLOCATE( nest_offl%v_top(0:1,nysv:nyn,nxl:nxr) ) 1007 ALLOCATE( nest_offl%w_top(0:1,nys:nyn,nxl:nxr) ) 1008 IF ( humidity ) ALLOCATE( nest_offl%q_top(0:1,nys:nyn,nxl:nxr) ) 1009 IF ( .NOT. neutral ) ALLOCATE( nest_offl%pt_top(0:1,nys:nyn,nxl:nxr) ) 1010 1011 ! 1012 !-- Read COSMO data at lateral and top boundaries 1094 1013 CALL netcdf_data_input_lsf 1095 1014 ! 1096 !-- Please note, at the moment INIFOR assumes only an equidistant vertical 1097 !-- grid. In case of vertical grid stretching, input of inital 3D data 1098 !-- need to be inter- and/or extrapolated. 1099 !-- Therefore, check if zw grid on file is identical to numeric zw grid. 1100 IF ( ANY( zu(1:nzt+1) /= force%zu_atmos(1:force%nzu) ) ) THEN 1101 ! 1102 !-- Also data at the boundaries need to be inter/extrapolated at both 1103 !-- time levels 1104 DO t = 0, 1 1105 IF ( force_bound_l ) THEN 1106 CALL netcdf_data_input_interpolate( force%u_left(t,:,:), & 1107 zu(1:nzt+1), & 1108 force%zu_atmos ) 1109 CALL netcdf_data_input_interpolate( force%v_left(t,:,:), & 1110 zu(1:nzt+1), & 1111 force%zu_atmos ) 1112 CALL netcdf_data_input_interpolate( force%w_left(t,:,:), & 1113 zw(1:nzt+1), & 1114 force%zw_atmos ) 1115 IF ( .NOT. neutral ) & 1116 CALL netcdf_data_input_interpolate( force%pt_left(t,:,:),& 1117 zu(1:nzt+1), & 1118 force%zu_atmos ) 1119 IF ( humidity ) & 1120 CALL netcdf_data_input_interpolate( force%q_left(t,:,:), & 1121 zu(1:nzt+1), & 1122 force%zu_atmos ) 1123 ENDIF 1124 IF ( force_bound_r ) THEN 1125 CALL netcdf_data_input_interpolate( force%u_right(t,:,:), & 1126 zu(1:nzt+1), & 1127 force%zu_atmos ) 1128 CALL netcdf_data_input_interpolate( force%v_right(t,:,:), & 1129 zu(1:nzt+1), & 1130 force%zu_atmos ) 1131 CALL netcdf_data_input_interpolate( force%w_right(t,:,:), & 1132 zw(1:nzt+1), & 1133 force%zw_atmos ) 1134 IF ( .NOT. neutral ) & 1135 CALL netcdf_data_input_interpolate( force%pt_right(t,:,:),& 1136 zu(1:nzt+1), & 1137 force%zu_atmos ) 1138 IF ( humidity ) & 1139 CALL netcdf_data_input_interpolate( force%q_right(t,:,:),& 1140 zu(1:nzt+1), & 1141 force%zu_atmos ) 1142 ENDIF 1143 IF ( force_bound_n ) THEN 1144 CALL netcdf_data_input_interpolate( force%u_north(t,:,:), & 1145 zu(1:nzt+1), & 1146 force%zu_atmos ) 1147 CALL netcdf_data_input_interpolate( force%v_north(t,:,:), & 1148 zu(1:nzt+1), & 1149 force%zu_atmos ) 1150 CALL netcdf_data_input_interpolate( force%w_north(t,:,:), & 1151 zw(1:nzt+1), & 1152 force%zw_atmos ) 1153 IF ( .NOT. neutral ) & 1154 CALL netcdf_data_input_interpolate( force%pt_north(t,:,:),& 1155 zu(1:nzt+1), & 1156 force%zu_atmos ) 1157 IF ( humidity ) & 1158 CALL netcdf_data_input_interpolate( force%q_north(t,:,:),& 1159 zu(1:nzt+1), & 1160 force%zu_atmos ) 1161 ENDIF 1162 IF ( force_bound_s ) THEN 1163 CALL netcdf_data_input_interpolate( force%u_south(t,:,:), & 1164 zu(1:nzt+1), & 1165 force%zu_atmos ) 1166 CALL netcdf_data_input_interpolate( force%v_south(t,:,:), & 1167 zu(1:nzt+1), & 1168 force%zu_atmos ) 1169 CALL netcdf_data_input_interpolate( force%w_south(t,:,:), & 1170 zw(1:nzt+1), & 1171 force%zw_atmos ) 1172 IF ( .NOT. neutral ) & 1173 CALL netcdf_data_input_interpolate( force%pt_south(t,:,:),& 1174 zu(1:nzt+1), & 1175 force%zu_atmos ) 1176 IF ( humidity ) & 1177 CALL netcdf_data_input_interpolate( force%q_south(t,:,:),& 1178 zu(1:nzt+1), & 1179 force%zu_atmos ) 1180 ENDIF 1181 ENDDO 1182 ENDIF 1183 1184 ! 1185 !-- Exchange ghost points 1186 CALL exchange_horiz( u, nbgp ) 1187 CALL exchange_horiz( v, nbgp ) 1188 CALL exchange_horiz( w, nbgp ) 1189 IF ( .NOT. neutral ) CALL exchange_horiz( pt, nbgp ) 1190 IF ( humidity ) CALL exchange_horiz( q, nbgp ) 1191 ! 1192 !-- At lateral boundaries, set also initial boundary conditions 1193 IF ( force_bound_l ) THEN 1194 u(:,:,nxl) = u(:,:,nxlu) 1195 v(:,:,nxl-1) = v(:,:,nxl) 1196 w(:,:,nxl-1) = w(:,:,nxl) 1197 IF ( .NOT. neutral ) pt(:,:,nxl-1) = pt(:,:,nxl) 1198 IF ( humidity ) q(:,:,nxl-1) = q(:,:,nxl) 1199 ENDIF 1200 IF ( force_bound_r ) THEN 1201 u(:,:,nxr+1) = u(:,:,nxr) 1202 v(:,:,nxr+1) = v(:,:,nxr) 1203 w(:,:,nxr+1) = w(:,:,nxr) 1204 IF ( .NOT. neutral ) pt(:,:,nxr+1) = pt(:,:,nxr) 1205 IF ( humidity ) q(:,:,nxr+1) = q(:,:,nxr) 1206 ENDIF 1207 IF ( force_bound_s ) THEN 1208 u(:,nys-1,:) = u(:,nys,:) 1209 v(:,nys,:) = v(:,nysv,:) 1210 w(:,nys-1,:) = w(:,nys,:) 1211 IF ( .NOT. neutral ) pt(:,nys-1,:) = pt(:,nys,:) 1212 IF ( humidity ) q(:,nys-1,:) = q(:,nys,:) 1213 ENDIF 1214 IF ( force_bound_n ) THEN 1215 u(:,nyn+1,:) = u(:,nyn,:) 1216 v(:,nyn+1,:) = v(:,nyn,:) 1217 w(:,nyn+1,:) = w(:,nyn,:) 1218 IF ( .NOT. neutral ) pt(:,nyn+1,:) = pt(:,nyn,:) 1219 IF ( humidity ) q(:,nyn+1,:) = q(:,nyn,:) 1220 ENDIF 1221 1015 !-- Write COSMO data at lateral and top boundaries 1016 CALL lsf_nesting_offline 1222 1017 ! 1223 1018 !-- After 3D data is initialized, ensure mass conservation 1224 CALL forcing_bc_mass_conservation1019 CALL lsf_nesting_offline_mass_conservation 1225 1020 ! 1226 1021 !-- Initialize surface pressure. Please note, time-dependent surface … … 1228 1023 !-- treatment of fluxes. 1229 1024 !-- For the moment, comment this out! 1230 ! surface_pressure = force%surface_pressure(0)1025 ! surface_pressure = nest_offl%surface_pressure(0) 1231 1026 1232 1027 ELSE
Note: See TracChangeset
for help on using the changeset viewer.