Ignore:
Timestamp:
May 18, 2024, 8:07:34 PM (7 months ago)
Author:
idelkadi
Message:

The addition of explicit loops with the "omp simd reduction" directive to solve the slowness problem linked to the "sum" command (svn4938 revesion) led to non-reproducibility in MPI and mixed MPI-OMP modes in the case of Tripleclouds and Mcica solvers.
We return to the svn4848 versions of the radiation_tripleclouds_*w.F90 and radiation_mcica_*w.F90 routines.

Location:
LMDZ6/trunk/libf/phylmd/ecrad/radiation
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • LMDZ6/trunk/libf/phylmd/ecrad/radiation/radiation_mcica_lw.F90

    r4853 r4946  
    1818!   2017-07-12  R. Hogan  Call fast adding method if only clouds scatter
    1919!   2017-10-23  R. Hogan  Renamed single-character variables
    20 
    21 #include "ecrad_config.h"
    2220
    2321module radiation_mcica_lw
     
    126124    ! Identify clear-sky layers
    127125    logical :: is_clear_sky_layer(nlev)
    128 
    129     ! Temporary storage for more efficient summation
    130 #ifdef DWD_REDUCTION_OPTIMIZATIONS
    131     real(jprb), dimension(nlev+1,2) :: sum_aux
    132 #else
    133     real(jprb) :: sum_up, sum_dn
    134 #endif
    135126
    136127    ! Index of the highest cloudy layer
     
    188179
    189180      ! Sum over g-points to compute broadband fluxes
    190 #ifdef DWD_REDUCTION_OPTIMIZATIONS
    191       sum_aux(:,:) = 0.0_jprb
    192       do jg = 1,ng
    193         do jlev = 1,nlev+1
    194           sum_aux(jlev,1) = sum_aux(jlev,1) + flux_up_clear(jg,jlev)
    195           sum_aux(jlev,2) = sum_aux(jlev,2) + flux_dn_clear(jg,jlev)
    196         end do
    197       end do
    198       flux%lw_up_clear(jcol,:) = sum_aux(:,1)
    199       flux%lw_dn_clear(jcol,:) = sum_aux(:,2)
    200 #else
    201       do jlev = 1,nlev+1
    202         sum_up = 0.0_jprb
    203         sum_dn = 0.0_jprb
    204         !$omp simd reduction(+:sum_up, sum_dn)
    205         do jg = 1,ng
    206           sum_up = sum_up + flux_up_clear(jg,jlev)
    207           sum_dn = sum_dn + flux_dn_clear(jg,jlev)
    208         end do
    209         flux%lw_up_clear(jcol,jlev) = sum_up
    210         flux%lw_dn_clear(jcol,jlev) = sum_dn
    211       end do
    212 #endif
    213 
     181      flux%lw_up_clear(jcol,:) = sum(flux_up_clear,1)
     182      flux%lw_dn_clear(jcol,:) = sum(flux_dn_clear,1)
    214183      ! Store surface spectral downwelling fluxes
    215184      flux%lw_dn_surf_clear_g(:,jcol) = flux_dn_clear(:,nlev+1)
     
    310279          else
    311280            ! Clear-sky layer: copy over clear-sky values
    312             do jg = 1,ng
    313               reflectance(jg,jlev) = ref_clear(jg,jlev)
    314               transmittance(jg,jlev) = trans_clear(jg,jlev)
    315               source_up(jg,jlev) = source_up_clear(jg,jlev)
    316               source_dn(jg,jlev) = source_dn_clear(jg,jlev)
    317             end do
     281            reflectance(:,jlev) = ref_clear(:,jlev)
     282            transmittance(:,jlev) = trans_clear(:,jlev)
     283            source_up(:,jlev) = source_up_clear(:,jlev)
     284            source_dn(:,jlev) = source_dn_clear(:,jlev)
    318285          end if
    319286        end do
     
    340307       
    341308        ! Store overcast broadband fluxes
    342 #ifdef DWD_REDUCTION_OPTIMIZATIONS
    343         sum_aux(:,:) = 0._jprb
    344         do jg = 1, ng
    345           do jlev = 1, nlev+1
    346             sum_aux(jlev,1) = sum_aux(jlev,1) + flux_up(jg,jlev)
    347             sum_aux(jlev,2) = sum_aux(jlev,2) + flux_dn(jg,jlev)
    348           end do
    349         end do
    350         flux%lw_up(jcol,:) = sum_aux(:,1)
    351         flux%lw_dn(jcol,:) = sum_aux(:,2)
    352 #else
    353         do jlev = 1,nlev+1
    354           sum_up = 0.0_jprb
    355           sum_dn = 0.0_jprb
    356           !$omp simd reduction(+:sum_up, sum_dn)
    357           do jg = 1,ng
    358             sum_up = sum_up + flux_up(jg,jlev)
    359             sum_dn = sum_dn + flux_dn(jg,jlev)
    360           end do
    361           flux%lw_up(jcol,jlev) = sum_up
    362           flux%lw_dn(jcol,jlev) = sum_dn
    363         end do
    364 #endif
     309        flux%lw_up(jcol,:) = sum(flux_up,1)
     310        flux%lw_dn(jcol,:) = sum(flux_dn,1)
    365311
    366312        ! Cloudy flux profiles currently assume completely overcast
    367313        ! skies; perform weighted average with clear-sky profile
    368         do jlev = 1,nlev+1
    369           flux%lw_up(jcol,jlev) =  total_cloud_cover *flux%lw_up(jcol,jlev) &
    370              &       + (1.0_jprb - total_cloud_cover)*flux%lw_up_clear(jcol,jlev)
    371           flux%lw_dn(jcol,jlev) =  total_cloud_cover *flux%lw_dn(jcol,jlev) &
    372              &       + (1.0_jprb - total_cloud_cover)*flux%lw_dn_clear(jcol,jlev)
    373         end do
     314        flux%lw_up(jcol,:) =  total_cloud_cover *flux%lw_up(jcol,:) &
     315             &  + (1.0_jprb - total_cloud_cover)*flux%lw_up_clear(jcol,:)
     316        flux%lw_dn(jcol,:) =  total_cloud_cover *flux%lw_dn(jcol,:) &
     317             &  + (1.0_jprb - total_cloud_cover)*flux%lw_dn_clear(jcol,:)
    374318        ! Store surface spectral downwelling fluxes
    375319        flux%lw_dn_surf_g(:,jcol) = total_cloud_cover*flux_dn(:,nlev+1) &
     
    391335        ! No cloud in profile and clear-sky fluxes already
    392336        ! calculated: copy them over
    393         do jlev = 1,nlev+1
    394           flux%lw_up(jcol,jlev) = flux%lw_up_clear(jcol,jlev)
    395           flux%lw_dn(jcol,jlev) = flux%lw_dn_clear(jcol,jlev)
    396         end do
     337        flux%lw_up(jcol,:) = flux%lw_up_clear(jcol,:)
     338        flux%lw_dn(jcol,:) = flux%lw_dn_clear(jcol,:)
    397339        flux%lw_dn_surf_g(:,jcol) = flux%lw_dn_surf_clear_g(:,jcol)
    398340        if (config%do_lw_derivatives) then
  • LMDZ6/trunk/libf/phylmd/ecrad/radiation/radiation_mcica_sw.F90

    r4853 r4946  
    1717!   2017-04-22  R. Hogan  Store surface fluxes at all g-points
    1818!   2017-10-23  R. Hogan  Renamed single-character variables
    19 
    20 #include "ecrad_config.h"
    2119
    2220module radiation_mcica_sw
     
    121119    ! Total cloud cover output from the cloud generator
    122120    real(jprb) :: total_cloud_cover
    123 
    124     ! Temporary storage for more efficient summation
    125 #ifdef DWD_REDUCTION_OPTIMIZATIONS
    126     real(jprb), dimension(nlev+1,3) :: sum_aux
    127 #else
    128     real(jprb) :: sum_up, sum_dn_diff, sum_dn_dir
    129 #endif
    130121
    131122    ! Number of g points
     
    184175       
    185176        ! Sum over g-points to compute and save clear-sky broadband
    186         ! fluxes. Note that the built-in "sum" function is very slow,
    187         ! and before being replaced by the alternatives below
    188         ! accounted for around 40% of the total cost of this routine.
    189 #ifdef DWD_REDUCTION_OPTIMIZATIONS
    190         ! Optimized summation for the NEC architecture
    191         sum_aux(:,:) = 0.0_jprb
    192         do jg = 1,ng
    193           do jlev = 1,nlev+1
    194             sum_aux(jlev,1) = sum_aux(jlev,1) + flux_up(jg,jlev)
    195             sum_aux(jlev,2) = sum_aux(jlev,2) + flux_dn_direct(jg,jlev)
    196             sum_aux(jlev,3) = sum_aux(jlev,3) + flux_dn_diffuse(jg,jlev)
    197           end do
    198         end do
    199         flux%sw_up_clear(jcol,:) = sum_aux(:,1)
    200         flux%sw_dn_clear(jcol,:) = sum_aux(:,2) + sum_aux(:,3)
     177        ! fluxes
     178        flux%sw_up_clear(jcol,:) = sum(flux_up,1)
    201179        if (allocated(flux%sw_dn_direct_clear)) then
    202           flux%sw_dn_direct_clear(jcol,:) = sum_aux(:,2)
     180          flux%sw_dn_direct_clear(jcol,:) &
     181               &  = sum(flux_dn_direct,1)
     182          flux%sw_dn_clear(jcol,:) = sum(flux_dn_diffuse,1) &
     183               &  + flux%sw_dn_direct_clear(jcol,:)
     184        else
     185          flux%sw_dn_clear(jcol,:) = sum(flux_dn_diffuse,1) &
     186               &  + sum(flux_dn_direct,1)
    203187        end if
    204 #else
    205         ! Optimized summation for the x86-64 architecture
    206         do jlev = 1,nlev+1
    207           sum_up      = 0.0_jprb
    208           sum_dn_diff = 0.0_jprb
    209           sum_dn_dir  = 0.0_jprb
    210           !$omp simd reduction(+:sum_up, sum_dn_diff, sum_dn_dir)
    211           do jg = 1,ng
    212             sum_up      = sum_up      + flux_up(jg,jlev)
    213             sum_dn_diff = sum_dn_diff + flux_dn_diffuse(jg,jlev)
    214             sum_dn_dir  = sum_dn_dir  + flux_dn_direct(jg,jlev)
    215           end do
    216           flux%sw_up_clear(jcol,jlev) = sum_up
    217           flux%sw_dn_clear(jcol,jlev) = sum_dn_diff + sum_dn_dir
    218           if (allocated(flux%sw_dn_direct_clear)) then
    219             flux%sw_dn_direct_clear(jcol,jlev) = sum_dn_dir
    220           end if
    221         end do
    222 #endif
    223        
    224188        ! Store spectral downwelling fluxes at surface
    225         do jg = 1,ng
    226           flux%sw_dn_diffuse_surf_clear_g(jg,jcol) = flux_dn_diffuse(jg,nlev+1)
    227           flux%sw_dn_direct_surf_clear_g(jg,jcol)  = flux_dn_direct(jg,nlev+1)
    228         end do
     189        flux%sw_dn_diffuse_surf_clear_g(:,jcol) = flux_dn_diffuse(:,nlev+1)
     190        flux%sw_dn_direct_surf_clear_g(:,jcol)  = flux_dn_direct(:,nlev+1)
    229191
    230192        ! Do cloudy-sky calculation
     
    287249            else
    288250              ! Clear-sky layer: copy over clear-sky values
    289               do jg = 1,ng
    290                 reflectance(jg,jlev) = ref_clear(jg,jlev)
    291                 transmittance(jg,jlev) = trans_clear(jg,jlev)
    292                 ref_dir(jg,jlev) = ref_dir_clear(jg,jlev)
    293                 trans_dir_diff(jg,jlev) = trans_dir_diff_clear(jg,jlev)
    294                 trans_dir_dir(jg,jlev) = trans_dir_dir_clear(jg,jlev)
    295               end do
     251              reflectance(:,jlev) = ref_clear(:,jlev)
     252              transmittance(:,jlev) = trans_clear(:,jlev)
     253              ref_dir(:,jlev) = ref_dir_clear(:,jlev)
     254              trans_dir_diff(:,jlev) = trans_dir_diff_clear(:,jlev)
     255              trans_dir_dir(:,jlev) = trans_dir_dir_clear(:,jlev)
    296256            end if
    297257          end do
     
    304264         
    305265          ! Store overcast broadband fluxes
    306 #ifdef DWD_REDUCTION_OPTIMIZATIONS
    307           sum_aux(:,:) = 0.0_jprb
    308           do jg = 1,ng
    309             do jlev = 1,nlev+1
    310               sum_aux(jlev,1) = sum_aux(jlev,1) + flux_up(jg,jlev)
    311               sum_aux(jlev,2) = sum_aux(jlev,2) + flux_dn_direct(jg,jlev)
    312               sum_aux(jlev,3) = sum_aux(jlev,3) + flux_dn_diffuse(jg,jlev)
    313             end do
    314           end do
    315           flux%sw_up(jcol,:) = sum_aux(:,1)
    316           flux%sw_dn(jcol,:) = sum_aux(:,2) + sum_aux(:,3)
     266          flux%sw_up(jcol,:) = sum(flux_up,1)
    317267          if (allocated(flux%sw_dn_direct)) then
    318             flux%sw_dn_direct(jcol,:) = sum_aux(:,2)
     268            flux%sw_dn_direct(jcol,:) = sum(flux_dn_direct,1)
     269            flux%sw_dn(jcol,:) = sum(flux_dn_diffuse,1) &
     270                 &  + flux%sw_dn_direct(jcol,:)
     271          else
     272            flux%sw_dn(jcol,:) = sum(flux_dn_diffuse,1) &
     273                 &  + sum(flux_dn_direct,1)
    319274          end if
    320 #else
    321           do jlev = 1,nlev+1
    322             sum_up      = 0.0_jprb
    323             sum_dn_diff = 0.0_jprb
    324             sum_dn_dir  = 0.0_jprb
    325             !$omp simd reduction(+:sum_up, sum_dn_diff, sum_dn_dir)
    326             do jg = 1,ng
    327               sum_up      = sum_up      + flux_up(jg,jlev)
    328               sum_dn_diff = sum_dn_diff + flux_dn_diffuse(jg,jlev)
    329               sum_dn_dir  = sum_dn_dir  + flux_dn_direct(jg,jlev)
    330             end do
    331             flux%sw_up(jcol,jlev) = sum_up
    332             flux%sw_dn(jcol,jlev) = sum_dn_diff + sum_dn_dir
    333             if (allocated(flux%sw_dn_direct)) then
    334               flux%sw_dn_direct(jcol,jlev) = sum_dn_dir
    335             end if
    336           end do
    337 #endif
    338          
     275
    339276          ! Cloudy flux profiles currently assume completely overcast
    340277          ! skies; perform weighted average with clear-sky profile
    341           do jlev = 1, nlev+1
    342             flux%sw_up(jcol,jlev) =  total_cloud_cover *flux%sw_up(jcol,jlev) &
    343                  &     + (1.0_jprb - total_cloud_cover)*flux%sw_up_clear(jcol,jlev)
    344             flux%sw_dn(jcol,jlev) =  total_cloud_cover *flux%sw_dn(jcol,jlev) &
    345                  &     + (1.0_jprb - total_cloud_cover)*flux%sw_dn_clear(jcol,jlev)
    346             if (allocated(flux%sw_dn_direct)) then
    347               flux%sw_dn_direct(jcol,jlev) = total_cloud_cover *flux%sw_dn_direct(jcol,jlev) &
    348                    &  + (1.0_jprb - total_cloud_cover)*flux%sw_dn_direct_clear(jcol,jlev)
    349             end if
    350           end do
     278          flux%sw_up(jcol,:) =  total_cloud_cover *flux%sw_up(jcol,:) &
     279               &  + (1.0_jprb - total_cloud_cover)*flux%sw_up_clear(jcol,:)
     280          flux%sw_dn(jcol,:) =  total_cloud_cover *flux%sw_dn(jcol,:) &
     281               &  + (1.0_jprb - total_cloud_cover)*flux%sw_dn_clear(jcol,:)
     282          if (allocated(flux%sw_dn_direct)) then
     283            flux%sw_dn_direct(jcol,:) = total_cloud_cover *flux%sw_dn_direct(jcol,:) &
     284                 &  + (1.0_jprb - total_cloud_cover)*flux%sw_dn_direct_clear(jcol,:)
     285          end if
    351286          ! Likewise for surface spectral fluxes
    352           do jg = 1,ng
    353             flux%sw_dn_diffuse_surf_g(jg,jcol) = flux_dn_diffuse(jg,nlev+1)
    354             flux%sw_dn_direct_surf_g(jg,jcol)  = flux_dn_direct(jg,nlev+1)
    355             flux%sw_dn_diffuse_surf_g(jg,jcol) = total_cloud_cover *flux%sw_dn_diffuse_surf_g(jg,jcol) &
    356                  &                 + (1.0_jprb - total_cloud_cover)*flux%sw_dn_diffuse_surf_clear_g(jg,jcol)
    357             flux%sw_dn_direct_surf_g(jg,jcol)  = total_cloud_cover *flux%sw_dn_direct_surf_g(jg,jcol) &
    358                  &                 + (1.0_jprb - total_cloud_cover)*flux%sw_dn_direct_surf_clear_g(jg,jcol)
    359           end do
    360 
     287          flux%sw_dn_diffuse_surf_g(:,jcol) = flux_dn_diffuse(:,nlev+1)
     288          flux%sw_dn_direct_surf_g(:,jcol)  = flux_dn_direct(:,nlev+1)
     289          flux%sw_dn_diffuse_surf_g(:,jcol) = total_cloud_cover *flux%sw_dn_diffuse_surf_g(:,jcol) &
     290               &     + (1.0_jprb - total_cloud_cover)*flux%sw_dn_diffuse_surf_clear_g(:,jcol)
     291          flux%sw_dn_direct_surf_g(:,jcol) = total_cloud_cover *flux%sw_dn_direct_surf_g(:,jcol) &
     292               &     + (1.0_jprb - total_cloud_cover)*flux%sw_dn_direct_surf_clear_g(:,jcol)
     293         
    361294        else
    362295          ! No cloud in profile and clear-sky fluxes already
    363296          ! calculated: copy them over
    364           do jlev = 1, nlev+1
    365             flux%sw_up(jcol,jlev) = flux%sw_up_clear(jcol,jlev)
    366             flux%sw_dn(jcol,jlev) = flux%sw_dn_clear(jcol,jlev)
    367             if (allocated(flux%sw_dn_direct)) then
    368               flux%sw_dn_direct(jcol,jlev) = flux%sw_dn_direct_clear(jcol,jlev)
    369             end if
    370           end do
    371           do jg = 1,ng
    372             flux%sw_dn_diffuse_surf_g(jg,jcol) = flux%sw_dn_diffuse_surf_clear_g(jg,jcol)
    373             flux%sw_dn_direct_surf_g(jg,jcol)  = flux%sw_dn_direct_surf_clear_g(jg,jcol)
    374           end do
     297          flux%sw_up(jcol,:) = flux%sw_up_clear(jcol,:)
     298          flux%sw_dn(jcol,:) = flux%sw_dn_clear(jcol,:)
     299          if (allocated(flux%sw_dn_direct)) then
     300            flux%sw_dn_direct(jcol,:) = flux%sw_dn_direct_clear(jcol,:)
     301          end if
     302          flux%sw_dn_diffuse_surf_g(:,jcol) = flux%sw_dn_diffuse_surf_clear_g(:,jcol)
     303          flux%sw_dn_direct_surf_g(:,jcol)  = flux%sw_dn_direct_surf_clear_g(:,jcol)
    375304
    376305        end if ! Cloud is present in profile
     
    378307      else
    379308        ! Set fluxes to zero if sun is below the horizon
    380         do jlev = 1, nlev+1
    381           flux%sw_up(jcol,jlev) = 0.0_jprb
    382           flux%sw_dn(jcol,jlev) = 0.0_jprb
    383           if (allocated(flux%sw_dn_direct)) then
    384             flux%sw_dn_direct(jcol,jlev) = 0.0_jprb
    385           end if
    386           flux%sw_up_clear(jcol,jlev) = 0.0_jprb
    387           flux%sw_dn_clear(jcol,jlev) = 0.0_jprb
    388           if (allocated(flux%sw_dn_direct_clear)) then
    389             flux%sw_dn_direct_clear(jcol,jlev) = 0.0_jprb
    390           end if
    391         end do
    392         do jg = 1,ng
    393           flux%sw_dn_diffuse_surf_g(jg,jcol) = 0.0_jprb
    394           flux%sw_dn_direct_surf_g(jg,jcol)  = 0.0_jprb
    395           flux%sw_dn_diffuse_surf_clear_g(jg,jcol) = 0.0_jprb
    396           flux%sw_dn_direct_surf_clear_g(jg,jcol)  = 0.0_jprb
    397         end do
     309        flux%sw_up(jcol,:) = 0.0_jprb
     310        flux%sw_dn(jcol,:) = 0.0_jprb
     311        if (allocated(flux%sw_dn_direct)) then
     312          flux%sw_dn_direct(jcol,:) = 0.0_jprb
     313        end if
     314        flux%sw_up_clear(jcol,:) = 0.0_jprb
     315        flux%sw_dn_clear(jcol,:) = 0.0_jprb
     316        if (allocated(flux%sw_dn_direct_clear)) then
     317          flux%sw_dn_direct_clear(jcol,:) = 0.0_jprb
     318        end if
     319        flux%sw_dn_diffuse_surf_g(:,jcol) = 0.0_jprb
     320        flux%sw_dn_direct_surf_g(:,jcol)  = 0.0_jprb
     321        flux%sw_dn_diffuse_surf_clear_g(:,jcol) = 0.0_jprb
     322        flux%sw_dn_direct_surf_clear_g(:,jcol)  = 0.0_jprb
    398323      end if ! Sun above horizon
    399324
  • LMDZ6/trunk/libf/phylmd/ecrad/radiation/radiation_tripleclouds_lw.F90

    r4853 r4946  
    170170    logical :: is_clear_sky_layer(0:nlev+1)
    171171
    172     ! Temporaries to speed up summations
    173     real(jprb) :: sum_dn, sum_up
    174    
    175172    ! Index of the highest cloudy layer
    176173    integer :: i_cloud_top
     
    264261      if (config%do_clear) then
    265262        ! Sum over g-points to compute broadband fluxes
    266         do jlev = 1,nlev+1
    267           sum_up = 0.0_jprb
    268           sum_dn = 0.0_jprb
    269           !$omp simd reduction(+:sum_up, sum_dn)
    270           do jg = 1,ng
    271             sum_up = sum_up + flux_up_clear(jg,jlev)
    272             sum_dn = sum_dn + flux_dn_clear(jg,jlev)
    273           end do
    274           flux%lw_up_clear(jcol,jlev) = sum_up
    275           flux%lw_dn_clear(jcol,jlev) = sum_dn
    276         end do
    277 
     263        flux%lw_up_clear(jcol,:) = sum(flux_up_clear,1)
     264        flux%lw_dn_clear(jcol,:) = sum(flux_dn_clear,1)
    278265        ! Store surface spectral downwelling fluxes / TOA upwelling
    279         do jg = 1,ng
    280           flux%lw_dn_surf_clear_g(jg,jcol) = flux_dn_clear(jg,nlev+1)
    281           flux%lw_up_toa_clear_g (jg,jcol) = flux_up_clear(jg,1)
    282         end do
     266        flux%lw_dn_surf_clear_g(:,jcol) = flux_dn_clear(:,nlev+1)
     267        flux%lw_up_toa_clear_g (:,jcol) = flux_up_clear(:,1)
    283268        ! Save the spectral fluxes if required
    284269        if (config%do_save_spectral_flux) then
     
    468453          end if
    469454        else
    470           sum_dn = 0.0_jprb
    471           !$omp simd reduction(+:sum_dn)
    472           do jg = 1,ng
    473             sum_dn = sum_dn + flux_dn_clear(jg,jlev)
    474           end do
    475           flux%lw_dn(jcol,jlev) = sum_dn
     455          flux%lw_dn(jcol,:) = sum(flux_dn_clear(:,jlev))
    476456          if (config%do_save_spectral_flux) then
    477457            call indexed_sum(flux_dn_clear(:,jlev), &
     
    490470           &  + total_albedo(:,1,i_cloud_top)*flux_dn_clear(:,i_cloud_top)
    491471      flux_up(:,2:) = 0.0_jprb
    492 
    493       sum_up = 0.0_jprb
    494       !$omp simd reduction(+:sum_up)
    495       do jg = 1,ng
    496         sum_up = sum_up + flux_up(jg,1)
    497       end do
    498       flux%lw_up(jcol,i_cloud_top) = sum_up
    499 
     472      flux%lw_up(jcol,i_cloud_top) = sum(flux_up(:,1))
    500473      if (config%do_save_spectral_flux) then
    501474        call indexed_sum(flux_up(:,1), &
     
    505478      do jlev = i_cloud_top-1,1,-1
    506479        flux_up(:,1) = trans_clear(:,jlev)*flux_up(:,1) + source_up_clear(:,jlev)
    507         sum_up = 0.0_jprb
    508         !$omp simd reduction(+:sum_up)
    509         do jg = 1,ng
    510           sum_up = sum_up + flux_up(jg,1)
    511         end do
    512         flux%lw_up(jcol,jlev) = sum_up
     480        flux%lw_up(jcol,jlev) = sum(flux_up(:,1))
    513481        if (config%do_save_spectral_flux) then
    514482          call indexed_sum(flux_up(:,1), &
     
    560528
    561529        ! Store the broadband fluxes
    562         sum_up = 0.0_jprb
    563         sum_dn = 0.0_jprb
    564         do jreg = 1,nregions
    565           !$omp simd reduction(+:sum_up, sum_dn)
    566           do jg = 1,ng
    567             sum_up = sum_up + flux_up(jg,jreg)
    568             sum_dn = sum_dn + flux_dn(jg,jreg)
    569           end do
    570         end do
    571         flux%lw_up(jcol,jlev+1) = sum_up
    572         flux%lw_dn(jcol,jlev+1) = sum_dn
     530        flux%lw_up(jcol,jlev+1) = sum(sum(flux_up,1))
     531        flux%lw_dn(jcol,jlev+1) = sum(sum(flux_dn,1))
    573532
    574533        ! Save the spectral fluxes if required
  • LMDZ6/trunk/libf/phylmd/ecrad/radiation/radiation_tripleclouds_lw.F90.or

    r4773 r4946  
    170170    logical :: is_clear_sky_layer(0:nlev+1)
    171171
     172    ! Temporaries to speed up summations
     173    real(jprb) :: sum_dn, sum_up
     174   
    172175    ! Index of the highest cloudy layer
    173176    integer :: i_cloud_top
     
    249252        call calc_ref_trans_lw(ng*nlev, &
    250253             &  od(:,:,jcol), ssa(:,:,jcol), g(:,:,jcol), &
    251              &  planck_hl(:,1:jlev,jcol), planck_hl(:,2:jlev+1,jcol), &
     254             &  planck_hl(:,1:nlev,jcol), planck_hl(:,2:nlev+1,jcol), &
    252255             &  ref_clear, trans_clear, &
    253256             &  source_up_clear, source_dn_clear)
     
    261264      if (config%do_clear) then
    262265        ! Sum over g-points to compute broadband fluxes
    263         flux%lw_up_clear(jcol,:) = sum(flux_up_clear,1)
    264         flux%lw_dn_clear(jcol,:) = sum(flux_dn_clear,1)
     266        do jlev = 1,nlev+1
     267          sum_up = 0.0_jprb
     268          sum_dn = 0.0_jprb
     269          !$omp simd reduction(+:sum_up, sum_dn)
     270          do jg = 1,ng
     271            sum_up = sum_up + flux_up_clear(jg,jlev)
     272            sum_dn = sum_dn + flux_dn_clear(jg,jlev)
     273          end do
     274          flux%lw_up_clear(jcol,jlev) = sum_up
     275          flux%lw_dn_clear(jcol,jlev) = sum_dn
     276        end do
     277
    265278        ! Store surface spectral downwelling fluxes / TOA upwelling
    266         flux%lw_dn_surf_clear_g(:,jcol) = flux_dn_clear(:,nlev+1)
    267         flux%lw_up_toa_clear_g (:,jcol) = flux_up_clear(:,1)
     279        do jg = 1,ng
     280          flux%lw_dn_surf_clear_g(jg,jcol) = flux_dn_clear(jg,nlev+1)
     281          flux%lw_up_toa_clear_g (jg,jcol) = flux_up_clear(jg,1)
     282        end do
    268283        ! Save the spectral fluxes if required
    269284        if (config%do_save_spectral_flux) then
     
    453468          end if
    454469        else
    455           flux%lw_dn(jcol,:) = sum(flux_dn_clear(:,jlev))
     470          sum_dn = 0.0_jprb
     471          !$omp simd reduction(+:sum_dn)
     472          do jg = 1,ng
     473            sum_dn = sum_dn + flux_dn_clear(jg,jlev)
     474          end do
     475          flux%lw_dn(jcol,jlev) = sum_dn
    456476          if (config%do_save_spectral_flux) then
    457477            call indexed_sum(flux_dn_clear(:,jlev), &
     
    470490           &  + total_albedo(:,1,i_cloud_top)*flux_dn_clear(:,i_cloud_top)
    471491      flux_up(:,2:) = 0.0_jprb
    472       flux%lw_up(jcol,i_cloud_top) = sum(flux_up(:,1))
     492
     493      sum_up = 0.0_jprb
     494      !$omp simd reduction(+:sum_up)
     495      do jg = 1,ng
     496        sum_up = sum_up + flux_up(jg,1)
     497      end do
     498      flux%lw_up(jcol,i_cloud_top) = sum_up
     499
    473500      if (config%do_save_spectral_flux) then
    474501        call indexed_sum(flux_up(:,1), &
     
    478505      do jlev = i_cloud_top-1,1,-1
    479506        flux_up(:,1) = trans_clear(:,jlev)*flux_up(:,1) + source_up_clear(:,jlev)
    480         flux%lw_up(jcol,jlev) = sum(flux_up(:,1))
     507        sum_up = 0.0_jprb
     508        !$omp simd reduction(+:sum_up)
     509        do jg = 1,ng
     510          sum_up = sum_up + flux_up(jg,1)
     511        end do
     512        flux%lw_up(jcol,jlev) = sum_up
    481513        if (config%do_save_spectral_flux) then
    482514          call indexed_sum(flux_up(:,1), &
     
    528560
    529561        ! Store the broadband fluxes
    530         flux%lw_up(jcol,jlev+1) = sum(sum(flux_up,1))
    531         flux%lw_dn(jcol,jlev+1) = sum(sum(flux_dn,1))
     562        sum_up = 0.0_jprb
     563        sum_dn = 0.0_jprb
     564        do jreg = 1,nregions
     565          !$omp simd reduction(+:sum_up, sum_dn)
     566          do jg = 1,ng
     567            sum_up = sum_up + flux_up(jg,jreg)
     568            sum_dn = sum_dn + flux_dn(jg,jreg)
     569          end do
     570        end do
     571        flux%lw_up(jcol,jlev+1) = sum_up
     572        flux%lw_dn(jcol,jlev+1) = sum_dn
    532573
    533574        ! Save the spectral fluxes if required
  • LMDZ6/trunk/libf/phylmd/ecrad/radiation/radiation_tripleclouds_sw.F90

    r4853 r4946  
    7474    ! Gas and aerosol optical depth, single-scattering albedo and
    7575    ! asymmetry factor at each shortwave g-point
    76     real(jprb), intent(in), dimension(config%n_g_sw,nlev,istartcol:iendcol) &
    77          &  :: od, ssa, g
     76!    real(jprb), intent(in), dimension(istartcol:iendcol,nlev,config%n_g_sw) :: &
     77    real(jprb), intent(in), dimension(config%n_g_sw,nlev,istartcol:iendcol) :: &
     78         &  od, ssa, g
    7879
    7980    ! Cloud and precipitation optical depth, single-scattering albedo and
    8081    ! asymmetry factor in each shortwave band
    81     real(jprb), intent(in), dimension(config%n_bands_sw,nlev,istartcol:iendcol) &
    82          &  :: od_cloud, ssa_cloud, g_cloud
     82    real(jprb), intent(in), dimension(config%n_bands_sw,nlev,istartcol:iendcol) :: &
     83         &  od_cloud, ssa_cloud, g_cloud
    8384
    8485    ! Optical depth, single scattering albedo and asymmetry factor in
     
    9192    ! flux into a plane perpendicular to the incoming radiation at
    9293    ! top-of-atmosphere in each of the shortwave g points
    93     real(jprb), intent(in), dimension(config%n_g_sw,istartcol:iendcol) &
    94          &  :: albedo_direct, albedo_diffuse, incoming_sw
     94    real(jprb), intent(in), dimension(config%n_g_sw,istartcol:iendcol) :: &
     95         &  albedo_direct, albedo_diffuse, incoming_sw
    9596
    9697    ! Output
     
    165166    real(jprb) :: scat_od, scat_od_cloud
    166167
    167     ! Temporaries to speed up summations
    168     real(jprb) :: sum_dn_diff, sum_dn_dir, sum_up
    169 
    170     ! Local cosine of solar zenith angle
    171168    real(jprb) :: mu0
    172169
     
    447444      end if
    448445     
    449       ! Store the TOA broadband fluxes, noting that there is no
    450       ! diffuse downwelling at TOA. The intrinsic "sum" command has
    451       ! been found to be very slow; better performance is found on
    452       ! x86-64 architecture with explicit loops and the "omp simd
    453       ! reduction" directive.
    454       sum_up     = 0.0_jprb
    455       sum_dn_dir = 0.0_jprb
    456       do jreg = 1,nregions
    457         !$omp simd reduction(+:sum_up, sum_dn_dir)
    458         do jg = 1,ng
    459           sum_up     = sum_up     + flux_up(jg,jreg)
    460           sum_dn_dir = sum_dn_dir + direct_dn(jg,jreg)
    461         end do
    462       end do
    463       flux%sw_up(jcol,1) = sum_up
    464       flux%sw_dn(jcol,1) = mu0 * sum_dn_dir
     446      ! Store the TOA broadband fluxes
     447      flux%sw_up(jcol,1) = sum(sum(flux_up,1))
     448      flux%sw_dn(jcol,1) = mu0 * sum(sum(direct_dn,1))
    465449      if (allocated(flux%sw_dn_direct)) then
    466450        flux%sw_dn_direct(jcol,1) = flux%sw_dn(jcol,1)
    467451      end if
    468452      if (config%do_clear) then
    469         sum_up     = 0.0_jprb
    470         sum_dn_dir = 0.0_jprb
    471         !$omp simd reduction(+:sum_up, sum_dn_dir)
    472         do jg = 1,ng
    473           sum_up     = sum_up     + flux_up_clear(jg)
    474           sum_dn_dir = sum_dn_dir + direct_dn_clear(jg)
    475         end do
    476         flux%sw_up_clear(jcol,1) = sum_up
    477         flux%sw_dn_clear(jcol,1) = mu0 * sum_dn_dir
     453        flux%sw_up_clear(jcol,1) = sum(flux_up_clear)
     454        flux%sw_dn_clear(jcol,1) = mu0 * sum(direct_dn_clear)
    478455        if (allocated(flux%sw_dn_direct_clear)) then
    479456          flux%sw_dn_direct_clear(jcol,1) = flux%sw_dn_clear(jcol,1)
     
    490467             &           config%i_spec_from_reordered_g_sw, &
    491468             &           flux%sw_dn_band(:,jcol,1))
    492         flux%sw_dn_band(:,jcol,1) = mu0 * flux%sw_dn_band(:,jcol,1)
     469        flux%sw_dn_band(:,jcol,1) = &
     470             &  mu0 * flux%sw_dn_band(:,jcol,1)
    493471        if (allocated(flux%sw_dn_direct_band)) then
    494472          flux%sw_dn_direct_band(:,jcol,1) = flux%sw_dn_band(:,jcol,1)
     
    571549               ! nothing to do
    572550
    573         ! Store the broadband fluxes. The intrinsic "sum" command has
    574         ! been found to be very slow; better performance is found on
    575         ! x86-64 architecture with explicit loops and the "omp simd
    576         ! reduction" directive.
    577         sum_up      = 0.0_jprb
    578         sum_dn_dir  = 0.0_jprb
    579         sum_dn_diff = 0.0_jprb
    580         do jreg = 1,nregions
    581           !$omp simd reduction(+:sum_up, sum_dn_diff, sum_dn_dir)
    582           do jg = 1,ng
    583             sum_up      = sum_up      + flux_up(jg,jreg)
    584             sum_dn_diff = sum_dn_diff + flux_dn(jg,jreg)
    585             sum_dn_dir  = sum_dn_dir  + direct_dn(jg,jreg)
    586           end do
    587         end do
    588         flux%sw_up(jcol,jlev+1) = sum_up
    589         flux%sw_dn(jcol,jlev+1) = mu0 * sum_dn_dir + sum_dn_diff
     551        ! Store the broadband fluxes
     552        flux%sw_up(jcol,jlev+1) = sum(sum(flux_up,1))
    590553        if (allocated(flux%sw_dn_direct)) then
    591           flux%sw_dn_direct(jcol,jlev+1) = mu0 * sum_dn_dir
     554          flux%sw_dn_direct(jcol,jlev+1) = mu0 * sum(sum(direct_dn,1))
     555          flux%sw_dn(jcol,jlev+1) &
     556               &  = flux%sw_dn_direct(jcol,jlev+1) + sum(sum(flux_dn,1))
     557        else
     558          flux%sw_dn(jcol,jlev+1) = mu0 * sum(sum(direct_dn,1)) + sum(sum(flux_dn,1))   
    592559        end if
    593560        if (config%do_clear) then
    594           sum_up      = 0.0_jprb
    595           sum_dn_dir  = 0.0_jprb
    596           sum_dn_diff = 0.0_jprb
    597           !$omp simd reduction(+:sum_up, sum_dn_diff, sum_dn_dir)
    598           do jg = 1,ng
    599             sum_up      = sum_up      + flux_up_clear(jg)
    600             sum_dn_diff = sum_dn_diff + flux_dn_clear(jg)
    601             sum_dn_dir  = sum_dn_dir  + direct_dn_clear(jg)
    602           end do
    603           flux%sw_up_clear(jcol,jlev+1) = sum_up
    604           flux%sw_dn_clear(jcol,jlev+1) = mu0 * sum_dn_dir + sum_dn_diff
     561          flux%sw_up_clear(jcol,jlev+1) = sum(flux_up_clear)
    605562          if (allocated(flux%sw_dn_direct_clear)) then
    606             flux%sw_dn_direct_clear(jcol,jlev+1) = mu0 * sum_dn_dir
     563            flux%sw_dn_direct_clear(jcol,jlev+1) = mu0 * sum(direct_dn_clear)
     564            flux%sw_dn_clear(jcol,jlev+1) &
     565                 &  = flux%sw_dn_direct_clear(jcol,jlev+1) + sum(flux_dn_clear)
     566          else
     567            flux%sw_dn_clear(jcol,jlev+1) = mu0 * sum(direct_dn_clear) &
     568                 &  + sum(flux_dn_clear)
    607569          end if
    608570        end if
     
    643605          end if
    644606        end if
     607
    645608      end do ! Final loop over levels
    646609     
Note: See TracChangeset for help on using the changeset viewer.