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.
limhdf.F90 in trunk/NEMO/LIM_SRC – NEMO

source: trunk/NEMO/LIM_SRC/limhdf.F90 @ 106

Last change on this file since 106 was 106, checked in by opalod, 20 years ago

CT : UPDATE067 : Add control indices nictl, njctl used in SUM function output to compare mono versus multi procs runs
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 6.8 KB
Line 
1MODULE limhdf
2   !!======================================================================
3   !!                    ***  MODULE limhdf   ***
4   !! LIM ice model : horizontal diffusion of sea-ice quantities
5   !!======================================================================
6#if defined key_ice_lim
7   !!----------------------------------------------------------------------
8   !!   'key_ice_lim'                                     LIM sea-ice model
9   !!----------------------------------------------------------------------
10   !!   lim_hdf  : diffusion trend on sea-ice variable
11   !!----------------------------------------------------------------------
12   !! * Modules used
13   USE dom_oce
14   USE ice_oce         ! ice variables
15   USE in_out_manager
16   USE ice
17!  USE limdyn
18   USE lbclnk
19   USE lib_mpp
20
21   IMPLICIT NONE
22   PRIVATE
23
24   !! * Routine accessibility
25   PUBLIC lim_hdf    ! called by lim_tra
26
27   !! * Module variables
28   LOGICAL  ::   linit = .TRUE.              ! ???
29   REAL(wp) ::   epsi04 = 1e-04              ! constant
30   REAL(wp), DIMENSION(jpi,jpj) ::   zfact   ! ???
31
32   !! * Substitution
33#  include "vectopt_loop_substitute.h90"
34   !!----------------------------------------------------------------------
35   !!   LIM 2.0 , UCL-LODYC-IPSL  (2003)
36   !!----------------------------------------------------------------------
37
38CONTAINS
39
40   SUBROUTINE lim_hdf( ptab )
41      !!-------------------------------------------------------------------
42      !!                  ***  ROUTINE lim_hdf  ***
43      !!
44      !! ** purpose :   Compute and add the diffusive trend on sea-ice
45      !!      variables
46      !!
47      !! ** method  :   Second order diffusive operator evaluated using a
48      !!      Cranck-Nicholson time Scheme.
49      !!
50      !! ** Action  :    update ptab with the diffusive contribution
51      !!
52      !! History :
53      !!        !  00-01 (LIM) Original code
54      !!        !  01-05 (G. Madec, R. Hordoir) opa norm
55      !!        !  02-08 (C. Ethe)  F90, free form
56      !!-------------------------------------------------------------------
57      ! * Arguments
58      REAL(wp), DIMENSION(jpi,jpj), INTENT( inout ) ::   &
59         ptab                 ! Field on which the diffusion is applied 
60      REAL(wp), DIMENSION(jpi,jpj) ::   &
61         ptab0                ! ???
62
63      ! * Local variables
64      INTEGER ::  ji, jj      ! dummy loop indices
65      INTEGER ::  &
66         its, iter            ! temporary integers
67      REAL(wp) ::  &
68         zalfa, zrlxint, zconv, zeps   ! temporary scalars
69      REAL(wp), DIMENSION(jpi,jpj) ::  & 
70         zrlx, zflu, zflv, &  ! temporary workspaces
71         zdiv0, zdiv          !    "           "
72      !!-------------------------------------------------------------------
73
74      ! Initialisation
75      ! ---------------   
76      ! Time integration parameters
77      zalfa = 0.5       ! =1.0/0.5/0.0 = implicit/Cranck-Nicholson/explicit
78      its   = 100       ! Maximum number of iteration
79      zeps  =  2. * epsi04
80
81      ! Arrays initialization
82      ptab0 (:, : ) = ptab(:,:)
83!bug  zflu (:,jpj) = 0.e0
84!bug  zflv (:,jpj) = 0.e0
85      zdiv0(:, 1 ) = 0.e0
86      zdiv0(:,jpj) = 0.e0
87      IF( .NOT.lk_vopt_loop ) THEN
88         zflu (jpi,:) = 0.e0   
89         zflv (jpi,:) = 0.e0
90         zdiv0(1,  :) = 0.e0
91         zdiv0(jpi,:) = 0.e0
92      ENDIF
93
94      ! Metric coefficient (compute at the first call and saved in
95      IF( linit ) THEN
96         DO jj = 2, jpjm1 
97            DO ji = fs_2 , fs_jpim1   ! vector opt.
98               zfact(ji,jj) = ( e2u(ji,jj) + e2u(ji-1,jj  ) + e1v(ji,jj) + e1v(ji,jj-1) ) &
99                  &          / ( e1t(ji,jj) * e2t(ji,jj) )
100            END DO
101         END DO
102         linit = .FALSE.
103      ENDIF
104
105
106      ! Sub-time step loop
107      zconv = 1.e0
108      iter  = 0
109
110      !                                                   !===================
111      DO WHILE ( ( zconv > zeps ) .AND. (iter <= its) )   ! Sub-time step loop
112         !                                                !===================
113         ! incrementation of the sub-time step number
114         iter = iter + 1
115
116         ! diffusive fluxes in U- and V- direction
117         DO jj = 1, jpjm1
118            DO ji = 1 , fs_jpim1   ! vector opt.
119               zflu(ji,jj) = pahu(ji,jj) * e2u(ji,jj) / e1u(ji,jj) * ( ptab(ji+1,jj) - ptab(ji,jj) )
120               zflv(ji,jj) = pahv(ji,jj) * e1v(ji,jj) / e2v(ji,jj) * ( ptab(ji,jj+1) - ptab(ji,jj) )
121            END DO
122         END DO
123
124         ! diffusive trend : divergence of the fluxes
125         DO jj= 2, jpjm1
126            DO ji = fs_2 , fs_jpim1   ! vector opt.
127               zdiv (ji,jj) = (  zflu(ji,jj) - zflu(ji-1,jj  )   &
128                  &            + zflv(ji,jj) - zflv(ji  ,jj-1)  ) / ( e1t (ji,jj) * e2t (ji,jj) )
129            END DO
130         END DO
131
132         ! save the first evaluation of the diffusive trend in zdiv0
133         IF( iter == 1 )   zdiv0(:,:) = zdiv(:,:)       
134
135         ! XXXX iterative evaluation?????
136         DO jj = 2, jpjm1
137            DO ji = fs_2 , fs_jpim1   ! vector opt.
138               zrlxint = (   ptab0(ji,jj)    &
139                  &       +  rdt_ice * (           zalfa   * ( zdiv(ji,jj) + zfact(ji,jj) * ptab(ji,jj) )   &
140                  &                      + ( 1.0 - zalfa ) *   zdiv0(ji,jj) )  )                             & 
141                  &    / ( 1.0 + zalfa * rdt_ice * zfact(ji,jj) )
142               zrlx(ji,jj) = ptab(ji,jj) + om * ( zrlxint - ptab(ji,jj) )
143            END DO
144         END DO
145
146         ! convergence test
147         zconv = 0.e0
148         DO jj = 2, jpjm1
149            DO ji = 2, jpim1
150               zconv = MAX( zconv, ABS( zrlx(ji,jj) - ptab(ji,jj) )  )
151            END DO
152         END DO
153         IF( lk_mpp )   CALL mpp_max( zconv )   ! max over the global domain
154
155         DO jj = 2, jpjm1
156            DO ji = 2 , jpim1
157               ptab(ji,jj) = zrlx(ji,jj)
158            END DO
159         END DO
160
161         ! lateral boundary condition on ptab
162         CALL lbc_lnk( ptab, 'T', 1. )
163         !                                         !==========================
164      END DO                                       ! end of sub-time step loop
165      !                                            !==========================
166
167      ptab(:,:) = ptab(:,:)
168
169      IF(l_ctl)   WRITE(numout,*) ' lim_hdf  : ', SUM( ptab(2:nictl,2:njctl)-ptab0(2:nictl,2:njctl) ), &
170      &                                        ' zconv= ', zconv, ' iter= ', iter
171
172   END SUBROUTINE lim_hdf
173
174#else
175   !!----------------------------------------------------------------------
176   !!   Default option          Dummy module           NO LIM sea-ice model
177   !!----------------------------------------------------------------------
178CONTAINS
179   SUBROUTINE lim_hdf         ! Empty routine
180   END SUBROUTINE lim_hdf
181#endif
182
183   !!======================================================================
184END MODULE limhdf
Note: See TracBrowser for help on using the repository browser.