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.
icecor.F90 in NEMO/branches/2019/dev_r11943_MERGE_2019/src/ICE – NEMO

source: NEMO/branches/2019/dev_r11943_MERGE_2019/src/ICE/icecor.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: 10.2 KB
Line 
1MODULE icecor
2   !!======================================================================
3   !!                     ***  MODULE  icecor  ***
4   !!   sea-ice: Corrections on sea-ice variables at the end of the time step
5   !!======================================================================
6   !! History :  3.0  !  2006-04  (M. Vancoppenolle) Original code
7   !!            3.5  !  2014-06  (C. Rousset)       Complete rewriting/cleaning
8   !!            4.0  !  2018     (many people)      SI3 [aka Sea Ice cube]
9   !!----------------------------------------------------------------------
10#if defined key_si3
11   !!----------------------------------------------------------------------
12   !!   'key_si3'                                       SI3 sea-ice model
13   !!----------------------------------------------------------------------
14   !!    ice_cor      : corrections on sea-ice variables
15   !!----------------------------------------------------------------------
16   USE dom_oce        ! ocean domain
17   USE phycst         ! physical constants
18   USE ice            ! sea-ice: variable
19   USE ice1D          ! sea-ice: thermodynamic variables
20   USE iceitd         ! sea-ice: rebining
21   USE icevar         ! sea-ice: operations
22   USE icectl         ! sea-ice: control prints
23   !
24   USE in_out_manager ! I/O manager
25   USE iom            ! I/O manager library
26   USE lib_mpp        ! MPP library
27   USE lib_fortran    ! fortran utilities (glob_sum + no signed zero)
28   USE lbclnk         ! lateral boundary conditions (or mpp links)
29   USE timing         ! Timing
30
31   IMPLICIT NONE
32   PRIVATE
33
34   PUBLIC   ice_cor   ! called by icestp.F90
35
36   !! * Substitutions
37#  include "vectopt_loop_substitute.h90"
38#  include "do_loop_substitute.h90"
39   !!----------------------------------------------------------------------
40   !! NEMO/ICE 4.0 , NEMO Consortium (2018)
41   !! $Id$
42   !! Software governed by the CeCILL license (see ./LICENSE)
43   !!----------------------------------------------------------------------
44CONTAINS
45
46   SUBROUTINE ice_cor( kt, kn )
47      !!----------------------------------------------------------------------
48      !!               ***  ROUTINE ice_cor  ***
49      !!               
50      !! ** Purpose :   Computes corrections on sea-ice global variables at
51      !!              the end of the dynamics (kn=1) and thermodynamics (kn=2)
52      !!----------------------------------------------------------------------
53      INTEGER, INTENT(in) ::   kt    ! number of iteration
54      INTEGER, INTENT(in) ::   kn    ! 1 = after dyn ; 2 = after thermo
55      !
56      INTEGER  ::   ji, jj, jk, jl   ! dummy loop indices
57      REAL(wp) ::   zsal, zzc
58      REAL(wp), DIMENSION(jpi,jpj) ::   zafx   ! concentration trends diag
59      !!----------------------------------------------------------------------
60      ! controls
61      IF( ln_timing    )   CALL timing_start('icecor')                                                             ! timing
62      IF( ln_icediachk )   CALL ice_cons_hsm(0, 'icecor', rdiag_v, rdiag_s, rdiag_t, rdiag_fv, rdiag_fs, rdiag_ft) ! conservation
63      IF( ln_icediachk )   CALL ice_cons2D  (0, 'icecor',  diag_v,  diag_s,  diag_t,  diag_fv,  diag_fs,  diag_ft) ! conservation
64      !
65      IF( kt == nit000 .AND. lwp .AND. kn == 2 ) THEN
66         WRITE(numout,*)
67         WRITE(numout,*) 'ice_cor:  correct sea ice variables if out of bounds ' 
68         WRITE(numout,*) '~~~~~~~'
69      ENDIF
70      !                             !-----------------------------------------------------
71      !                             !  ice thickness must exceed himin (for temp. diff.) !
72      !                             !-----------------------------------------------------
73      WHERE( a_i(:,:,:) >= epsi20 )   ;   h_i(:,:,:) = v_i(:,:,:) / a_i(:,:,:)
74      ELSEWHERE                       ;   h_i(:,:,:) = 0._wp
75      END WHERE
76      WHERE( h_i(:,:,:) < rn_himin )      a_i(:,:,:) = a_i(:,:,:) * h_i(:,:,:) / rn_himin
77      !
78      !                             !-----------------------------------------------------
79      !                             !  ice concentration should not exceed amax          !
80      !                             !-----------------------------------------------------
81      at_i(:,:) = SUM( a_i(:,:,:), dim=3 )
82      DO jl = 1, jpl
83         WHERE( at_i(:,:) > rn_amax_2d(:,:) )   a_i(:,:,jl) = a_i(:,:,jl) * rn_amax_2d(:,:) / at_i(:,:)
84      END DO
85   
86      !                             !-----------------------------------------------------
87      IF ( nn_icesal == 2 ) THEN    !  salinity must stay in bounds [Simin,Simax]        !
88         !                          !-----------------------------------------------------
89         zzc = rhoi * r1_rdtice
90         DO jl = 1, jpl
91            DO_2D_11_11
92               zsal = sv_i(ji,jj,jl)
93               sv_i(ji,jj,jl) = MIN(  MAX( rn_simin*v_i(ji,jj,jl) , sv_i(ji,jj,jl) ) , rn_simax*v_i(ji,jj,jl)  )
94               sfx_res(ji,jj) = sfx_res(ji,jj) - ( sv_i(ji,jj,jl) - zsal ) * zzc   ! associated salt flux
95            END_2D
96         END DO
97      ENDIF
98      !                             !-----------------------------------------------------
99      !                             !  Rebin categories with thickness out of bounds     !
100      !                             !-----------------------------------------------------
101      IF ( jpl > 1 )   CALL ice_itd_reb( kt )
102
103      !                             !-----------------------------------------------------
104      CALL ice_var_zapsmall         !  Zap small values                                  !
105      !                             !-----------------------------------------------------
106
107      !                             !-----------------------------------------------------
108      IF( kn == 2 ) THEN            !  Ice drift case: Corrections to avoid wrong values !
109         DO_2D_00_00
110            IF ( at_i(ji,jj) == 0._wp ) THEN    ! what to do if there is no ice
111               IF ( at_i(ji+1,jj) == 0._wp )   u_ice(ji  ,jj) = 0._wp   ! right side
112               IF ( at_i(ji-1,jj) == 0._wp )   u_ice(ji-1,jj) = 0._wp   ! left side
113               IF ( at_i(ji,jj+1) == 0._wp )   v_ice(ji,jj  ) = 0._wp   ! upper side
114               IF ( at_i(ji,jj-1) == 0._wp )   v_ice(ji,jj-1) = 0._wp   ! bottom side
115            ENDIF
116         END_2D
117         CALL lbc_lnk_multi( 'icecor', u_ice, 'U', -1., v_ice, 'V', -1. )
118      ENDIF
119
120      !                             !-----------------------------------------------------
121      SELECT CASE( kn )             !  Diagnostics                                       !
122      !                             !-----------------------------------------------------
123      CASE( 1 )                        !--- dyn trend diagnostics
124         !
125         IF( ln_icediachk .OR. iom_use('hfxdhc') ) THEN
126            diag_heat(:,:) = - SUM(SUM( e_i (:,:,1:nlay_i,:) - e_i_b (:,:,1:nlay_i,:), dim=4 ), dim=3 ) * r1_rdtice &      ! W.m-2
127               &             - SUM(SUM( e_s (:,:,1:nlay_s,:) - e_s_b (:,:,1:nlay_s,:), dim=4 ), dim=3 ) * r1_rdtice
128            diag_sice(:,:) =   SUM(     sv_i(:,:,:)          - sv_i_b(:,:,:)                  , dim=3 ) * r1_rdtice * rhoi
129            diag_vice(:,:) =   SUM(     v_i (:,:,:)          - v_i_b (:,:,:)                  , dim=3 ) * r1_rdtice * rhoi
130            diag_vsnw(:,:) =   SUM(     v_s (:,:,:)          - v_s_b (:,:,:)                  , dim=3 ) * r1_rdtice * rhos
131         ENDIF
132         !                       ! concentration tendency (dynamics)
133         IF( iom_use('afxdyn') .OR. iom_use('afxthd') .OR. iom_use('afxtot') ) THEN
134            zafx(:,:) = SUM( a_i(:,:,:) - a_i_b(:,:,:), dim=3 ) * r1_rdtice 
135            CALL iom_put( 'afxdyn' , zafx )
136         ENDIF
137         !
138      CASE( 2 )                        !--- thermo trend diagnostics & ice aging
139         !
140         oa_i(:,:,:) = oa_i(:,:,:) + a_i(:,:,:) * rdt_ice   ! ice natural aging incrementation
141         !
142         IF( ln_icediachk .OR. iom_use('hfxdhc') ) THEN
143            diag_heat(:,:) = diag_heat(:,:) &
144               &             - SUM(SUM( e_i (:,:,1:nlay_i,:) - e_i_b (:,:,1:nlay_i,:), dim=4 ), dim=3 ) * r1_rdtice &
145               &             - SUM(SUM( e_s (:,:,1:nlay_s,:) - e_s_b (:,:,1:nlay_s,:), dim=4 ), dim=3 ) * r1_rdtice
146            diag_sice(:,:) = diag_sice(:,:) &
147               &             + SUM(     sv_i(:,:,:)          - sv_i_b(:,:,:)                  , dim=3 ) * r1_rdtice * rhoi
148            diag_vice(:,:) = diag_vice(:,:) &
149               &             + SUM(     v_i (:,:,:)          - v_i_b (:,:,:)                  , dim=3 ) * r1_rdtice * rhoi
150            diag_vsnw(:,:) = diag_vsnw(:,:) &
151               &             + SUM(     v_s (:,:,:)          - v_s_b (:,:,:)                  , dim=3 ) * r1_rdtice * rhos
152            CALL iom_put ( 'hfxdhc' , diag_heat ) 
153         ENDIF
154         !                       ! concentration tendency (total + thermo)
155         IF( iom_use('afxdyn') .OR. iom_use('afxthd') .OR. iom_use('afxtot') ) THEN
156            zafx(:,:) = zafx(:,:) + SUM( a_i(:,:,:) - a_i_b(:,:,:), dim=3 ) * r1_rdtice
157            CALL iom_put( 'afxthd' , SUM( a_i(:,:,:) - a_i_b(:,:,:), dim=3 ) * r1_rdtice )
158            CALL iom_put( 'afxtot' , zafx )
159         ENDIF
160         !
161      END SELECT
162      !
163      ! controls
164      IF( sn_cfctl%l_prtctl ) &
165         &                 CALL ice_prt3D   ('icecor')                                                             ! prints
166      IF( ln_icectl .AND. kn == 2 ) &
167         &                 CALL ice_prt     ( kt, iiceprt, jiceprt, 2, ' - Final state - ' )                       ! prints
168      IF( ln_icediachk )   CALL ice_cons_hsm(1, 'icecor', rdiag_v, rdiag_s, rdiag_t, rdiag_fv, rdiag_fs, rdiag_ft) ! conservation
169      IF( ln_icediachk )   CALL ice_cons2D  (1, 'icecor',  diag_v,  diag_s,  diag_t,  diag_fv,  diag_fs,  diag_ft) ! conservation
170      IF( ln_timing    )   CALL timing_stop ('icecor')                                                             ! timing
171      !
172   END SUBROUTINE ice_cor
173
174#else
175   !!----------------------------------------------------------------------
176   !!   Default option           Dummy module         NO SI3 sea-ice model
177   !!----------------------------------------------------------------------
178#endif
179
180   !!======================================================================
181END MODULE icecor
Note: See TracBrowser for help on using the repository browser.