source: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/iceadv_prather.F90 @ 8506

Last change on this file since 8506 was 8504, checked in by clem, 3 years ago

changes in style - part4 - clarify ice advection routines

File size: 37.3 KB
Line 
1MODULE iceadv_prather 
2   !!======================================================================
3   !!                       ***  MODULE iceadv_prather   ***
4   !! LIM sea-ice model : sea-ice advection => Prather scheme
5   !!======================================================================
6   !! History :  LIM  ! 2008-03 (M. Vancoppenolle)  LIM-3 from LIM-2 code
7   !!            3.2  ! 2009-06 (F. Dupont)  correct a error in the North fold b.c.
8   !!            4.0  ! 2011-02 (G. Madec) dynamical allocation
9   !!--------------------------------------------------------------------
10#if defined key_lim3
11   !!----------------------------------------------------------------------
12   !!   'key_lim3'                                       LIM3 sea-ice model
13   !!----------------------------------------------------------------------
14   !!   ice_adv_prather     : advection of sea ice using Prather scheme
15   !!----------------------------------------------------------------------
16   USE dom_oce        ! ocean domain
17   USE ice            ! sea-ice variables
18   USE sbc_oce , ONLY : nn_fsbc   ! frequency of sea-ice call
19   !
20   USE lbclnk         ! lateral boundary condition - MPP exchanges
21   USE in_out_manager ! I/O manager
22   USE prtctl         ! Print control
23   USE lib_mpp        ! MPP library
24   USE lib_fortran    ! to use key_nosignedzero
25
26   IMPLICIT NONE
27   PRIVATE
28
29   PUBLIC   ice_adv_prather   ! called by iceadv
30
31   !! * Substitutions
32#  include "vectopt_loop_substitute.h90"
33   !!----------------------------------------------------------------------
34   !! NEMO/ICE 4.0 , NEMO Consortium (2017)
35   !! $Id: iceadv.F90 6746 2016-06-27 17:20:57Z clem $
36   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
37   !!----------------------------------------------------------------------
38CONTAINS
39
40   SUBROUTINE ice_adv_prather( kt, pu_ice, pv_ice,  &
41      &                        pato_i, pv_i, pv_s, psmv_i, poa_i, pa_i, pa_ip, pv_ip, pe_s, pe_i )
42      !!----------------------------------------------------------------------
43      !!                **  routine ice_adv_prather  **
44      !! 
45      !! ** purpose :   Computes and adds the advection trend to sea-ice
46      !!
47      !! ** method  :   Uses Prather second order scheme that advects tracers
48      !!                but also their quadratic forms. The method preserves
49      !!                tracer structures by conserving second order moments.
50      !!
51      !! Reference:  Prather, 1986, JGR, 91, D6. 6671-6681.
52      !!----------------------------------------------------------------------
53      INTEGER                     , INTENT(in   ) ::   kt         ! time step
54      REAL(wp), DIMENSION(:,:)    , INTENT(in   ) ::   pu_ice     ! ice i-velocity
55      REAL(wp), DIMENSION(:,:)    , INTENT(in   ) ::   pv_ice     ! ice j-velocity
56      REAL(wp), DIMENSION(:,:)    , INTENT(inout) ::   pato_i     ! open water area
57      REAL(wp), DIMENSION(:,:,:)  , INTENT(inout) ::   pv_i       ! ice volume
58      REAL(wp), DIMENSION(:,:,:)  , INTENT(inout) ::   pv_s       ! snw volume
59      REAL(wp), DIMENSION(:,:,:)  , INTENT(inout) ::   psmv_i     ! salt content
60      REAL(wp), DIMENSION(:,:,:)  , INTENT(inout) ::   poa_i      ! age content
61      REAL(wp), DIMENSION(:,:,:)  , INTENT(inout) ::   pa_i       ! ice concentration
62      REAL(wp), DIMENSION(:,:,:)  , INTENT(inout) ::   pa_ip      ! melt pond fraction
63      REAL(wp), DIMENSION(:,:,:)  , INTENT(inout) ::   pv_ip      ! melt pond volume
64      REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) ::   pe_s       ! snw heat content
65      REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) ::   pe_i       ! ice heat content
66      !
67      INTEGER  ::   jk, jl, jt              ! dummy loop indices
68      INTEGER  ::   initad                  ! number of sub-timestep for the advection
69      REAL(wp) ::   zcfl , zusnit           !   -      -
70      REAL(wp), ALLOCATABLE, DIMENSION(:,:)     ::   zarea
71      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:)   ::   z0opw
72      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:)   ::   z0ice, z0snw, z0ai, z0es , z0smi , z0oi
73      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:)   ::   z0ap , z0vp
74      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) ::   z0ei
75      !!----------------------------------------------------------------------
76      !
77      IF( kt == nit000 .AND. lwp ) THEN
78         WRITE(numout,*)
79         WRITE(numout,*) 'ice_adv_prather: Prather advection scheme'
80         WRITE(numout,*) '~~~~~~~~~~~~~~~'
81      ENDIF
82      !
83      ALLOCATE( zarea(jpi,jpj)     , z0opw(jpi,jpj, 1 ) ,                                           &
84         &      z0ice(jpi,jpj,jpl) , z0snw(jpi,jpj,jpl) , z0ai(jpi,jpj,jpl) , z0es(jpi,jpj,jpl) ,   &
85         &      z0smi(jpi,jpj,jpl) , z0oi (jpi,jpj,jpl) , z0ap(jpi,jpj,jpl) , z0vp(jpi,jpj,jpl) ,   &
86         &      z0ei (jpi,jpj,nlay_i,jpl)                                                           )
87      !
88      ! --- If ice drift field is too fast, use an appropriate time step for advection (CFL test for stability) --- !       
89      zcfl  =            MAXVAL( ABS( pu_ice(:,:) ) * rdt_ice * r1_e1u(:,:) )
90      zcfl  = MAX( zcfl, MAXVAL( ABS( pv_ice(:,:) ) * rdt_ice * r1_e2v(:,:) ) )
91      IF( lk_mpp )   CALL mpp_max( zcfl )
92     
93      IF( zcfl > 0.5 ) THEN   ;   initad = 2   ;   zusnit = 0.5_wp
94      ELSE                    ;   initad = 1   ;   zusnit = 1.0_wp
95      ENDIF
96     
97      zarea(:,:) = e1e2t(:,:)
98      !-------------------------
99      ! transported fields                                       
100      !-------------------------
101      z0opw(:,:,1) = pato_i(:,:) * e1e2t(:,:)              ! Open water area
102      DO jl = 1, jpl
103         z0snw(:,:,jl) = pv_s  (:,:,  jl) * e1e2t(:,:)     ! Snow volume
104         z0ice(:,:,jl) = pv_i  (:,:,  jl) * e1e2t(:,:)     ! Ice  volume
105         z0ai (:,:,jl) = pa_i  (:,:,  jl) * e1e2t(:,:)     ! Ice area
106         z0smi(:,:,jl) = psmv_i(:,:,  jl) * e1e2t(:,:)     ! Salt content
107         z0oi (:,:,jl) = poa_i (:,:,  jl) * e1e2t(:,:)     ! Age content
108         z0es (:,:,jl) = pe_s  (:,:,1,jl) * e1e2t(:,:)     ! Snow heat content
109         DO jk = 1, nlay_i
110            z0ei(:,:,jk,jl) = pe_i(:,:,jk,jl) * e1e2t(:,:) ! Ice  heat content
111         END DO
112         IF ( nn_pnd_scheme > 0 ) THEN
113            z0ap(:,:,jl)  = pa_ip(:,:,jl) * e1e2t(:,:)     ! Melt pond fraction
114            z0vp(:,:,jl)  = pv_ip(:,:,jl) * e1e2t(:,:)     ! Melt pond volume
115         ENDIF
116      END DO
117
118      !                                                    !--------------------------------------------!
119      IF( MOD( ( kt - 1) / nn_fsbc , 2 ) == 0 ) THEN       !==  odd ice time step:  adv_x then adv_y  ==!
120         !                                                 !--------------------------------------------!
121         DO jt = 1, initad
122            CALL adv_x( zusnit, pu_ice, 1._wp, zarea, z0opw (:,:,1), sxopw(:,:),   &             !--- ice open water area
123               &                                      sxxopw(:,:)  , syopw(:,:), syyopw(:,:), sxyopw(:,:)  )
124            CALL adv_y( zusnit, pv_ice, 0._wp, zarea, z0opw (:,:,1), sxopw(:,:),   &
125               &                                      sxxopw(:,:)  , syopw(:,:), syyopw(:,:), sxyopw(:,:)  )
126            DO jl = 1, jpl
127               CALL adv_x( zusnit, pu_ice, 1._wp, zarea, z0ice (:,:,jl), sxice(:,:,jl),   &    !--- ice volume  ---
128                  &                                      sxxice(:,:,jl), syice(:,:,jl), syyice(:,:,jl), sxyice(:,:,jl)  )
129               CALL adv_y( zusnit, pv_ice, 0._wp, zarea, z0ice (:,:,jl), sxice(:,:,jl),   &
130                  &                                      sxxice(:,:,jl), syice(:,:,jl), syyice(:,:,jl), sxyice(:,:,jl)  )
131               CALL adv_x( zusnit, pu_ice, 1._wp, zarea, z0snw (:,:,jl), sxsn (:,:,jl),   &    !--- snow volume  ---
132                  &                                      sxxsn (:,:,jl), sysn (:,:,jl), syysn (:,:,jl), sxysn (:,:,jl)  )
133               CALL adv_y( zusnit, pv_ice, 0._wp, zarea, z0snw (:,:,jl), sxsn (:,:,jl),   &
134                  &                                      sxxsn (:,:,jl), sysn (:,:,jl), syysn (:,:,jl), sxysn (:,:,jl)  )
135               CALL adv_x( zusnit, pu_ice, 1._wp, zarea, z0smi (:,:,jl), sxsal(:,:,jl),   &    !--- ice salinity ---
136                  &                                      sxxsal(:,:,jl), sysal(:,:,jl), syysal(:,:,jl), sxysal(:,:,jl)  )
137               CALL adv_y( zusnit, pv_ice, 0._wp, zarea, z0smi (:,:,jl), sxsal(:,:,jl),   &
138                  &                                      sxxsal(:,:,jl), sysal(:,:,jl), syysal(:,:,jl), sxysal(:,:,jl)  )
139               CALL adv_x( zusnit, pu_ice, 1._wp, zarea, z0oi  (:,:,jl), sxage(:,:,jl),   &    !--- ice age      ---     
140                  &                                      sxxage(:,:,jl), syage(:,:,jl), syyage(:,:,jl), sxyage(:,:,jl)  )
141               CALL adv_y( zusnit, pv_ice, 0._wp, zarea, z0oi  (:,:,jl), sxage(:,:,jl),   &
142                  &                                      sxxage(:,:,jl), syage(:,:,jl), syyage(:,:,jl), sxyage(:,:,jl)  )
143               CALL adv_x( zusnit, pu_ice, 1._wp, zarea, z0ai  (:,:,jl), sxa  (:,:,jl),   &    !--- ice concentrations ---
144                  &                                      sxxa  (:,:,jl), sya  (:,:,jl), syya  (:,:,jl), sxya  (:,:,jl)  )
145               CALL adv_y( zusnit, pv_ice, 0._wp, zarea, z0ai  (:,:,jl), sxa  (:,:,jl),   & 
146                  &                                      sxxa  (:,:,jl), sya  (:,:,jl), syya  (:,:,jl), sxya  (:,:,jl)  )
147               CALL adv_x( zusnit, pu_ice, 1._wp, zarea, z0es  (:,:,jl), sxc0 (:,:,jl),   &    !--- snow heat contents ---
148                  &                                      sxxc0 (:,:,jl), syc0 (:,:,jl), syyc0 (:,:,jl), sxyc0 (:,:,jl)  )
149               CALL adv_y( zusnit, pv_ice, 0._wp, zarea, z0es  (:,:,jl), sxc0 (:,:,jl),   &
150                  &                                      sxxc0 (:,:,jl), syc0 (:,:,jl), syyc0 (:,:,jl), sxyc0 (:,:,jl)  )
151               DO jk = 1, nlay_i                                                               !--- ice heat contents ---
152                  CALL adv_x( zusnit, pu_ice, 1._wp, zarea, z0ei(:,:,jk,jl), sxe (:,:,jk,jl),   & 
153                     &                                      sxxe(:,:,jk,jl), sye (:,:,jk,jl),   &
154                     &                                      syye(:,:,jk,jl), sxye(:,:,jk,jl) )
155                  CALL adv_y( zusnit, pv_ice, 0._wp, zarea, z0ei(:,:,jk,jl), sxe (:,:,jk,jl),   & 
156                     &                                      sxxe(:,:,jk,jl), sye (:,:,jk,jl),   &
157                     &                                      syye(:,:,jk,jl), sxye(:,:,jk,jl) )
158               END DO
159               IF ( nn_pnd_scheme > 0 ) THEN
160                  CALL adv_x( zusnit, pu_ice, 1._wp, zarea, z0ap  (:,:,jl), sxap (:,:,jl),   &    !--- melt pond fraction --
161                     &                                      sxxap (:,:,jl), syap (:,:,jl), syyap (:,:,jl), sxyap (:,:,jl)  )
162                  CALL adv_y( zusnit, pv_ice, 0._wp, zarea, z0ap  (:,:,jl), sxap (:,:,jl),   & 
163                     &                                      sxxap (:,:,jl), syap (:,:,jl), syyap (:,:,jl), sxyap (:,:,jl)  )
164                  CALL adv_x( zusnit, pu_ice, 1._wp, zarea, z0vp  (:,:,jl), sxvp (:,:,jl),   &    !--- melt pond volume   --
165                     &                                      sxxvp (:,:,jl), syvp (:,:,jl), syyvp (:,:,jl), sxyvp (:,:,jl)  )
166                  CALL adv_y( zusnit, pv_ice, 0._wp, zarea, z0vp  (:,:,jl), sxvp (:,:,jl),   & 
167                     &                                      sxxvp (:,:,jl), syvp (:,:,jl), syyvp (:,:,jl), sxyvp (:,:,jl)  )
168               ENDIF
169            END DO
170         END DO
171      !                                                    !--------------------------------------------!
172      ELSE                                                 !== even ice time step:  adv_y then adv_x  ==!
173         !                                                 !--------------------------------------------!
174         DO jt = 1, initad
175            CALL adv_y( zusnit, pv_ice, 1._wp, zarea, z0opw (:,:,1), sxopw(:,:),   &             !--- ice open water area
176               &                                      sxxopw(:,:)  , syopw(:,:), syyopw(:,:), sxyopw(:,:)  )
177            CALL adv_x( zusnit, pu_ice, 0._wp, zarea, z0opw (:,:,1), sxopw(:,:),   &
178               &                                      sxxopw(:,:)  , syopw(:,:), syyopw(:,:), sxyopw(:,:)  )
179            DO jl = 1, jpl
180               CALL adv_y( zusnit, pv_ice, 1._wp, zarea, z0ice (:,:,jl), sxice(:,:,jl),   &    !--- ice volume  ---
181                  &                                      sxxice(:,:,jl), syice(:,:,jl), syyice(:,:,jl), sxyice(:,:,jl)  )
182               CALL adv_x( zusnit, pu_ice, 0._wp, zarea, z0ice (:,:,jl), sxice(:,:,jl),   &
183                  &                                      sxxice(:,:,jl), syice(:,:,jl), syyice(:,:,jl), sxyice(:,:,jl)  )
184               CALL adv_y( zusnit, pv_ice, 1._wp, zarea, z0snw (:,:,jl), sxsn (:,:,jl),   &    !--- snow volume  ---
185                  &                                      sxxsn (:,:,jl), sysn (:,:,jl), syysn (:,:,jl), sxysn (:,:,jl)  )
186               CALL adv_x( zusnit, pu_ice, 0._wp, zarea, z0snw (:,:,jl), sxsn (:,:,jl),   &
187                  &                                      sxxsn (:,:,jl), sysn (:,:,jl), syysn (:,:,jl), sxysn (:,:,jl)  )
188               CALL adv_y( zusnit, pv_ice, 1._wp, zarea, z0smi (:,:,jl), sxsal(:,:,jl),   &    !--- ice salinity ---
189                  &                                      sxxsal(:,:,jl), sysal(:,:,jl), syysal(:,:,jl), sxysal(:,:,jl)  )
190               CALL adv_x( zusnit, pu_ice, 0._wp, zarea, z0smi (:,:,jl), sxsal(:,:,jl),   &
191                  &                                      sxxsal(:,:,jl), sysal(:,:,jl), syysal(:,:,jl), sxysal(:,:,jl)  )
192               CALL adv_y( zusnit, pv_ice, 1._wp, zarea, z0oi  (:,:,jl), sxage(:,:,jl),   &   !--- ice age      ---
193                  &                                      sxxage(:,:,jl), syage(:,:,jl), syyage(:,:,jl), sxyage(:,:,jl)  )
194               CALL adv_x( zusnit, pu_ice, 0._wp, zarea, z0oi  (:,:,jl), sxage(:,:,jl),   &
195                  &                                      sxxage(:,:,jl), syage(:,:,jl), syyage(:,:,jl), sxyage(:,:,jl)  )
196               CALL adv_y( zusnit, pv_ice, 1._wp, zarea, z0ai  (:,:,jl), sxa  (:,:,jl),   &   !--- ice concentrations ---
197                  &                                      sxxa  (:,:,jl), sya  (:,:,jl), syya  (:,:,jl), sxya  (:,:,jl)  )
198               CALL adv_x( zusnit, pu_ice, 0._wp, zarea, z0ai  (:,:,jl), sxa  (:,:,jl),   &
199                  &                                      sxxa  (:,:,jl), sya  (:,:,jl), syya  (:,:,jl), sxya  (:,:,jl)  )
200               CALL adv_y( zusnit, pv_ice, 1._wp, zarea, z0es  (:,:,jl), sxc0 (:,:,jl),   &  !--- snow heat contents ---
201                  &                                      sxxc0 (:,:,jl), syc0 (:,:,jl), syyc0 (:,:,jl), sxyc0 (:,:,jl)  )
202               CALL adv_x( zusnit, pu_ice, 0._wp, zarea, z0es  (:,:,jl), sxc0 (:,:,jl),   &
203                  &                                      sxxc0 (:,:,jl), syc0 (:,:,jl), syyc0 (:,:,jl), sxyc0 (:,:,jl)  )
204               DO jk = 1, nlay_i                                                             !--- ice heat contents ---
205                  CALL adv_y( zusnit, pv_ice, 1._wp, zarea, z0ei(:,:,jk,jl), sxe (:,:,jk,jl),   & 
206                     &                                      sxxe(:,:,jk,jl), sye (:,:,jk,jl),   &
207                     &                                      syye(:,:,jk,jl), sxye(:,:,jk,jl) )
208                  CALL adv_x( zusnit, pu_ice, 0._wp, zarea, z0ei(:,:,jk,jl), sxe (:,:,jk,jl),   & 
209                     &                                      sxxe(:,:,jk,jl), sye (:,:,jk,jl),   &
210                     &                                      syye(:,:,jk,jl), sxye(:,:,jk,jl) )
211               END DO
212               IF ( nn_pnd_scheme > 0 ) THEN
213                  CALL adv_y( zusnit, pv_ice, 1._wp, zarea, z0ap  (:,:,jl), sxap (:,:,jl),   &   !--- melt pond fraction ---
214                     &                                      sxxap (:,:,jl), syap (:,:,jl), syyap (:,:,jl), sxyap (:,:,jl)  )
215                  CALL adv_x( zusnit, pu_ice, 0._wp, zarea, z0ap  (:,:,jl), sxap (:,:,jl),   &
216                     &                                      sxxap (:,:,jl), syap (:,:,jl), syyap (:,:,jl), sxyap (:,:,jl)  )
217                  CALL adv_y( zusnit, pv_ice, 1._wp, zarea, z0vp  (:,:,jl), sxvp (:,:,jl),   &   !--- melt pond volume   ---
218                     &                                      sxxvp (:,:,jl), syvp (:,:,jl), syyvp (:,:,jl), sxyvp (:,:,jl)  )
219                  CALL adv_x( zusnit, pu_ice, 0._wp, zarea, z0vp  (:,:,jl), sxvp (:,:,jl),   &
220                     &                                      sxxvp (:,:,jl), syvp (:,:,jl), syyvp (:,:,jl), sxyvp (:,:,jl)  )
221               ENDIF
222            END DO
223         END DO
224      ENDIF
225
226      !-------------------------------------------
227      ! Recover the properties from their contents
228      !-------------------------------------------
229      pato_i(:,:) = z0opw(:,:,1) * r1_e1e2t(:,:)
230      DO jl = 1, jpl
231         pv_i  (:,:,  jl) = z0ice(:,:,jl) * r1_e1e2t(:,:)
232         pv_s  (:,:,  jl) = z0snw(:,:,jl) * r1_e1e2t(:,:)
233         psmv_i(:,:,  jl) = z0smi(:,:,jl) * r1_e1e2t(:,:)
234         poa_i (:,:,  jl) = z0oi (:,:,jl) * r1_e1e2t(:,:)
235         pa_i  (:,:,  jl) = z0ai (:,:,jl) * r1_e1e2t(:,:)
236         pe_s  (:,:,1,jl) = z0es (:,:,jl) * r1_e1e2t(:,:)
237         DO jk = 1, nlay_i
238            pe_i(:,:,jk,jl) = z0ei(:,:,jk,jl) * r1_e1e2t(:,:)
239         END DO
240         ! MV MP 2016
241         IF ( nn_pnd_scheme > 0 ) THEN
242            pa_ip  (:,:,jl)   = z0ap (:,:,jl) * r1_e1e2t(:,:)
243            pv_ip  (:,:,jl)   = z0vp (:,:,jl) * r1_e1e2t(:,:)
244         ENDIF
245         ! END MV MP 2016
246      END DO
247      !
248      DEALLOCATE( zarea , z0opw , z0ice, z0snw , z0ai , z0es , z0smi , z0oi , z0ap , z0vp , z0ei )
249      !
250   END SUBROUTINE ice_adv_prather
251   
252   SUBROUTINE adv_x( pdf, put , pcrh, psm , ps0 ,   &
253      &              psx, psxx, psy , psyy, psxy )
254      !!----------------------------------------------------------------------
255      !!                **  routine adv_x  **
256      !! 
257      !! ** purpose :   Computes and adds the advection trend to sea-ice
258      !!                variable on x axis
259      !!----------------------------------------------------------------------
260      REAL(wp)                    , INTENT(in   ) ::   pdf                ! reduction factor for the time step
261      REAL(wp)                    , INTENT(in   ) ::   pcrh               ! call adv_x then adv_y (=1) or the opposite (=0)
262      REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) ::   put                ! i-direction ice velocity at U-point [m/s]
263      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   psm                ! area
264      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   ps0                ! field to be advected
265      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   psx , psy          ! 1st moments
266      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   psxx, psyy, psxy   ! 2nd moments
267      !!
268      INTEGER  ::   ji, jj                               ! dummy loop indices
269      REAL(wp) ::   zs1max, zrdt, zslpmax, ztemp         ! local scalars
270      REAL(wp) ::   zs1new, zalf , zalfq , zbt           !   -      -
271      REAL(wp) ::   zs2new, zalf1, zalf1q, zbt1          !   -      -
272      REAL(wp), DIMENSION(jpi,jpj) ::   zf0 , zfx  , zfy   , zbet   ! 2D workspace
273      REAL(wp), DIMENSION(jpi,jpj) ::   zfm , zfxx , zfyy  , zfxy   !  -      -
274      REAL(wp), DIMENSION(jpi,jpj) ::   zalg, zalg1, zalg1q         !  -      -
275      !-----------------------------------------------------------------------
276
277      ! Limitation of moments.                                           
278
279      zrdt = rdt_ice * pdf      ! If ice drift field is too fast, use an appropriate time step for advection.
280
281      DO jj = 1, jpj
282         DO ji = 1, jpi
283            zslpmax = MAX( 0._wp, ps0(ji,jj) )
284            zs1max  = 1.5 * zslpmax
285            zs1new  = MIN( zs1max, MAX( -zs1max, psx(ji,jj) ) )
286            zs2new  = MIN(  2.0 * zslpmax - 0.3334 * ABS( zs1new ),      &
287               &            MAX( ABS( zs1new ) - zslpmax, psxx(ji,jj) )  )
288            rswitch = ( 1.0 - MAX( 0._wp, SIGN( 1._wp, -zslpmax) ) ) * tmask(ji,jj,1)   ! Case of empty boxes & Apply mask
289
290            ps0 (ji,jj) = zslpmax 
291            psx (ji,jj) = zs1new      * rswitch
292            psxx(ji,jj) = zs2new      * rswitch
293            psy (ji,jj) = psy (ji,jj) * rswitch
294            psyy(ji,jj) = psyy(ji,jj) * rswitch
295            psxy(ji,jj) = MIN( zslpmax, MAX( -zslpmax, psxy(ji,jj) ) ) * rswitch
296         END DO
297      END DO
298
299      !  Initialize volumes of boxes  (=area if adv_x first called, =psm otherwise)                                     
300      psm (:,:)  = MAX( pcrh * e1e2t(:,:) + ( 1.0 - pcrh ) * psm(:,:) , epsi20 )
301
302      !  Calculate fluxes and moments between boxes i<-->i+1             
303      DO jj = 1, jpj                      !  Flux from i to i+1 WHEN u GT 0
304         DO ji = 1, jpi
305            zbet(ji,jj)  =  MAX( 0._wp, SIGN( 1._wp, put(ji,jj) ) )
306            zalf         =  MAX( 0._wp, put(ji,jj) ) * zrdt * e2u(ji,jj) / psm(ji,jj)
307            zalfq        =  zalf * zalf
308            zalf1        =  1.0 - zalf
309            zalf1q       =  zalf1 * zalf1
310            !
311            zfm (ji,jj)  =  zalf  *   psm (ji,jj)
312            zf0 (ji,jj)  =  zalf  * ( ps0 (ji,jj) + zalf1 * ( psx(ji,jj) + (zalf1 - zalf) * psxx(ji,jj) )  )
313            zfx (ji,jj)  =  zalfq * ( psx (ji,jj) + 3.0 * zalf1 * psxx(ji,jj) )
314            zfxx(ji,jj)  =  zalf  *   psxx(ji,jj) * zalfq
315            zfy (ji,jj)  =  zalf  * ( psy (ji,jj) + zalf1 * psxy(ji,jj) )
316            zfxy(ji,jj)  =  zalfq *   psxy(ji,jj)
317            zfyy(ji,jj)  =  zalf  *   psyy(ji,jj)
318
319            !  Readjust moments remaining in the box.
320            psm (ji,jj)  =  psm (ji,jj) - zfm(ji,jj)
321            ps0 (ji,jj)  =  ps0 (ji,jj) - zf0(ji,jj)
322            psx (ji,jj)  =  zalf1q * ( psx(ji,jj) - 3.0 * zalf * psxx(ji,jj) )
323            psxx(ji,jj)  =  zalf1  * zalf1q * psxx(ji,jj)
324            psy (ji,jj)  =  psy (ji,jj) - zfy(ji,jj)
325            psyy(ji,jj)  =  psyy(ji,jj) - zfyy(ji,jj)
326            psxy(ji,jj)  =  zalf1q * psxy(ji,jj)
327         END DO
328      END DO
329
330      DO jj = 1, jpjm1                      !  Flux from i+1 to i when u LT 0.
331         DO ji = 1, fs_jpim1
332            zalf          = MAX( 0._wp, -put(ji,jj) ) * zrdt * e2u(ji,jj) / psm(ji+1,jj) 
333            zalg  (ji,jj) = zalf
334            zalfq         = zalf * zalf
335            zalf1         = 1.0 - zalf
336            zalg1 (ji,jj) = zalf1
337            zalf1q        = zalf1 * zalf1
338            zalg1q(ji,jj) = zalf1q
339            !
340            zfm   (ji,jj) = zfm (ji,jj) + zalf  *   psm (ji+1,jj)
341            zf0   (ji,jj) = zf0 (ji,jj) + zalf  * ( ps0 (ji+1,jj) - zalf1 * ( psx(ji+1,jj) - (zalf1 - zalf ) * psxx(ji+1,jj) ) )
342            zfx   (ji,jj) = zfx (ji,jj) + zalfq * ( psx (ji+1,jj) - 3.0 * zalf1 * psxx(ji+1,jj) )
343            zfxx  (ji,jj) = zfxx(ji,jj) + zalf  *   psxx(ji+1,jj) * zalfq
344            zfy   (ji,jj) = zfy (ji,jj) + zalf  * ( psy (ji+1,jj) - zalf1 * psxy(ji+1,jj) )
345            zfxy  (ji,jj) = zfxy(ji,jj) + zalfq *   psxy(ji+1,jj)
346            zfyy  (ji,jj) = zfyy(ji,jj) + zalf  *   psyy(ji+1,jj)
347         END DO
348      END DO
349
350      DO jj = 2, jpjm1                     !  Readjust moments remaining in the box.
351         DO ji = fs_2, fs_jpim1
352            zbt  =       zbet(ji-1,jj)
353            zbt1 = 1.0 - zbet(ji-1,jj)
354            !
355            psm (ji,jj) = zbt * psm(ji,jj) + zbt1 * ( psm(ji,jj) - zfm(ji-1,jj) )
356            ps0 (ji,jj) = zbt * ps0(ji,jj) + zbt1 * ( ps0(ji,jj) - zf0(ji-1,jj) )
357            psx (ji,jj) = zalg1q(ji-1,jj) * ( psx(ji,jj) + 3.0 * zalg(ji-1,jj) * psxx(ji,jj) )
358            psxx(ji,jj) = zalg1 (ji-1,jj) * zalg1q(ji-1,jj) * psxx(ji,jj)
359            psy (ji,jj) = zbt * psy (ji,jj) + zbt1 * ( psy (ji,jj) - zfy (ji-1,jj) )
360            psyy(ji,jj) = zbt * psyy(ji,jj) + zbt1 * ( psyy(ji,jj) - zfyy(ji-1,jj) )
361            psxy(ji,jj) = zalg1q(ji-1,jj) * psxy(ji,jj)
362         END DO
363      END DO
364
365      !   Put the temporary moments into appropriate neighboring boxes.   
366      DO jj = 2, jpjm1                     !   Flux from i to i+1 IF u GT 0.
367         DO ji = fs_2, fs_jpim1
368            zbt  =       zbet(ji-1,jj)
369            zbt1 = 1.0 - zbet(ji-1,jj)
370            psm(ji,jj)  = zbt * ( psm(ji,jj) + zfm(ji-1,jj) ) + zbt1 * psm(ji,jj)
371            zalf        = zbt * zfm(ji-1,jj) / psm(ji,jj)
372            zalf1       = 1.0 - zalf
373            ztemp       = zalf * ps0(ji,jj) - zalf1 * zf0(ji-1,jj)
374            !
375            ps0 (ji,jj) = zbt * ( ps0(ji,jj) + zf0(ji-1,jj) ) + zbt1 * ps0(ji,jj)
376            psx (ji,jj) = zbt * ( zalf * zfx(ji-1,jj) + zalf1 * psx(ji,jj) + 3.0 * ztemp ) + zbt1 * psx(ji,jj)
377            psxx(ji,jj) = zbt * ( zalf * zalf * zfxx(ji-1,jj) + zalf1 * zalf1 * psxx(ji,jj)                               &
378               &                + 5.0 * ( zalf * zalf1 * ( psx (ji,jj) - zfx(ji-1,jj) ) - ( zalf1 - zalf ) * ztemp )  )   &
379               &                                                + zbt1 * psxx(ji,jj)
380            psxy(ji,jj) = zbt * ( zalf * zfxy(ji-1,jj) + zalf1 * psxy(ji,jj)             &
381               &                + 3.0 * (- zalf1*zfy(ji-1,jj)  + zalf * psy(ji,jj) ) )   &
382               &                                                + zbt1 * psxy(ji,jj)
383            psy (ji,jj) = zbt * ( psy (ji,jj) + zfy (ji-1,jj) ) + zbt1 * psy (ji,jj)
384            psyy(ji,jj) = zbt * ( psyy(ji,jj) + zfyy(ji-1,jj) ) + zbt1 * psyy(ji,jj)
385         END DO
386      END DO
387
388      DO jj = 2, jpjm1                     !  Flux from i+1 to i IF u LT 0.
389         DO ji = fs_2, fs_jpim1
390            zbt  =       zbet(ji,jj)
391            zbt1 = 1.0 - zbet(ji,jj)
392            psm(ji,jj)  = zbt * psm(ji,jj)  + zbt1 * ( psm(ji,jj) + zfm(ji,jj) )
393            zalf        = zbt1 * zfm(ji,jj) / psm(ji,jj)
394            zalf1       = 1.0 - zalf
395            ztemp       = - zalf * ps0(ji,jj) + zalf1 * zf0(ji,jj)
396            !
397            ps0(ji,jj)  = zbt * ps0 (ji,jj) + zbt1 * ( ps0(ji,jj) + zf0(ji,jj) )
398            psx(ji,jj)  = zbt * psx (ji,jj) + zbt1 * ( zalf * zfx(ji,jj) + zalf1 * psx(ji,jj) + 3.0 * ztemp )
399            psxx(ji,jj) = zbt * psxx(ji,jj) + zbt1 * ( zalf * zalf * zfxx(ji,jj)  + zalf1 * zalf1 * psxx(ji,jj)  &
400               &                                      + 5.0 *( zalf * zalf1 * ( - psx(ji,jj) + zfx(ji,jj) )      &
401               &                                      + ( zalf1 - zalf ) * ztemp ) )
402            psxy(ji,jj) = zbt * psxy(ji,jj) + zbt1 * (  zalf * zfxy(ji,jj) + zalf1 * psxy(ji,jj)  &
403               &                                      + 3.0 * ( zalf1 * zfy(ji,jj) - zalf * psy(ji,jj) )  )
404            psy(ji,jj)  = zbt * psy (ji,jj)  + zbt1 * ( psy (ji,jj) + zfy (ji,jj) )
405            psyy(ji,jj) = zbt * psyy(ji,jj)  + zbt1 * ( psyy(ji,jj) + zfyy(ji,jj) )
406         END DO
407      END DO
408
409      !-- Lateral boundary conditions
410      CALL lbc_lnk_multi( psm , 'T',  1., ps0 , 'T',  1.   &
411         &              , psx , 'T', -1., psy , 'T', -1.   &   ! caution gradient ==> the sign changes
412         &              , psxx, 'T',  1., psyy, 'T',  1.   &
413         &              , psxy, 'T',  1. )
414
415      IF(ln_ctl) THEN
416         CALL prt_ctl(tab2d_1=psm  , clinfo1=' adv_x: psm  :', tab2d_2=ps0 , clinfo2=' ps0  : ')
417         CALL prt_ctl(tab2d_1=psx  , clinfo1=' adv_x: psx  :', tab2d_2=psxx, clinfo2=' psxx : ')
418         CALL prt_ctl(tab2d_1=psy  , clinfo1=' adv_x: psy  :', tab2d_2=psyy, clinfo2=' psyy : ')
419         CALL prt_ctl(tab2d_1=psxy , clinfo1=' adv_x: psxy :')
420      ENDIF
421      !
422   END SUBROUTINE adv_x
423
424
425   SUBROUTINE adv_y( pdf, pvt , pcrh, psm , ps0 ,   &
426      &              psx, psxx, psy , psyy, psxy )
427      !!---------------------------------------------------------------------
428      !!                **  routine adv_y  **
429      !!           
430      !! ** purpose :   Computes and adds the advection trend to sea-ice
431      !!                variable on y axis
432      !!---------------------------------------------------------------------
433      REAL(wp)                    , INTENT(in   ) ::   pdf                ! reduction factor for the time step
434      REAL(wp)                    , INTENT(in   ) ::   pcrh               ! call adv_x then adv_y (=1) or the opposite (=0)
435      REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) ::   pvt                ! j-direction ice velocity at V-point [m/s]
436      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   psm                ! area
437      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   ps0                ! field to be advected
438      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   psx , psy          ! 1st moments
439      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   psxx, psyy, psxy   ! 2nd moments
440      !!
441      INTEGER  ::   ji, jj                               ! dummy loop indices
442      REAL(wp) ::   zs1max, zrdt, zslpmax, ztemp         ! temporary scalars
443      REAL(wp) ::   zs1new, zalf , zalfq , zbt           !    -         -
444      REAL(wp) ::   zs2new, zalf1, zalf1q, zbt1          !    -         -
445      REAL(wp), DIMENSION(jpi,jpj) ::   zf0, zfx , zfy , zbet   ! 2D workspace
446      REAL(wp), DIMENSION(jpi,jpj) ::   zfm, zfxx, zfyy, zfxy   !  -      -
447      REAL(wp), DIMENSION(jpi,jpj) ::   zalg, zalg1, zalg1q     !  -      -
448      !---------------------------------------------------------------------
449
450      ! Limitation of moments.
451
452      zrdt = rdt_ice * pdf ! If ice drift field is too fast, use an appropriate time step for advection.
453
454      DO jj = 1, jpj
455         DO ji = 1, jpi
456            zslpmax = MAX( 0._wp, ps0(ji,jj) )
457            zs1max  = 1.5 * zslpmax
458            zs1new  = MIN( zs1max, MAX( -zs1max, psy(ji,jj) ) )
459            zs2new  = MIN(  ( 2.0 * zslpmax - 0.3334 * ABS( zs1new ) ),   &
460               &             MAX( ABS( zs1new )-zslpmax, psyy(ji,jj) )  )
461            rswitch = ( 1.0 - MAX( 0._wp, SIGN( 1._wp, -zslpmax) ) ) * tmask(ji,jj,1)   ! Case of empty boxes & Apply mask
462            !
463            ps0 (ji,jj) = zslpmax 
464            psx (ji,jj) = psx (ji,jj) * rswitch
465            psxx(ji,jj) = psxx(ji,jj) * rswitch
466            psy (ji,jj) = zs1new * rswitch
467            psyy(ji,jj) = zs2new * rswitch
468            psxy(ji,jj) = MIN( zslpmax, MAX( -zslpmax, psxy(ji,jj) ) ) * rswitch
469         END DO
470      END DO
471
472      !  Initialize volumes of boxes (=area if adv_x first called, =psm otherwise)
473      psm(:,:)  = MAX(  pcrh * e1e2t(:,:) + ( 1.0 - pcrh ) * psm(:,:) , epsi20  )
474
475      !  Calculate fluxes and moments between boxes j<-->j+1             
476      DO jj = 1, jpj                     !  Flux from j to j+1 WHEN v GT 0   
477         DO ji = 1, jpi
478            zbet(ji,jj)  =  MAX( 0._wp, SIGN( 1._wp, pvt(ji,jj) ) )
479            zalf         =  MAX( 0._wp, pvt(ji,jj) ) * zrdt * e1v(ji,jj) / psm(ji,jj)
480            zalfq        =  zalf * zalf
481            zalf1        =  1.0 - zalf
482            zalf1q       =  zalf1 * zalf1
483            !
484            zfm (ji,jj)  =  zalf  * psm(ji,jj)
485            zf0 (ji,jj)  =  zalf  * ( ps0(ji,jj) + zalf1 * ( psy(ji,jj)  + (zalf1-zalf) * psyy(ji,jj)  ) ) 
486            zfy (ji,jj)  =  zalfq *( psy(ji,jj) + 3.0*zalf1*psyy(ji,jj) )
487            zfyy(ji,jj)  =  zalf  * zalfq * psyy(ji,jj)
488            zfx (ji,jj)  =  zalf  * ( psx(ji,jj) + zalf1 * psxy(ji,jj) )
489            zfxy(ji,jj)  =  zalfq * psxy(ji,jj)
490            zfxx(ji,jj)  =  zalf  * psxx(ji,jj)
491            !
492            !  Readjust moments remaining in the box.
493            psm (ji,jj)  =  psm (ji,jj) - zfm(ji,jj)
494            ps0 (ji,jj)  =  ps0 (ji,jj) - zf0(ji,jj)
495            psy (ji,jj)  =  zalf1q * ( psy(ji,jj) -3.0 * zalf * psyy(ji,jj) )
496            psyy(ji,jj)  =  zalf1 * zalf1q * psyy(ji,jj)
497            psx (ji,jj)  =  psx (ji,jj) - zfx(ji,jj)
498            psxx(ji,jj)  =  psxx(ji,jj) - zfxx(ji,jj)
499            psxy(ji,jj)  =  zalf1q * psxy(ji,jj)
500         END DO
501      END DO
502      !
503      DO jj = 1, jpjm1                   !  Flux from j+1 to j when v LT 0.
504         DO ji = 1, jpi
505            zalf          = ( MAX(0._wp, -pvt(ji,jj) ) * zrdt * e1v(ji,jj) ) / psm(ji,jj+1) 
506            zalg  (ji,jj) = zalf
507            zalfq         = zalf * zalf
508            zalf1         = 1.0 - zalf
509            zalg1 (ji,jj) = zalf1
510            zalf1q        = zalf1 * zalf1
511            zalg1q(ji,jj) = zalf1q
512            !
513            zfm   (ji,jj) = zfm (ji,jj) + zalf  *   psm (ji,jj+1)
514            zf0   (ji,jj) = zf0 (ji,jj) + zalf  * ( ps0 (ji,jj+1) - zalf1 * (psy(ji,jj+1) - (zalf1 - zalf ) * psyy(ji,jj+1) ) )
515            zfy   (ji,jj) = zfy (ji,jj) + zalfq * ( psy (ji,jj+1) - 3.0 * zalf1 * psyy(ji,jj+1) )
516            zfyy  (ji,jj) = zfyy(ji,jj) + zalf  *   psyy(ji,jj+1) * zalfq
517            zfx   (ji,jj) = zfx (ji,jj) + zalf  * ( psx (ji,jj+1) - zalf1 * psxy(ji,jj+1) )
518            zfxy  (ji,jj) = zfxy(ji,jj) + zalfq *   psxy(ji,jj+1)
519            zfxx  (ji,jj) = zfxx(ji,jj) + zalf  *   psxx(ji,jj+1)
520         END DO
521      END DO
522
523      !  Readjust moments remaining in the box.
524      DO jj = 2, jpj
525         DO ji = 1, jpi
526            zbt  =         zbet(ji,jj-1)
527            zbt1 = ( 1.0 - zbet(ji,jj-1) )
528            !
529            psm (ji,jj) = zbt * psm(ji,jj) + zbt1 * ( psm(ji,jj) - zfm(ji,jj-1) )
530            ps0 (ji,jj) = zbt * ps0(ji,jj) + zbt1 * ( ps0(ji,jj) - zf0(ji,jj-1) )
531            psy (ji,jj) = zalg1q(ji,jj-1) * ( psy(ji,jj) + 3.0 * zalg(ji,jj-1) * psyy(ji,jj) )
532            psyy(ji,jj) = zalg1 (ji,jj-1) * zalg1q(ji,jj-1) * psyy(ji,jj)
533            psx (ji,jj) = zbt * psx (ji,jj) + zbt1 * ( psx (ji,jj) - zfx (ji,jj-1) )
534            psxx(ji,jj) = zbt * psxx(ji,jj) + zbt1 * ( psxx(ji,jj) - zfxx(ji,jj-1) )
535            psxy(ji,jj) = zalg1q(ji,jj-1) * psxy(ji,jj)
536         END DO
537      END DO
538
539      !   Put the temporary moments into appropriate neighboring boxes.   
540      DO jj = 2, jpjm1                    !   Flux from j to j+1 IF v GT 0.
541         DO ji = 1, jpi
542            zbt  =         zbet(ji,jj-1)
543            zbt1 = ( 1.0 - zbet(ji,jj-1) )
544            psm(ji,jj)  = zbt * ( psm(ji,jj) + zfm(ji,jj-1) ) + zbt1 * psm(ji,jj) 
545            zalf        = zbt * zfm(ji,jj-1) / psm(ji,jj) 
546            zalf1       = 1.0 - zalf
547            ztemp       = zalf * ps0(ji,jj) - zalf1 * zf0(ji,jj-1)
548            !
549            ps0(ji,jj)  = zbt  * ( ps0(ji,jj) + zf0(ji,jj-1) ) + zbt1 * ps0(ji,jj)
550            psy(ji,jj)  = zbt  * ( zalf * zfy(ji,jj-1) + zalf1 * psy(ji,jj) + 3.0 * ztemp )   &
551               &                                               + zbt1 * psy(ji,jj) 
552            psyy(ji,jj) = zbt  * ( zalf * zalf * zfyy(ji,jj-1) + zalf1 * zalf1 * psyy(ji,jj)                             &
553               &                 + 5.0 * ( zalf * zalf1 * ( psy(ji,jj) - zfy(ji,jj-1) ) - ( zalf1 - zalf ) * ztemp ) )   & 
554               &                                               + zbt1 * psyy(ji,jj)
555            psxy(ji,jj) = zbt  * (  zalf * zfxy(ji,jj-1) + zalf1 * psxy(ji,jj)               &
556               &                  + 3.0 * (- zalf1 * zfx(ji,jj-1) + zalf * psx(ji,jj) )  )   &
557               &                                                + zbt1 * psxy(ji,jj)
558            psx (ji,jj) = zbt * ( psx (ji,jj) + zfx (ji,jj-1) ) + zbt1 * psx (ji,jj)
559            psxx(ji,jj) = zbt * ( psxx(ji,jj) + zfxx(ji,jj-1) ) + zbt1 * psxx(ji,jj)
560         END DO
561      END DO
562
563      DO jj = 2, jpjm1                   !  Flux from j+1 to j IF v LT 0.
564         DO ji = 1, jpi
565            zbt  =         zbet(ji,jj)
566            zbt1 = ( 1.0 - zbet(ji,jj) )
567            psm(ji,jj)  = zbt * psm(ji,jj) + zbt1 * ( psm(ji,jj) + zfm(ji,jj) )
568            zalf        = zbt1 * zfm(ji,jj) / psm(ji,jj)
569            zalf1       = 1.0 - zalf
570            ztemp       = - zalf * ps0 (ji,jj) + zalf1 * zf0(ji,jj)
571            ps0 (ji,jj) =   zbt  * ps0 (ji,jj) + zbt1  * ( ps0(ji,jj) + zf0(ji,jj) )
572            psy (ji,jj) =   zbt  * psy (ji,jj) + zbt1  * ( zalf * zfy(ji,jj) + zalf1 * psy(ji,jj) + 3.0 * ztemp )
573            psyy(ji,jj) =   zbt  * psyy(ji,jj) + zbt1 * (  zalf * zalf * zfyy(ji,jj) + zalf1 * zalf1 * psyy(ji,jj)   &
574               &                                         + 5.0 *( zalf *zalf1 *( -psy(ji,jj) + zfy(ji,jj) )          &
575               &                                         + ( zalf1 - zalf ) * ztemp )                                )
576            psxy(ji,jj) =   zbt  * psxy(ji,jj) + zbt1 * (  zalf * zfxy(ji,jj) + zalf1 * psxy(ji,jj)       &
577               &                                         + 3.0 * ( zalf1 * zfx(ji,jj) - zalf * psx(ji,jj) )  )
578            psx (ji,jj) =   zbt  * psx (ji,jj) + zbt1 * ( psx (ji,jj) + zfx (ji,jj) )
579            psxx(ji,jj) =   zbt  * psxx(ji,jj) + zbt1 * ( psxx(ji,jj) + zfxx(ji,jj) )
580         END DO
581      END DO
582
583      !-- Lateral boundary conditions
584      CALL lbc_lnk_multi( psm , 'T',  1.,  ps0 , 'T',  1.   &
585         &              , psx , 'T', -1.,  psy , 'T', -1.   &   ! caution gradient ==> the sign changes
586         &              , psxx, 'T',  1.,  psyy, 'T',  1.   &
587         &              , psxy, 'T',  1. )
588
589      IF(ln_ctl) THEN
590         CALL prt_ctl(tab2d_1=psm  , clinfo1=' adv_y: psm  :', tab2d_2=ps0 , clinfo2=' ps0  : ')
591         CALL prt_ctl(tab2d_1=psx  , clinfo1=' adv_y: psx  :', tab2d_2=psxx, clinfo2=' psxx : ')
592         CALL prt_ctl(tab2d_1=psy  , clinfo1=' adv_y: psy  :', tab2d_2=psyy, clinfo2=' psyy : ')
593         CALL prt_ctl(tab2d_1=psxy , clinfo1=' adv_y: psxy :')
594      ENDIF
595      !
596   END SUBROUTINE adv_y
597
598#else
599   !!----------------------------------------------------------------------
600   !!   Default option            Dummy module         NO LIM sea-ice model
601   !!----------------------------------------------------------------------
602#endif
603
604   !!======================================================================
605END MODULE iceadv_prather
Note: See TracBrowser for help on using the repository browser.