Changeset 4168 for palm/trunk/SOURCE/radiation_model_mod.f90

Timestamp:

Aug 16, 2019 1:50:17 PM (5 years ago)

Author:

suehring

Message:

Replace get_topography_top_index functions by pre-calculated arrays in order to save computational resources

File:

• palm/trunk/SOURCE/radiation_model_mod.f90 (modified) (20 diffs)

palm/trunk/SOURCE/radiation_model_mod.f90

@@ -28 +28 @@
 ! -----------------
 ! $Id$
-! Changed behaviour of masked output over surface to follow terrain and ignore 
-! buildings (J.Resler, T.Gronemeier)
+! Replace function get_topography_top_index by topo_top_ind
 ! 
 ! 4157 2019-08-14 09:19:12Z suehring
@@ -691 +690 @@
     USE indices,                                                               &
         ONLY:  nnx, nny, nx, nxl, nxlg, nxr, nxrg, ny, nyn, nyng, nys, nysg,   &
-               nzb, nzt
+               nzb, nzt, topo_top_ind
 
     USE, INTRINSIC :: iso_c_binding
@@ -737 +736 @@
 
     USE surface_mod,                                                           &
-        ONLY:  get_topography_top_index, get_topography_top_index_ji,          &
-               ind_pav_green, ind_veg_wall, ind_wat_win,                       &
+        ONLY:  ind_pav_green, ind_veg_wall, ind_wat_win,                       &
                surf_lsm_h, surf_lsm_v, surf_type, surf_usm_h, surf_usm_v,      &
                vertical_surfaces_exist
@@ -745 +743 @@
 
     CHARACTER(10) :: radiation_scheme = 'clear-sky' ! 'constant', 'clear-sky', or 'rrtmg'
-
-    REAL(wp), PARAMETER ::  fill_value = -9999.0_wp       !< value for the _FillValue attribute
 
 !
@@ -3583 +3579 @@
 !
 !--       Determine minimum topography top index. 
-          k_topo_l = MINVAL( get_topography_top_index( 's' ) )
+          k_topo_l = MINVAL( topo_top_ind(nys:nyn,nxl:nxr,0) )
 #if defined( __parallel )
           CALL MPI_ALLREDUCE( k_topo_l, k_topo, 1, MPI_INTEGER, MPI_MIN, &
@@ -3773 +3769 @@
              DO  i = nxl, nxr
                 DO  j = nys, nyn
-                   k_topo_l = get_topography_top_index_ji( j, i, 's' )
+                   k_topo_l = topo_top_ind(j,i,0)
                    DO k = k_topo_l+1, nzt+1
                       rad_lw_hr(k,j,i)     = rrtm_lwhr(0,k-k_topo_l)  * d_hours_day
@@ -4066 +4062 @@
 !
 !--             Obtain topography top index (lower bound of RRTMG)
-                k_topo = get_topography_top_index_ji( j, i, 's' )
+                k_topo = topo_top_ind(j,i,0)
 
                 IF ( lw_radiation )  THEN
@@ -5688 +5684 @@
 !
 !--          Following expression equals former kk = k - nzb_s_inner(j,i)
-             kk = k - get_topography_top_index_ji( j, i, 's' )  !- lad arrays are defined flat
+             kk = k - topo_top_ind(j,i,0)  !- lad arrays are defined flat
              pc_heating_rate(kk, j, i) = (pcbinsw(ipcgb) + pcbinlw(ipcgb)) &
                  * pchf_prep(k) * pt(k, j, i) !-- = dT/dt
@@ -5701 +5697 @@
                  j = pcbl(iy, ipcgb)
                  k = pcbl(iz, ipcgb)
-                 kk = k - get_topography_top_index_ji( j, i, 's' )  !- lad arrays are defined flat
+                 kk = k - topo_top_ind(j,i,0)  !- lad arrays are defined flat
                  CALL pcm_calc_transpiration_rate( i, j, k, kk, pcbinsw(ipcgb), pcbinlw(ipcgb), &
                                                    pc_transpiration_rate(kk,j,i), pc_latent_rate(kk,j,i) )
@@ -6217 +6213 @@
 !--    removed later). The following contruct finds the lowest / largest index
 !--    for any upward-facing wall (see bit 12).
-       nzubl = MINVAL( get_topography_top_index( 's' ) )
-       nzutl = MAXVAL( get_topography_top_index( 's' ) )
+       nzubl = MINVAL( topo_top_ind(nys:nyn,nxl:nxr,0) )
+       nzutl = MAXVAL( topo_top_ind(nys:nyn,nxl:nxr,0) )
 
        nzubl = MAX( nzubl, nzb )
@@ -6235 +6231 @@
 !
 !--                Find topography top index
-                   k_topo = get_topography_top_index_ji( j, i, 's' )
+                   k_topo = topo_top_ind(j,i,0)
 
                    DO k = nzt+1, 0, -1
@@ -6358 +6354 @@
 !
 !--                Find topography top index
-                   k_topo = get_topography_top_index_ji( j, i, 's' )
+                   k_topo = topo_top_ind(j,i,0)
 
                    DO k = k_topo + 1, pct(j,i)
@@ -6793 +6789 @@
         ALLOCATE( nzterrl_l((nyn-nys+1)*(nxr-nxl+1)) )
         nzterrl(nys:nyn,nxl:nxr) => nzterrl_l(1:(nyn-nys+1)*(nxr-nxl+1))
-        nzterrl = get_topography_top_index( 's' )
+        nzterrl = topo_top_ind(nys:nyn,nxl:nxr,0)
         CALL MPI_AllGather( nzterrl_l, nnx*nny, MPI_INTEGER, &
                             nzterr, nnx*nny, MPI_INTEGER, comm2d, ierr )
@@ -6803 +6799 @@
         DEALLOCATE(nzterrl_l)
 #else
-        nzterr = RESHAPE( get_topography_top_index( 's' ), (/(nx+1)*(ny+1)/) )
+        nzterr = RESHAPE( topo_top_ind(nys:nyn,nxl:nxr,0), (/(nx+1)*(ny+1)/) )
 #endif
         IF ( plant_canopy )  THEN
@@ -6884 +6880 @@
             DO i = nxl, nxr
                 DO j = nys, nyn
-                    k = get_topography_top_index_ji( j, i, 's' )
+                    k = topo_top_ind(j,i,0)
 
                     sub_lad(k:nz_plant_t, j, i) = lad_s(0:nz_plant_t-k, j, i)
@@ -10143 +10139 @@
     LOGICAL      ::  two_d !< flag parameter that indicates 2D variables (horizontal cross sections)
 
+    REAL(wp) ::  fill_value = -999.0_wp    !< value for the _FillValue attribute
+
     REAL(wp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) ::  local_pf !< 
 
@@ -10540 +10538 @@
     LOGICAL      ::  found       !<
 
+    REAL(wp)     ::  fill_value = -999.0_wp    !< value for the _FillValue attribute
+
     REAL(sp), DIMENSION(nxl:nxr,nys:nyn,nzb_do:nzt_do) ::  local_pf !<
 
@@ -11107 +11107 @@
     CHARACTER (LEN=*) ::  variable   !< 
 
+    CHARACTER(LEN=5) ::  grid        !< flag to distinquish between staggered grids
+
     INTEGER(iwp) ::  av              !< 
     INTEGER(iwp) ::  i               !< 
     INTEGER(iwp) ::  j               !< 
-    INTEGER(iwp) ::  k               !<
-    INTEGER(iwp) ::  im              !< loop index for masked variables
-    INTEGER(iwp) ::  jm              !< loop index for masked variables
-    INTEGER(iwp) ::  kk              !< 
+    INTEGER(iwp) ::  k               !< 
     INTEGER(iwp) ::  mid             !< masked output running index
-    INTEGER(iwp) ::  ktt             !< k index of highest terrain surface
+    INTEGER(iwp) ::  topo_top_index  !< k index of highest horizontal surface
 
     LOGICAL ::  found                !< true if output array was found
@@ -11127 +11126 @@
     REAL(wp), DIMENSION(:,:,:), POINTER ::  to_be_resorted  !< points to array which needs to be resorted for output
 
+
+    found    = .TRUE.
+    grid     = 's'
     resorted = .FALSE.
-    found = .TRUE.
 
     SELECT CASE ( TRIM( variable ) )
@@ -11215 +11216 @@
              DO  j = 1, mask_size_l(mid,2)
 !
-!--             Get k index of the highest terraing surface
-                im = mask_i(mid,i)
-                jm = mask_j(mid,j)
-                ktt = MINLOC( MERGE( 1, 0, BTEST( wall_flags_0(:,jm,im), 5 )), &
-                              DIM = 1 ) - 1
+!--             Get k index of highest horizontal surface
+                topo_top_index = topo_top_ind(mask_j(mid,j), &
+                                              mask_i(mid,i),   &
+                                              0 )
+!
+!--             Save output array
                 DO  k = 1, mask_size_l(mid,3)
-                   kk = MIN( ktt+mask_k(mid,k), nzt+1 )
-!
-!--                Set value if not in building
-                   IF ( BTEST( wall_flags_0(kk,jm,im), 6 ) )  THEN
-                      local_pf(i,j,k) = fill_value
-                   ELSE
-                      local_pf(i,j,k) = to_be_resorted(kk,jm,im)
-                   ENDIF
+                   local_pf(i,j,k) = to_be_resorted(                         &
+                                          MIN( topo_top_index+mask_k(mid,k), &
+                                               nzt+1 ),                      &
+                                          mask_j(mid,j),                     &
+                                          mask_i(mid,i)                     )
                 ENDDO
              ENDDO