Home


About us
PALM
PALM-4U
Installation
Release notes
Documentation
Getting Help
Gallery
Developers area

Changeset 3359


Ignore:
Timestamp:
Oct 16, 2018 7:36:26 PM (6 years ago)
Author:
knoop
Message:

Restructured loops and ifs in production_e to ease vectorization on GPUs

File:
1 edited

Legend:

Unmodified
Added
Removed

  • palm/trunk/SOURCE/turbulence_closure_mod.f90

    r3300 r3359  
    2525! -----------------
    2626! $Id$
     27! Restructured loops and ifs in production_e to ease vectorization on GPUs
     28!
     29! 3300 2018-10-02 14:16:54Z gronemeier
    2730! - removed global array wall_flags_0_global, hence reduced accuracy of l_wall
    2831!   calculation
    … …  
    25702573    INTEGER(iwp) ::  surf_e  !< end index of surface elements at given i-j position
    25712574    INTEGER(iwp) ::  surf_s  !< start index of surface elements at given i-j position
     2575    INTEGER(iwp) ::  flag_nr !< number of masking flag
    25722576
    25732577    REAL(wp)     ::  def         !<
    … …  
    25842588    REAL(wp)     ::  wsvs        !< momentum flux w"v"
    25852589
    2586     REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) ::  dudx  !< Gradient of u-component in x-direction
    2587     REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) ::  dudy  !< Gradient of u-component in y-direction
    2588     REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) ::  dudz  !< Gradient of u-component in z-direction
    2589     REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) ::  dvdx  !< Gradient of v-component in x-direction
    2590     REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) ::  dvdy  !< Gradient of v-component in y-direction
    2591     REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) ::  dvdz  !< Gradient of v-component in z-direction
    2592     REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) ::  dwdx  !< Gradient of w-component in x-direction
    2593     REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) ::  dwdy  !< Gradient of w-component in y-direction
    2594     REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) ::  dwdz  !< Gradient of w-component in z-direction
    2595 
     2590    REAL(wp), DIMENSION(nzb+1:nzt) ::  dudx  !< Gradient of u-component in x-direction
     2591    REAL(wp), DIMENSION(nzb+1:nzt) ::  dudy  !< Gradient of u-component in y-direction
     2592    REAL(wp), DIMENSION(nzb+1:nzt) ::  dudz  !< Gradient of u-component in z-direction
     2593    REAL(wp), DIMENSION(nzb+1:nzt) ::  dvdx  !< Gradient of v-component in x-direction
     2594    REAL(wp), DIMENSION(nzb+1:nzt) ::  dvdy  !< Gradient of v-component in y-direction
     2595    REAL(wp), DIMENSION(nzb+1:nzt) ::  dvdz  !< Gradient of v-component in z-direction
     2596    REAL(wp), DIMENSION(nzb+1:nzt) ::  dwdx  !< Gradient of w-component in x-direction
     2597    REAL(wp), DIMENSION(nzb+1:nzt) ::  dwdy  !< Gradient of w-component in y-direction
     2598    REAL(wp), DIMENSION(nzb+1:nzt) ::  dwdz  !< Gradient of w-component in z-direction
     2599
     2600
     2601
     2602!
     2603!-- Calculate TKE production by shear. Calculate gradients at all grid
     2604!-- points first, gradients at surface-bounded grid points will be
     2605!-- overwritten further below.
    25962606    DO  i = nxl, nxr
    2597 
    2598        IF ( constant_flux_layer )  THEN
    2599 
    2600 !
    2601 !--       Calculate TKE production by shear. Calculate gradients at all grid
    2602 !--       points first, gradients at surface-bounded grid points will be
    2603 !--       overwritten further below.
    2604           DO  j = nys, nyn
    2605              DO  k = nzb+1, nzt
    2606 
    2607                 dudx(k,j) =           ( u(k,j,i+1) - u(k,j,i)     ) * ddx
    2608                 dudy(k,j) = 0.25_wp * ( u(k,j+1,i) + u(k,j+1,i+1) -            &
    2609                                         u(k,j-1,i) - u(k,j-1,i+1) ) * ddy
    2610                 dudz(k,j) = 0.5_wp  * ( u(k+1,j,i) + u(k+1,j,i+1) -            &
    2611                                         u(k-1,j,i) - u(k-1,j,i+1) ) *          &
    2612                                                          dd2zu(k)
    2613    
    2614                 dvdx(k,j) = 0.25_wp * ( v(k,j,i+1) + v(k,j+1,i+1) -            &
    2615                                         v(k,j,i-1) - v(k,j+1,i-1) ) * ddx
    2616                 dvdy(k,j) =           ( v(k,j+1,i) - v(k,j,i)     ) * ddy
    2617                 dvdz(k,j) = 0.5_wp  * ( v(k+1,j,i) + v(k+1,j+1,i) -            &
    2618                                         v(k-1,j,i) - v(k-1,j+1,i) ) *          &
    2619                                                          dd2zu(k)
    2620 
    2621                 dwdx(k,j) = 0.25_wp * ( w(k,j,i+1) + w(k-1,j,i+1) -            &
    2622                                         w(k,j,i-1) - w(k-1,j,i-1) ) * ddx
    2623                 dwdy(k,j) = 0.25_wp * ( w(k,j+1,i) + w(k-1,j+1,i) -            &
    2624                                         w(k,j-1,i) - w(k-1,j-1,i) ) * ddy
    2625                 dwdz(k,j) =           ( w(k,j,i)   - w(k-1,j,i)   ) * ddzw(k)
    2626 
    2627              ENDDO
     2607       DO  j = nys, nyn
     2608          DO  k = nzb+1, nzt
     2609
     2610             dudx(k) =           ( u(k,j,i+1) - u(k,j,i)     ) * ddx
     2611             dudy(k) = 0.25_wp * ( u(k,j+1,i) + u(k,j+1,i+1) -                 &
     2612                                   u(k,j-1,i) - u(k,j-1,i+1) ) * ddy
     2613             dudz(k) = 0.5_wp  * ( u(k+1,j,i) + u(k+1,j,i+1) -                 &
     2614                                   u(k-1,j,i) - u(k-1,j,i+1) ) * dd2zu(k)
     2615
     2616             dvdx(k) = 0.25_wp * ( v(k,j,i+1) + v(k,j+1,i+1) -                 &
     2617                                   v(k,j,i-1) - v(k,j+1,i-1) ) * ddx
     2618             dvdy(k) =           ( v(k,j+1,i) - v(k,j,i)     ) * ddy
     2619             dvdz(k) = 0.5_wp  * ( v(k+1,j,i) + v(k+1,j+1,i) -                 &
     2620                                     v(k-1,j,i) - v(k-1,j+1,i) ) * dd2zu(k)
     2621
     2622             dwdx(k) = 0.25_wp * ( w(k,j,i+1) + w(k-1,j,i+1) -                 &
     2623                                   w(k,j,i-1) - w(k-1,j,i-1) ) * ddx
     2624             dwdy(k) = 0.25_wp * ( w(k,j+1,i) + w(k-1,j+1,i) -                 &
     2625                                   w(k,j-1,i) - w(k-1,j-1,i) ) * ddy
     2626             dwdz(k) =           ( w(k,j,i)   - w(k-1,j,i)   ) * ddzw(k)
     2627
    26282628          ENDDO
    26292629
    2630 !
    2631 !--       Position beneath wall
    2632 !--       (2) - Will allways be executed.
    2633 !--       'bottom and wall: use u_0,v_0 and wall functions'
    2634           DO  j = nys, nyn
     2630
     2631          flag_nr = 29
     2632
     2633
     2634          IF ( constant_flux_layer )  THEN
     2635!
     2636
     2637             flag_nr = 0
     2638
     2639!--          Position beneath wall
     2640!--          (2) - Will allways be executed.
     2641!--          'bottom and wall: use u_0,v_0 and wall functions'
    26352642!
    26362643!--          Compute gradients at north- and south-facing surfaces.
    2637 !--          First, for default surfaces, then for urban surfaces. 
     2644!--          First, for default surfaces, then for urban surfaces.
    26382645!--          Note, so far no natural vertical surfaces implemented
    26392646             DO  l = 0, 1
    … …  
    26442651                   usvs        = surf_def_v(l)%mom_flux_tke(0,m)
    26452652                   wsvs        = surf_def_v(l)%mom_flux_tke(1,m)
    2646    
     2653
    26472654                   km_neutral = kappa * ( usvs**2 + wsvs**2 )**0.25_wp         &
    26482655                                   * 0.5_wp * dy
    … …  
    26502657!--                -1.0 for right-facing wall, 1.0 for left-facing wall
    26512658                   sign_dir = MERGE( 1.0_wp, -1.0_wp,                          &
    2652                                      BTEST( wall_flags_0(k,j-1,i), 0 ) )
    2653                    dudy(k,j) = sign_dir * usvs / ( km_neutral + 1E-10_wp )
    2654                    dwdy(k,j) = sign_dir * wsvs / ( km_neutral + 1E-10_wp )
     2659                                     BTEST( wall_flags_0(k,j-1,i), flag_nr ) )
     2660                   dudy(k) = sign_dir * usvs / ( km_neutral + 1E-10_wp )
     2661                   dwdy(k) = sign_dir * wsvs / ( km_neutral + 1E-10_wp )
    26552662                ENDDO
    26562663!
    … …  
    26622669                   usvs        = surf_lsm_v(l)%mom_flux_tke(0,m)
    26632670                   wsvs        = surf_lsm_v(l)%mom_flux_tke(1,m)
    2664    
     2671
    26652672                   km_neutral = kappa * ( usvs**2 + wsvs**2 )**0.25_wp         &
    26662673                                   * 0.5_wp * dy
    … …  
    26682675!--                -1.0 for right-facing wall, 1.0 for left-facing wall
    26692676                   sign_dir = MERGE( 1.0_wp, -1.0_wp,                          &
    2670                                      BTEST( wall_flags_0(k,j-1,i), 0 ) )
    2671                    dudy(k,j) = sign_dir * usvs / ( km_neutral + 1E-10_wp )
    2672                    dwdy(k,j) = sign_dir * wsvs / ( km_neutral + 1E-10_wp )
    2673                 ENDDO 
     2677                                     BTEST( wall_flags_0(k,j-1,i), flag_nr ) )
     2678                   dudy(k) = sign_dir * usvs / ( km_neutral + 1E-10_wp )
     2679                   dwdy(k) = sign_dir * wsvs / ( km_neutral + 1E-10_wp )
     2680                ENDDO
    26742681!
    26752682!--             Urban surfaces
    … …  
    26802687                   usvs        = surf_usm_v(l)%mom_flux_tke(0,m)
    26812688                   wsvs        = surf_usm_v(l)%mom_flux_tke(1,m)
    2682    
     2689
    26832690                   km_neutral = kappa * ( usvs**2 + wsvs**2 )**0.25_wp         &
    26842691                                   * 0.5_wp * dy
    … …  
    26862693!--                -1.0 for right-facing wall, 1.0 for left-facing wall
    26872694                   sign_dir = MERGE( 1.0_wp, -1.0_wp,                          &
    2688                                      BTEST( wall_flags_0(k,j-1,i), 0 ) )
    2689                    dudy(k,j) = sign_dir * usvs / ( km_neutral + 1E-10_wp )
    2690                    dwdy(k,j) = sign_dir * wsvs / ( km_neutral + 1E-10_wp )
    2691                 ENDDO 
     2695                                     BTEST( wall_flags_0(k,j-1,i), flag_nr ) )
     2696                   dudy(k) = sign_dir * usvs / ( km_neutral + 1E-10_wp )
     2697                   dwdy(k) = sign_dir * wsvs / ( km_neutral + 1E-10_wp )
     2698                ENDDO
    26922699             ENDDO
    26932700!
    … …  
    27062713!--                -1.0 for right-facing wall, 1.0 for left-facing wall
    27072714                   sign_dir = MERGE( 1.0_wp, -1.0_wp,                          &
    2708                                      BTEST( wall_flags_0(k,j,i-1), 0 ) )
    2709                    dvdx(k,j) = sign_dir * vsus / ( km_neutral + 1E-10_wp )
    2710                    dwdx(k,j) = sign_dir * wsus / ( km_neutral + 1E-10_wp )
    2711                 ENDDO 
    2712 !
    2713 !--             Natural surfaces                   
     2715                                     BTEST( wall_flags_0(k,j,i-1), flag_nr ) )
     2716                   dvdx(k) = sign_dir * vsus / ( km_neutral + 1E-10_wp )
     2717                   dwdx(k) = sign_dir * wsus / ( km_neutral + 1E-10_wp )
     2718                ENDDO
     2719!
     2720!--             Natural surfaces
    27142721                surf_s = surf_lsm_v(l)%start_index(j,i)
    27152722                surf_e = surf_lsm_v(l)%end_index(j,i)
    … …  
    27242731!--                -1.0 for right-facing wall, 1.0 for left-facing wall
    27252732                   sign_dir = MERGE( 1.0_wp, -1.0_wp,                          &
    2726                                      BTEST( wall_flags_0(k,j,i-1), 0 ) )
    2727                    dvdx(k,j) = sign_dir * vsus / ( km_neutral + 1E-10_wp )
    2728                    dwdx(k,j) = sign_dir * wsus / ( km_neutral + 1E-10_wp )
    2729                 ENDDO   
    2730 !
    2731 !--             Urban surfaces                   
     2733                                     BTEST( wall_flags_0(k,j,i-1), flag_nr ) )
     2734                   dvdx(k) = sign_dir * vsus / ( km_neutral + 1E-10_wp )
     2735                   dwdx(k) = sign_dir * wsus / ( km_neutral + 1E-10_wp )
     2736                ENDDO
     2737!
     2738!--             Urban surfaces
    27322739                surf_s = surf_usm_v(l)%start_index(j,i)
    27332740                surf_e = surf_usm_v(l)%end_index(j,i)
    … …  
    27422749!--                -1.0 for right-facing wall, 1.0 for left-facing wall
    27432750                   sign_dir = MERGE( 1.0_wp, -1.0_wp,                          &
    2744                                      BTEST( wall_flags_0(k,j,i-1), 0 ) )
    2745                    dvdx(k,j) = sign_dir * vsus / ( km_neutral + 1E-10_wp )
    2746                    dwdx(k,j) = sign_dir * wsus / ( km_neutral + 1E-10_wp )
    2747                 ENDDO 
     2751                                     BTEST( wall_flags_0(k,j,i-1), flag_nr ) )
     2752                   dvdx(k) = sign_dir * vsus / ( km_neutral + 1E-10_wp )
     2753                   dwdx(k) = sign_dir * wsus / ( km_neutral + 1E-10_wp )
     2754                ENDDO
    27482755             ENDDO
    27492756!
    … …  
    27542761                k = surf_def_h(0)%k(m)
    27552762!
    2756 !--             Please note, actually, an interpolation of u_0 and v_0 
    2757 !--             onto the grid center would be required. However, this 
     2763!--             Please note, actually, an interpolation of u_0 and v_0
     2764!--             onto the grid center would be required. However, this
    27582765!--             would require several data transfers between 2D-grid and
    2759 !--             wall type. The effect of this missing interpolation is 
     2766!--             wall type. The effect of this missing interpolation is
    27602767!--             negligible. (See also production_e_init).
    2761                 dudz(k,j) = ( u(k+1,j,i) - surf_def_h(0)%u_0(m) ) * dd2zu(k)   
    2762                 dvdz(k,j) = ( v(k+1,j,i) - surf_def_h(0)%v_0(m) ) * dd2zu(k)
    2763          
     2768                dudz(k) = ( u(k+1,j,i) - surf_def_h(0)%u_0(m) ) * dd2zu(k)
     2769                dvdz(k) = ( v(k+1,j,i) - surf_def_h(0)%v_0(m) ) * dd2zu(k)
     2770
    27642771             ENDDO
    27652772!
    … …  
    27702777                k = surf_lsm_h%k(m)
    27712778
    2772                 dudz(k,j) = ( u(k+1,j,i) - surf_lsm_h%u_0(m) ) * dd2zu(k)   
    2773                 dvdz(k,j) = ( v(k+1,j,i) - surf_lsm_h%v_0(m) ) * dd2zu(k)
    2774          
     2779                dudz(k) = ( u(k+1,j,i) - surf_lsm_h%u_0(m) ) * dd2zu(k)
     2780                dvdz(k) = ( v(k+1,j,i) - surf_lsm_h%v_0(m) ) * dd2zu(k)
     2781
    27752782             ENDDO
    27762783!
    … …  
    27812788                k = surf_usm_h%k(m)
    27822789
    2783                 dudz(k,j) = ( u(k+1,j,i) - surf_usm_h%u_0(m) ) * dd2zu(k)   
    2784                 dvdz(k,j) = ( v(k+1,j,i) - surf_usm_h%v_0(m) ) * dd2zu(k)
    2785          
     2790                dudz(k) = ( u(k+1,j,i) - surf_usm_h%u_0(m) ) * dd2zu(k)
     2791                dvdz(k) = ( v(k+1,j,i) - surf_usm_h%v_0(m) ) * dd2zu(k)
     2792
    27862793             ENDDO
    27872794!
    2788 !--          Compute gradients at downward-facing walls, only for 
     2795!--          Compute gradients at downward-facing walls, only for
    27892796!--          non-natural default surfaces
    27902797             surf_s = surf_def_h(1)%start_index(j,i)
    … …  
    27932800                k = surf_def_h(1)%k(m)
    27942801
    2795                 dudz(k,j) = ( surf_def_h(1)%u_0(m) - u(k-1,j,i) ) * dd2zu(k)   
    2796                 dvdz(k,j) = ( surf_def_h(1)%v_0(m) - v(k-1,j,i) ) * dd2zu(k)
     2802                dudz(k) = ( surf_def_h(1)%u_0(m) - u(k-1,j,i) ) * dd2zu(k)
     2803                dvdz(k) = ( surf_def_h(1)%v_0(m) - v(k-1,j,i) ) * dd2zu(k)
    27972804
    27982805             ENDDO
     2806
     2807
     2808          ENDIF
     2809
     2810
     2811          DO  k = nzb+1, nzt
     2812
     2813             def = 2.0_wp * ( dudx(k)**2 + dvdy(k)**2 + dwdz(k)**2 ) +         &
     2814                              dudy(k)**2 + dvdx(k)**2 + dwdx(k)**2 +           &
     2815                              dwdy(k)**2 + dudz(k)**2 + dvdz(k)**2 +           &
     2816                   2.0_wp * ( dvdx(k)*dudy(k) + dwdx(k)*dudz(k) +              &
     2817                              dwdy(k)*dvdz(k) )
     2818
     2819             IF ( def < 0.0_wp )  def = 0.0_wp
     2820
     2821             flag = MERGE( 1.0_wp, 0.0_wp, BTEST(wall_flags_0(k,j,i),flag_nr) )
     2822
     2823             tend(k,j,i) = tend(k,j,i) + km(k,j,i) * def * flag
     2824
    27992825          ENDDO
    28002826
    2801           DO  j = nys, nyn
    2802              DO  k = nzb+1, nzt
    2803 
    2804                 def = 2.0_wp * ( dudx(k,j)**2 + dvdy(k,j)**2 + dwdz(k,j)**2 ) + &
    2805                                  dudy(k,j)**2 + dvdx(k,j)**2 + dwdx(k,j)**2 +   &
    2806                                  dwdy(k,j)**2 + dudz(k,j)**2 + dvdz(k,j)**2 +   &
    2807                       2.0_wp * ( dvdx(k,j)*dudy(k,j) + dwdx(k,j)*dudz(k,j)  +   &
    2808                                  dwdy(k,j)*dvdz(k,j) )
    2809 
    2810                 IF ( def < 0.0_wp )  def = 0.0_wp
    2811 
    2812                 flag  = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_0(k,j,i), 0 ) )
    2813 
    2814                 tend(k,j,i) = tend(k,j,i) + km(k,j,i) * def * flag
    2815 
    2816              ENDDO
    2817           ENDDO
    2818 
    2819        ELSE
    2820 
    2821           DO  j = nys, nyn
    2822 !
    2823 !--          Calculate TKE production by shear. Here, no additional
    2824 !--          wall-bounded code is considered.
    2825 !--          Why?
    2826              DO  k = nzb+1, nzt
    2827 
    2828                 dudx(k,j)  =           ( u(k,j,i+1) - u(k,j,i)     ) * ddx
    2829                 dudy(k,j)  = 0.25_wp * ( u(k,j+1,i) + u(k,j+1,i+1) -           &
    2830                                          u(k,j-1,i) - u(k,j-1,i+1) ) * ddy
    2831                 dudz(k,j)  = 0.5_wp  * ( u(k+1,j,i) + u(k+1,j,i+1) -           &
    2832                                          u(k-1,j,i) - u(k-1,j,i+1) ) *         &
    2833                                                                      dd2zu(k)
    2834 
    2835                 dvdx(k,j)  = 0.25_wp * ( v(k,j,i+1) + v(k,j+1,i+1) -           &
    2836                                          v(k,j,i-1) - v(k,j+1,i-1) ) * ddx
    2837                 dvdy(k,j)  =           ( v(k,j+1,i) - v(k,j,i)     ) * ddy
    2838                 dvdz(k,j)  = 0.5_wp  * ( v(k+1,j,i) + v(k+1,j+1,i) -           &
    2839                                          v(k-1,j,i) - v(k-1,j+1,i) ) *         &
    2840                                                                      dd2zu(k)
    2841 
    2842                 dwdx(k,j)  = 0.25_wp * ( w(k,j,i+1) + w(k-1,j,i+1) -           &
    2843                                          w(k,j,i-1) - w(k-1,j,i-1) ) * ddx
    2844                 dwdy(k,j)  = 0.25_wp * ( w(k,j+1,i) + w(k-1,j+1,i) -           &
    2845                                          w(k,j-1,i) - w(k-1,j-1,i) ) * ddy
    2846                 dwdz(k,j)  =           ( w(k,j,i)   - w(k-1,j,i)   ) *         &
    2847                                                                      ddzw(k)
    2848    
    2849                 def = 2.0_wp * (                                               &
    2850                              dudx(k,j)**2 + dvdy(k,j)**2 + dwdz(k,j)**2        &
    2851                                ) +                                             &
    2852                              dudy(k,j)**2 + dvdx(k,j)**2 + dwdx(k,j)**2 +      &
    2853                              dwdy(k,j)**2 + dudz(k,j)**2 + dvdz(k,j)**2 +      &
    2854                       2.0_wp * (                                               &
    2855                              dvdx(k,j)*dudy(k,j) + dwdx(k,j)*dudz(k,j)  +      &
    2856                              dwdy(k,j)*dvdz(k,j)                               &
    2857                                )
    2858 
    2859                 IF ( def < 0.0_wp )  def = 0.0_wp
    2860 
    2861                 flag  = MERGE( 1.0_wp, 0.0_wp,                                 &
    2862                                BTEST( wall_flags_0(k,j,i), 29 ) )
    2863                 tend(k,j,i) = tend(k,j,i) + km(k,j,i) * def * flag
    2864    
    2865              ENDDO
    2866           ENDDO
    2867 
    2868        ENDIF
    2869 
    2870 !
    2871 !--    If required, calculate TKE production by buoyancy
    2872        IF ( .NOT. neutral )  THEN
    2873 
    2874           IF ( .NOT. humidity )  THEN
    2875 
    2876              IF ( ocean_mode )  THEN
    2877 !
    2878 !--             So far in the ocean no special treatment of density flux
    2879 !--             in the bottom and top surface layer
     2827
     2828       ENDDO
     2829    ENDDO
     2830
     2831!
     2832!-- If required, calculate TKE production by buoyancy
     2833    IF ( .NOT. neutral )  THEN
     2834
     2835       IF ( .NOT. humidity )  THEN
     2836
     2837          IF ( ocean_mode )  THEN
     2838!
     2839!--          So far in the ocean no special treatment of density flux
     2840!--          in the bottom and top surface layer
     2841             DO  i = nxl, nxr
    28802842                DO  j = nys, nyn
    28812843                   DO  k = nzb+1, nzt
    … …  
    28912853                                MERGE( 1.0_wp, 0.0_wp,                         &
    28922854                                       BTEST( wall_flags_0(k,j,i), 9 )         &
    2893                                      ) 
     2855                                     )
    28942856                   ENDDO
    28952857!
    … …  
    29212883                                             use_single_reference_value ) *    &
    29222884                                      drho_air_zw(k) *                         &
    2923                                       surf_def_h(2)%shf(m) 
     2885                                      surf_def_h(2)%shf(m)
    29242886                      ENDDO
    29252887                   ENDIF
    29262888
    29272889                ENDDO
    2928 
    2929              ELSE
    2930 
     2890             ENDDO
     2891
     2892          ELSE ! or IF ( .NOT. ocean_mode )  THEN
     2893
     2894             DO  i = nxl, nxr
    29312895                DO  j = nys, nyn
     2896
    29322897                   DO  k = nzb+1, nzt
    29332898!
    … …  
    29452910                                MERGE( 1.0_wp, 0.0_wp,                         &
    29462911                                       BTEST( wall_flags_0(k,j,i), 9 )         &
    2947                                      ) 
     2912                                     )
    29482913                   ENDDO
    29492914
    … …  
    29602925                                        use_single_reference_value )           &
    29612926                                                   * drho_air_zw(k-1)          &
    2962                                                    * surf_def_h(l)%shf(m)       
    2963                          ENDDO     
     2927                                                   * surf_def_h(l)%shf(m)
     2928                         ENDDO
    29642929                      ENDDO
    29652930!
    … …  
    29732938                                        use_single_reference_value )           &
    29742939                                                   * drho_air_zw(k-1)          &
    2975                                                    * surf_lsm_h%shf(m)   
     2940                                                   * surf_lsm_h%shf(m)
    29762941                      ENDDO
    29772942!
    … …  
    29852950                                        use_single_reference_value )           &
    29862951                                                   * drho_air_zw(k-1)          &
    2987                                                    * surf_usm_h%shf(m)   
    2988                       ENDDO                         
     2952                                                   * surf_usm_h%shf(m)
     2953                      ENDDO
    29892954                   ENDIF
    29902955
    … …  
    29982963                                        use_single_reference_value )           &
    29992964                                                   * drho_air_zw(k)            &
    3000                                                    * surf_def_h(2)%shf(m) 
     2965                                                   * surf_def_h(2)%shf(m)
    30012966                      ENDDO
    30022967                   ENDIF
     2968
    30032969                ENDDO
    3004 
    3005              ENDIF
    3006 
    3007           ELSE
    3008 
     2970             ENDDO
     2971
     2972          ENDIF ! from IF ( .NOT. ocean_mode )
     2973
     2974       ELSE ! or IF ( humidity )  THEN
     2975
     2976          DO  i = nxl, nxr
    30092977             DO  j = nys, nyn
    30102978
    … …  
    30593027                      k1 = 1.0_wp + 0.61_wp * q(k,j,i) - ql(k,j,i)
    30603028                      k2 = 0.61_wp * pt(k,j,i)
    3061                       tend(k,j,i) = tend(k,j,i) -                              & 
     3029                      tend(k,j,i) = tend(k,j,i) -                              &
    30623030                                    kh(k,j,i) * g /                            &
    30633031                                 MERGE( vpt_reference, vpt(k,j,i),             &
    … …  
    30773045                ENDDO
    30783046
    3079              ENDDO
    3080 
    3081              IF ( use_surface_fluxes )  THEN
    3082 
    3083                 DO  j = nys, nyn
     3047                IF ( use_surface_fluxes )  THEN
    30843048!
    30853049!--                Treat horizontal default surfaces
    … …  
    31613125                   surf_e = surf_usm_h%end_index(j,i)
    31623126                   DO  m = surf_s, surf_e
    3163                       k = surf_usm_h%k(m)
     3127                      k = surf_lsm_h%k(m)
    31643128
    31653129                      IF ( .NOT. bulk_cloud_model .AND. .NOT. cloud_droplets ) THEN
    … …  
    31933157                   ENDDO
    31943158
    3195                 ENDDO
    3196 
    3197              ENDIF
    3198 
    3199              IF ( use_top_fluxes )  THEN
    3200 
    3201                 DO  j = nys, nyn
     3159                ENDIF ! from IF ( use_surface_fluxes )  THEN
     3160
     3161                IF ( use_top_fluxes )  THEN
    32023162
    32033163                   surf_s = surf_def_h(2)%start_index(j,i)
    … …  
    32373197                   ENDDO
    32383198
    3239                 ENDDO
    3240 
    3241              ENDIF
    3242 
    3243           ENDIF
     3199                ENDIF ! from IF ( use_top_fluxes )  THEN
     3200
     3201             ENDDO
     3202          ENDDO
    32443203
    32453204       ENDIF
    32463205
    3247     ENDDO
     3206    ENDIF
    32483207
    32493208 END SUBROUTINE production_e
Note: See TracChangeset for help on using the changeset viewer.
                                                                                                                                                                                                                                                                                                                                                                               
  | Impressum | ©Leibniz Universität Hannover |