New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
crsfld.F90 in NEMO/branches/2019/dev_r11943_MERGE_2019/src/OCE/CRS – NEMO

source: NEMO/branches/2019/dev_r11943_MERGE_2019/src/OCE/CRS/crsfld.F90 @ 12340

Last change on this file since 12340 was 12340, checked in by acc, 4 years ago

Branch 2019/dev_r11943_MERGE_2019. This commit introduces basic do loop macro
substitution to the 2019 option 1, merge branch. These changes have been SETTE
tested. The only addition is the do_loop_substitute.h90 file in the OCE directory but
the macros defined therein are used throughout the code to replace identifiable, 2D-
and 3D- nested loop opening and closing statements with single-line alternatives. Code
indents are also adjusted accordingly.

The following explanation is taken from comments in the new header file:

This header file contains preprocessor definitions and macros used in the do-loop
substitutions introduced between version 4.0 and 4.2. The primary aim of these macros
is to assist in future applications of tiling to improve performance. This is expected
to be achieved by alternative versions of these macros in selected locations. The
initial introduction of these macros simply replaces all identifiable nested 2D- and
3D-loops with single line statements (and adjusts indenting accordingly). Do loops
are identifiable if they comform to either:

DO jk = ....

DO jj = .... DO jj = ...

DO ji = .... DO ji = ...
. OR .
. .

END DO END DO

END DO END DO

END DO

and white-space variants thereof.

Additionally, only loops with recognised jj and ji loops limits are treated; these are:
Lower limits of 1, 2 or fs_2
Upper limits of jpi, jpim1 or fs_jpim1 (for ji) or jpj, jpjm1 or fs_jpjm1 (for jj)

The macro naming convention takes the form: DO_2D_BT_LR where:

B is the Bottom offset from the PE's inner domain;
T is the Top offset from the PE's inner domain;
L is the Left offset from the PE's inner domain;
R is the Right offset from the PE's inner domain

So, given an inner domain of 2,jpim1 and 2,jpjm1, a typical example would replace:

DO jj = 2, jpj

DO ji = 1, jpim1
.
.

END DO

END DO

with:

DO_2D_01_10
.
.
END_2D

similar conventions apply to the 3D loops macros. jk loop limits are retained
through macro arguments and are not restricted. This includes the possibility of
strides for which an extra set of DO_3DS macros are defined.

In the example definition below the inner PE domain is defined by start indices of
(kIs, kJs) and end indices of (kIe, KJe)

#define DO_2D_00_00 DO jj = kJs, kJe ; DO ji = kIs, kIe
#define END_2D END DO ; END DO

TO DO:


Only conventional nested loops have been identified and replaced by this step. There are constructs such as:

DO jk = 2, jpkm1

z2d(:,:) = z2d(:,:) + e3w(:,:,jk,Kmm) * z3d(:,:,jk) * wmask(:,:,jk)

END DO

which may need to be considered.

  • Property svn:keywords set to Id
File size: 11.4 KB
Line 
1MODULE crsfld
2   !!======================================================================
3   !!                     ***  MODULE  crsdfld  ***
4   !!  Ocean coarsening :  coarse ocean fields
5   !!=====================================================================
6   !!   2012-07  (J. Simeon, C. Calone, G. Madec, C. Ethe)
7   !!----------------------------------------------------------------------
8
9   !!----------------------------------------------------------------------
10   !!   crs_fld       : create the standard output files for coarse grid and prep
11   !!                       other variables needed to be passed to TOP
12   !!----------------------------------------------------------------------
13   USE crs
14   USE crsdom
15   USE crslbclnk
16   USE oce             ! ocean dynamics and tracers
17   USE dom_oce         ! ocean space and time domain
18   USE sbc_oce         ! Surface boundary condition: ocean fields
19   USE zdf_oce         ! vertical  physics: ocean fields
20   USE ldftra          ! ocean active tracers: lateral diffusivity & EIV coefficients
21   USE zdfddm          ! vertical  physics: double diffusion
22   !
23   USE in_out_manager  ! I/O manager
24   USE iom             !
25   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
26   USE timing          ! preformance summary
27
28   IMPLICIT NONE
29   PRIVATE
30
31   PUBLIC   crs_fld                 ! routines called by step.F90
32
33   !! * Substitutions
34#  include "vectopt_loop_substitute.h90"
35#  include "do_loop_substitute.h90"
36   !!----------------------------------------------------------------------
37   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
38   !! $Id$
39   !! Software governed by the CeCILL license (see ./LICENSE)
40   !!----------------------------------------------------------------------
41CONTAINS
42
43   SUBROUTINE crs_fld( kt, Kmm )
44      !!---------------------------------------------------------------------
45      !!                  ***  ROUTINE crs_fld  ***
46      !!                   
47      !! ** Purpose :   Basic output of coarsened dynamics and tracer fields
48      !!      NETCDF format is used by default
49      !!      1. Accumulate in time the dimensionally-weighted fields
50      !!      2. At time of output, rescale [1] by dimension and time
51      !!         to yield the spatial and temporal average.
52      !!  See. sbcmod.F90
53      !!
54      !! ** Method  : 
55      !!----------------------------------------------------------------------
56      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
57      INTEGER, INTENT(in) ::   Kmm  ! time level index
58      !
59      INTEGER  ::   ji, jj, jk        ! dummy loop indices
60      REAL(wp) ::   z2dcrsu, z2dcrsv  ! local scalars
61      REAL(wp) ::   zztmp             !   -      -
62      !
63      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   ze3t, ze3u, ze3v, ze3w   ! 3D workspace for e3
64      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zt  , zs  , z3d
65      REAL(wp), DIMENSION(jpi_crs,jpj_crs,jpk) ::   zt_crs, zs_crs 
66      !!----------------------------------------------------------------------
67      !
68      IF( ln_timing )   CALL timing_start('crs_fld')
69
70      ! Depth work arrrays
71      ze3t(:,:,:) = e3t(:,:,:,Kmm)
72      ze3u(:,:,:) = e3u(:,:,:,Kmm)
73      ze3v(:,:,:) = e3v(:,:,:,Kmm)
74      ze3w(:,:,:) = e3w(:,:,:,Kmm)
75
76      IF( kt == nit000  ) THEN
77         tsn_crs  (:,:,:,:) = 0._wp    ! temp/sal  array, now
78         un_crs   (:,:,:  ) = 0._wp    ! u-velocity
79         vn_crs   (:,:,:  ) = 0._wp    ! v-velocity
80         wn_crs   (:,:,:  ) = 0._wp    ! w
81         avs_crs  (:,:,:  ) = 0._wp    ! avt
82         hdivn_crs(:,:,:  ) = 0._wp    ! hdiv
83         sshn_crs (:,:    ) = 0._wp    ! ssh
84         utau_crs (:,:    ) = 0._wp    ! taux
85         vtau_crs (:,:    ) = 0._wp    ! tauy
86         wndm_crs (:,:    ) = 0._wp    ! wind speed
87         qsr_crs  (:,:    ) = 0._wp    ! qsr
88         emp_crs  (:,:    ) = 0._wp    ! emp
89         emp_b_crs(:,:    ) = 0._wp    ! emp
90         rnf_crs  (:,:    ) = 0._wp    ! runoff
91         fr_i_crs (:,:    ) = 0._wp    ! ice cover
92      ENDIF
93
94      CALL iom_swap( "nemo_crs" )    ! swap on the coarse grid
95
96      ! 2. Coarsen fields at each time step
97      ! --------------------------------------------------------
98
99      !  Temperature
100      zt(:,:,:) = ts(:,:,:,jp_tem,Kmm)  ;      zt_crs(:,:,:) = 0._wp
101      CALL crs_dom_ope( zt, 'VOL', 'T', tmask, zt_crs, p_e12=e1e2t, p_e3=ze3t, psgn=1.0 )
102      tsn_crs(:,:,:,jp_tem) = zt_crs(:,:,:)
103
104      CALL iom_put( "toce", tsn_crs(:,:,:,jp_tem) )    ! temp
105      CALL iom_put( "sst" , tsn_crs(:,:,1,jp_tem) )    ! sst
106
107     
108      !  Salinity
109      zs(:,:,:) = ts(:,:,:,jp_sal,Kmm)  ;      zs_crs(:,:,:) = 0._wp
110      CALL crs_dom_ope( zs, 'VOL', 'T', tmask, zs_crs, p_e12=e1e2t, p_e3=ze3t, psgn=1.0 )
111      tsn_crs(:,:,:,jp_sal) = zt_crs(:,:,:)
112
113      CALL iom_put( "soce" , tsn_crs(:,:,:,jp_sal) )    ! sal
114      CALL iom_put( "sss"  , tsn_crs(:,:,1,jp_sal) )    ! sss
115
116      !  U-velocity
117      CALL crs_dom_ope( uu(:,:,:,Kmm), 'SUM', 'U', umask, un_crs, p_e12=e2u, p_e3=ze3u, p_surf_crs=e2e3u_msk, psgn=-1.0 )
118      !
119      zt(:,:,:) = 0._wp     ;    zs(:,:,:) = 0._wp  ;   zt_crs(:,:,:) = 0._wp   ;    zs_crs(:,:,:) = 0._wp
120      DO_3D_00_00( 1, jpkm1 )
121         zt(ji,jj,jk)  = uu(ji,jj,jk,Kmm) * 0.5 * ( ts(ji,jj,jk,jp_tem,Kmm) + ts(ji+1,jj,jk,jp_tem,Kmm) ) 
122         zs(ji,jj,jk)  = uu(ji,jj,jk,Kmm) * 0.5 * ( ts(ji,jj,jk,jp_sal,Kmm) + ts(ji+1,jj,jk,jp_sal,Kmm) ) 
123      END_3D
124      CALL crs_dom_ope( zt, 'SUM', 'U', umask, zt_crs, p_e12=e2u, p_e3=ze3u, p_surf_crs=e2e3u_msk, psgn=-1.0 )
125      CALL crs_dom_ope( zs, 'SUM', 'U', umask, zs_crs, p_e12=e2u, p_e3=ze3u, p_surf_crs=e2e3u_msk, psgn=-1.0 )
126
127      CALL iom_put( "uoce"  , un_crs )   ! i-current
128      CALL iom_put( "uocet" , zt_crs )   ! uT
129      CALL iom_put( "uoces" , zs_crs )   ! uS
130
131      !  V-velocity
132      CALL crs_dom_ope( vv(:,:,:,Kmm), 'SUM', 'V', vmask, vn_crs, p_e12=e1v, p_e3=ze3v, p_surf_crs=e1e3v_msk, psgn=-1.0 )
133      !                                                                                 
134      zt(:,:,:) = 0._wp     ;    zs(:,:,:) = 0._wp  ;   zt_crs(:,:,:) = 0._wp   ;    zs_crs(:,:,:) = 0._wp
135      DO_3D_00_00( 1, jpkm1 )
136         zt(ji,jj,jk)  = vv(ji,jj,jk,Kmm) * 0.5 * ( ts(ji,jj,jk,jp_tem,Kmm) + ts(ji,jj+1,jk,jp_tem,Kmm) ) 
137         zs(ji,jj,jk)  = vv(ji,jj,jk,Kmm) * 0.5 * ( ts(ji,jj,jk,jp_sal,Kmm) + ts(ji,jj+1,jk,jp_sal,Kmm) ) 
138      END_3D
139      CALL crs_dom_ope( zt, 'SUM', 'V', vmask, zt_crs, p_e12=e1v, p_e3=ze3v, p_surf_crs=e1e3v_msk, psgn=-1.0 )
140      CALL crs_dom_ope( zs, 'SUM', 'V', vmask, zs_crs, p_e12=e1v, p_e3=ze3v, p_surf_crs=e1e3v_msk, psgn=-1.0 )
141 
142      CALL iom_put( "voce"  , vn_crs )   ! i-current
143      CALL iom_put( "vocet" , zt_crs )   ! vT
144      CALL iom_put( "voces" , zs_crs )   ! vS
145
146      IF( iom_use( "eken") ) THEN     !      kinetic energy
147         z3d(:,:,jk) = 0._wp 
148         DO_3D_00_00( 1, jpkm1 )
149            zztmp  = r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm)
150            z3d(ji,jj,jk) = 0.25_wp * zztmp * (                                    &
151               &            uu(ji-1,jj,jk,Kmm)**2 * e2u(ji-1,jj) * e3u(ji-1,jj,jk,Kmm)   &
152               &          + uu(ji  ,jj,jk,Kmm)**2 * e2u(ji  ,jj) * e3u(ji  ,jj,jk,Kmm)   &
153               &          + vv(ji,jj-1,jk,Kmm)**2 * e1v(ji,jj-1) * e3v(ji,jj-1,jk,Kmm)   &
154               &          + vv(ji,jj  ,jk,Kmm)**2 * e1v(ji,jj  ) * e3v(ji,jj  ,jk,Kmm)   )
155         END_3D
156         CALL lbc_lnk( 'crsfld', z3d, 'T', 1. )
157         !
158         CALL crs_dom_ope( z3d, 'VOL', 'T', tmask, zt_crs, p_e12=e1e2t, p_e3=ze3t, psgn=1.0 )
159         CALL iom_put( "eken", zt_crs )
160      ENDIF
161      !  Horizontal divergence ( following OCE/DYN/divhor.F90 )
162      DO jk = 1, jpkm1
163         DO ji = 2, jpi_crsm1
164            DO jj = 2, jpj_crsm1
165               IF( tmask_crs(ji,jj,jk ) > 0 ) THEN
166                   z2dcrsu =  ( un_crs(ji  ,jj  ,jk) * crs_surfu_wgt(ji  ,jj  ,jk) ) &
167                      &     - ( un_crs(ji-1,jj  ,jk) * crs_surfu_wgt(ji-1,jj  ,jk) )
168                   z2dcrsv =  ( vn_crs(ji  ,jj  ,jk) * crs_surfv_wgt(ji  ,jj  ,jk) ) &
169                      &     - ( vn_crs(ji  ,jj-1,jk) * crs_surfv_wgt(ji  ,jj-1,jk) )
170                   !
171                   hdivn_crs(ji,jj,jk) = ( z2dcrsu + z2dcrsv ) / crs_volt_wgt(ji,jj,jk) 
172               ENDIF
173            END DO
174         END DO
175      END DO
176      CALL crs_lbc_lnk( hdivn_crs, 'T', 1.0 )
177      !
178      CALL iom_put( "hdiv", hdivn_crs ) 
179
180
181      !  W-velocity
182      IF( ln_crs_wn ) THEN
183         CALL crs_dom_ope( ww, 'SUM', 'W', tmask, wn_crs, p_e12=e1e2t, p_surf_crs=e1e2w_msk, psgn=1.0 )
184       !  CALL crs_dom_ope( ww, 'VOL', 'W', tmask, wn_crs, p_e12=e1e2t, p_e3=ze3w )
185      ELSE
186        wn_crs(:,:,jpk) = 0._wp
187        DO jk = jpkm1, 1, -1
188           wn_crs(:,:,jk) = wn_crs(:,:,jk+1) - e3t_crs(:,:,jk) * hdivn_crs(:,:,jk)
189        ENDDO
190      ENDIF
191      CALL iom_put( "woce", wn_crs  )   ! vertical velocity
192      !  free memory
193
194      !  avs
195      SELECT CASE ( nn_crs_kz )
196         CASE ( 0 )
197            CALL crs_dom_ope( avt, 'VOL', 'W', tmask, avt_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 )
198            CALL crs_dom_ope( avs, 'VOL', 'W', tmask, avs_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 )
199         CASE ( 1 )
200            CALL crs_dom_ope( avt, 'MAX', 'W', tmask, avt_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 )
201            CALL crs_dom_ope( avs, 'MAX', 'W', tmask, avs_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 )
202         CASE ( 2 )
203            CALL crs_dom_ope( avt, 'MIN', 'W', tmask, avt_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 )
204            CALL crs_dom_ope( avs, 'MIN', 'W', tmask, avs_crs, p_e12=e1e2t, p_e3=ze3w, psgn=1.0 )
205      END SELECT
206      !
207      CALL iom_put( "avt", avt_crs )   !  Kz on T
208      CALL iom_put( "avs", avs_crs )   !  Kz on S
209     
210      !  sbc fields 
211      CALL crs_dom_ope( ssh(:,:,Kmm) , 'VOL', 'T', tmask, sshn_crs , p_e12=e1e2t, p_e3=ze3t           , psgn=1.0 ) 
212      CALL crs_dom_ope( utau , 'SUM', 'U', umask, utau_crs , p_e12=e2u  , p_surf_crs=e2u_crs  , psgn=1.0 )
213      CALL crs_dom_ope( vtau , 'SUM', 'V', vmask, vtau_crs , p_e12=e1v  , p_surf_crs=e1v_crs  , psgn=1.0 )
214      CALL crs_dom_ope( wndm , 'SUM', 'T', tmask, wndm_crs , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 )
215      CALL crs_dom_ope( rnf  , 'MAX', 'T', tmask, rnf_crs                                     , psgn=1.0 )
216      CALL crs_dom_ope( qsr  , 'SUM', 'T', tmask, qsr_crs  , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 )
217      CALL crs_dom_ope( emp_b, 'SUM', 'T', tmask, emp_b_crs, p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 )
218      CALL crs_dom_ope( emp  , 'SUM', 'T', tmask, emp_crs  , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 )
219      CALL crs_dom_ope( sfx  , 'SUM', 'T', tmask, sfx_crs  , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 )
220      CALL crs_dom_ope( fr_i , 'SUM', 'T', tmask, fr_i_crs , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 )
221
222      CALL iom_put( "ssh"      , sshn_crs )   ! ssh output
223      CALL iom_put( "utau"     , utau_crs )   ! i-tau output
224      CALL iom_put( "vtau"     , vtau_crs )   ! j-tau output
225      CALL iom_put( "wspd"     , wndm_crs )   ! wind speed output
226      CALL iom_put( "runoffs"  , rnf_crs  )   ! runoff output
227      CALL iom_put( "qsr"      , qsr_crs  )   ! qsr output
228      CALL iom_put( "empmr"    , emp_crs  )   ! water flux output
229      CALL iom_put( "saltflx"  , sfx_crs  )   ! salt flux output
230      CALL iom_put( "ice_cover", fr_i_crs )   ! ice cover output
231
232      !
233      CALL iom_swap( "nemo" )     ! return back on high-resolution grid
234      !
235      IF( ln_timing )   CALL timing_stop('crs_fld')
236      !
237   END SUBROUTINE crs_fld
238
239   !!======================================================================
240END MODULE crsfld
Note: See TracBrowser for help on using the repository browser.