- Timestamp:
- Feb 21, 2023, 3:26:41 PM (17 months ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
LMDZ6/branches/LMDZ_ECRad/libf/phylmd/ecrad/radiation_flux.F90
r3908 r4444 100 100 end type flux_type 101 101 102 ! Added for DWD (2020) 103 #ifdef __SX__ 104 logical, parameter :: use_indexed_sum_vec = .true. 105 #else 106 logical, parameter :: use_indexed_sum_vec = .false. 107 #endif 108 102 109 contains 103 110 … … 132 139 if (config%n_spec_lw == 0) then 133 140 write(nulerr,'(a)') '*** Error: number of LW spectral points to save not yet defined ' & 134 & // 'so cannot allocate dspectral flux arrays'141 & // 'so cannot allocate spectral flux arrays' 135 142 call radiation_abort() 136 143 end if … … 321 328 end if 322 329 330 if (allocated(this%lw_dn_surf_g)) deallocate(this%lw_dn_surf_g) 331 if (allocated(this%lw_dn_surf_clear_g)) deallocate(this%lw_dn_surf_clear_g) 332 if (allocated(this%sw_dn_diffuse_surf_g)) deallocate(this%sw_dn_diffuse_surf_g) 333 if (allocated(this%sw_dn_direct_surf_g)) deallocate(this%sw_dn_direct_surf_g) 334 if (allocated(this%sw_dn_diffuse_surf_clear_g)) deallocate(this%sw_dn_diffuse_surf_clear_g) 335 if (allocated(this%sw_dn_direct_surf_clear_g)) deallocate(this%sw_dn_direct_surf_clear_g) 336 323 337 if (lhook) call dr_hook('radiation_flux:deallocate',1,hook_handle) 324 338 … … 349 363 if (config%do_sw .and. config%do_surface_sw_spectral_flux) then 350 364 351 do jcol = istartcol,iendcol 352 call indexed_sum(this%sw_dn_direct_surf_g(:,jcol), & 353 & config%i_band_from_reordered_g_sw, & 354 & this%sw_dn_direct_surf_band(:,jcol)) 355 call indexed_sum(this%sw_dn_diffuse_surf_g(:,jcol), & 356 & config%i_band_from_reordered_g_sw, & 357 & this%sw_dn_surf_band(:,jcol)) 358 this%sw_dn_surf_band(:,jcol) & 359 & = this%sw_dn_surf_band(:,jcol) & 360 & + this%sw_dn_direct_surf_band(:,jcol) 361 end do 365 if (use_indexed_sum_vec) then 366 call indexed_sum_vec(this%sw_dn_direct_surf_g, & 367 & config%i_band_from_reordered_g_sw, & 368 & this%sw_dn_direct_surf_band, istartcol, iendcol) 369 call indexed_sum_vec(this%sw_dn_diffuse_surf_g, & 370 & config%i_band_from_reordered_g_sw, & 371 & this%sw_dn_surf_band, istartcol, iendcol) 372 do jcol = istartcol,iendcol 373 this%sw_dn_surf_band(:,jcol) & 374 & = this%sw_dn_surf_band(:,jcol) & 375 & + this%sw_dn_direct_surf_band(:,jcol) 376 end do 377 else 378 do jcol = istartcol,iendcol 379 call indexed_sum(this%sw_dn_direct_surf_g(:,jcol), & 380 & config%i_band_from_reordered_g_sw, & 381 & this%sw_dn_direct_surf_band(:,jcol)) 382 call indexed_sum(this%sw_dn_diffuse_surf_g(:,jcol), & 383 & config%i_band_from_reordered_g_sw, & 384 & this%sw_dn_surf_band(:,jcol)) 385 this%sw_dn_surf_band(:,jcol) & 386 & = this%sw_dn_surf_band(:,jcol) & 387 & + this%sw_dn_direct_surf_band(:,jcol) 388 end do 389 end if 362 390 363 391 if (config%do_clear) then 364 do jcol = istartcol,iendcol 365 call indexed_sum(this%sw_dn_direct_surf_clear_g(:,jcol), & 366 & config%i_band_from_reordered_g_sw, & 367 & this%sw_dn_direct_surf_clear_band(:,jcol)) 368 call indexed_sum(this%sw_dn_diffuse_surf_clear_g(:,jcol), & 369 & config%i_band_from_reordered_g_sw, & 370 & this%sw_dn_surf_clear_band(:,jcol)) 371 this%sw_dn_surf_clear_band(:,jcol) & 372 & = this%sw_dn_surf_clear_band(:,jcol) & 373 & + this%sw_dn_direct_surf_clear_band(:,jcol) 374 end do 392 if (use_indexed_sum_vec) then 393 call indexed_sum_vec(this%sw_dn_direct_surf_clear_g, & 394 & config%i_band_from_reordered_g_sw, & 395 & this%sw_dn_direct_surf_clear_band, istartcol, iendcol) 396 call indexed_sum_vec(this%sw_dn_diffuse_surf_clear_g, & 397 & config%i_band_from_reordered_g_sw, & 398 & this%sw_dn_surf_clear_band, istartcol, iendcol) 399 do jcol = istartcol,iendcol 400 this%sw_dn_surf_clear_band(:,jcol) & 401 & = this%sw_dn_surf_clear_band(:,jcol) & 402 & + this%sw_dn_direct_surf_clear_band(:,jcol) 403 end do 404 else 405 do jcol = istartcol,iendcol 406 call indexed_sum(this%sw_dn_direct_surf_clear_g(:,jcol), & 407 & config%i_band_from_reordered_g_sw, & 408 & this%sw_dn_direct_surf_clear_band(:,jcol)) 409 call indexed_sum(this%sw_dn_diffuse_surf_clear_g(:,jcol), & 410 & config%i_band_from_reordered_g_sw, & 411 & this%sw_dn_surf_clear_band(:,jcol)) 412 this%sw_dn_surf_clear_band(:,jcol) & 413 & = this%sw_dn_surf_clear_band(:,jcol) & 414 & + this%sw_dn_direct_surf_clear_band(:,jcol) 415 end do 416 end if 375 417 end if 376 418 … … 383 425 this%sw_dn_direct_surf_canopy (:,istartcol:iendcol) = this%sw_dn_direct_surf_g (:,istartcol:iendcol) 384 426 else if (config%do_nearest_spectral_sw_albedo) then 385 do jcol = istartcol,iendcol 386 call indexed_sum(this%sw_dn_direct_surf_g(:,jcol), & 387 & config%i_albedo_from_band_sw(config%i_band_from_reordered_g_sw), & 388 & this%sw_dn_direct_surf_canopy(:,jcol)) 389 call indexed_sum(this%sw_dn_diffuse_surf_g(:,jcol), & 390 & config%i_albedo_from_band_sw(config%i_band_from_reordered_g_sw), & 391 & this%sw_dn_diffuse_surf_canopy(:,jcol)) 392 end do 427 if (use_indexed_sum_vec) then 428 call indexed_sum_vec(this%sw_dn_direct_surf_g, & 429 & config%i_albedo_from_band_sw(config%i_band_from_reordered_g_sw), & 430 & this%sw_dn_direct_surf_canopy, istartcol, iendcol) 431 call indexed_sum_vec(this%sw_dn_diffuse_surf_g, & 432 & config%i_albedo_from_band_sw(config%i_band_from_reordered_g_sw), & 433 & this%sw_dn_diffuse_surf_canopy, istartcol, iendcol) 434 else 435 do jcol = istartcol,iendcol 436 call indexed_sum(this%sw_dn_direct_surf_g(:,jcol), & 437 & config%i_albedo_from_band_sw(config%i_band_from_reordered_g_sw), & 438 & this%sw_dn_direct_surf_canopy(:,jcol)) 439 call indexed_sum(this%sw_dn_diffuse_surf_g(:,jcol), & 440 & config%i_albedo_from_band_sw(config%i_band_from_reordered_g_sw), & 441 & this%sw_dn_diffuse_surf_canopy(:,jcol)) 442 end do 443 end if 393 444 else 394 445 ! More accurate calculations using weights, but requires … … 425 476 this%lw_dn_surf_canopy(:,istartcol:iendcol) = this%lw_dn_surf_g(:,istartcol:iendcol) 426 477 else if (config%do_nearest_spectral_lw_emiss) then 427 do jcol = istartcol,iendcol 428 call indexed_sum(this%lw_dn_surf_g(:,jcol), & 429 & config%i_emiss_from_band_lw(config%i_band_from_reordered_g_lw), & 430 & this%lw_dn_surf_canopy(:,jcol)) 431 end do 478 if (use_indexed_sum_vec) then 479 call indexed_sum_vec(this%lw_dn_surf_g, & 480 & config%i_emiss_from_band_lw(config%i_band_from_reordered_g_lw), & 481 & this%lw_dn_surf_canopy, istartcol, iendcol) 482 else 483 do jcol = istartcol,iendcol 484 call indexed_sum(this%lw_dn_surf_g(:,jcol), & 485 & config%i_emiss_from_band_lw(config%i_band_from_reordered_g_lw), & 486 & this%lw_dn_surf_canopy(:,jcol)) 487 end do 488 end if 432 489 else 433 490 ! Compute fluxes in each longwave emissivity interval using 434 491 ! weights; first sum over g points to get the values in bands 435 do jcol = istartcol,iendcol 436 call indexed_sum(this%lw_dn_surf_g(:,jcol), & 437 & config%i_band_from_reordered_g_lw, & 438 & lw_dn_surf_band(:,jcol)) 439 end do 492 if (use_indexed_sum_vec) then 493 call indexed_sum_vec(this%lw_dn_surf_g, & 494 & config%i_band_from_reordered_g_lw, & 495 & lw_dn_surf_band, istartcol, iendcol) 496 else 497 do jcol = istartcol,iendcol 498 call indexed_sum(this%lw_dn_surf_g(:,jcol), & 499 & config%i_band_from_reordered_g_lw, & 500 & lw_dn_surf_band(:,jcol)) 501 end do 502 end if 440 503 nalbedoband = size(config%lw_emiss_weights,1) 441 504 this%lw_dn_surf_canopy(:,istartcol:iendcol) = 0.0_jprb … … 465 528 466 529 use yomhook, only : lhook, dr_hook 467 use radiation_c onfig,only : out_of_bounds_2d530 use radiation_check, only : out_of_bounds_2d 468 531 469 532 class(flux_type), intent(inout) :: this … … 497 560 function heating_rate_out_of_physical_bounds(this, nlev, istartcol, iendcol, pressure_hl) result(is_bad) 498 561 499 use radiation_c onfig, only : out_of_bounds_2d562 use radiation_check, only : out_of_bounds_2d 500 563 use radiation_constants, only : AccelDueToGravity 501 564 … … 581 644 end subroutine indexed_sum 582 645 646 !--------------------------------------------------------------------- 647 ! Vectorized version of "add_indexed_sum" 648 subroutine indexed_sum_vec(source, ind, dest, ist, iend) 649 650 real(jprb), intent(in) :: source(:,:) 651 integer, intent(in) :: ind(:) 652 real(jprb), intent(out) :: dest(:,:) 653 integer, intent(in) :: ist, iend 654 655 integer :: ig, jg, jc 656 657 dest = 0.0 658 659 do jg = lbound(source,1), ubound(source,1) 660 ig = ind(jg) 661 do jc = ist, iend 662 dest(ig,jc) = dest(ig,jc) + source(jg,jc) 663 end do 664 end do 665 666 end subroutine indexed_sum_vec 583 667 584 668 !---------------------------------------------------------------------
Note: See TracChangeset
for help on using the changeset viewer.