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 @ 247

Last change on this file since 247 was 247, checked in by opalod, 19 years ago

CL : Add CVS Header and CeCILL licence information
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 6.9 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-LOCEAN-IPSL (2005)
36   !! $Header$
37   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt
38   !!----------------------------------------------------------------------
39
40CONTAINS
41
42   SUBROUTINE lim_hdf( ptab )
43      !!-------------------------------------------------------------------
44      !!                  ***  ROUTINE lim_hdf  ***
45      !!
46      !! ** purpose :   Compute and add the diffusive trend on sea-ice
47      !!      variables
48      !!
49      !! ** method  :   Second order diffusive operator evaluated using a
50      !!      Cranck-Nicholson time Scheme.
51      !!
52      !! ** Action  :    update ptab with the diffusive contribution
53      !!
54      !! History :
55      !!        !  00-01 (LIM) Original code
56      !!        !  01-05 (G. Madec, R. Hordoir) opa norm
57      !!        !  02-08 (C. Ethe)  F90, free form
58      !!-------------------------------------------------------------------
59      ! * Arguments
60      REAL(wp), DIMENSION(jpi,jpj), INTENT( inout ) ::   &
61         ptab                 ! Field on which the diffusion is applied 
62      REAL(wp), DIMENSION(jpi,jpj) ::   &
63         ptab0                ! ???
64
65      ! * Local variables
66      INTEGER ::  ji, jj      ! dummy loop indices
67      INTEGER ::  &
68         its, iter            ! temporary integers
69      REAL(wp) ::  &
70         zalfa, zrlxint, zconv, zeps   ! temporary scalars
71      REAL(wp), DIMENSION(jpi,jpj) ::  & 
72         zrlx, zflu, zflv, &  ! temporary workspaces
73         zdiv0, zdiv          !    "           "
74      !!-------------------------------------------------------------------
75
76      ! Initialisation
77      ! ---------------   
78      ! Time integration parameters
79      zalfa = 0.5       ! =1.0/0.5/0.0 = implicit/Cranck-Nicholson/explicit
80      its   = 100       ! Maximum number of iteration
81      zeps  =  2. * epsi04
82
83      ! Arrays initialization
84      ptab0 (:, : ) = ptab(:,:)
85!bug  zflu (:,jpj) = 0.e0
86!bug  zflv (:,jpj) = 0.e0
87      zdiv0(:, 1 ) = 0.e0
88      zdiv0(:,jpj) = 0.e0
89      IF( .NOT.lk_vopt_loop ) THEN
90         zflu (jpi,:) = 0.e0   
91         zflv (jpi,:) = 0.e0
92         zdiv0(1,  :) = 0.e0
93         zdiv0(jpi,:) = 0.e0
94      ENDIF
95
96      ! Metric coefficient (compute at the first call and saved in
97      IF( linit ) THEN
98         DO jj = 2, jpjm1 
99            DO ji = fs_2 , fs_jpim1   ! vector opt.
100               zfact(ji,jj) = ( e2u(ji,jj) + e2u(ji-1,jj  ) + e1v(ji,jj) + e1v(ji,jj-1) ) &
101                  &          / ( e1t(ji,jj) * e2t(ji,jj) )
102            END DO
103         END DO
104         linit = .FALSE.
105      ENDIF
106
107
108      ! Sub-time step loop
109      zconv = 1.e0
110      iter  = 0
111
112      !                                                   !===================
113      DO WHILE ( ( zconv > zeps ) .AND. (iter <= its) )   ! Sub-time step loop
114         !                                                !===================
115         ! incrementation of the sub-time step number
116         iter = iter + 1
117
118         ! diffusive fluxes in U- and V- direction
119         DO jj = 1, jpjm1
120            DO ji = 1 , fs_jpim1   ! vector opt.
121               zflu(ji,jj) = pahu(ji,jj) * e2u(ji,jj) / e1u(ji,jj) * ( ptab(ji+1,jj) - ptab(ji,jj) )
122               zflv(ji,jj) = pahv(ji,jj) * e1v(ji,jj) / e2v(ji,jj) * ( ptab(ji,jj+1) - ptab(ji,jj) )
123            END DO
124         END DO
125
126         ! diffusive trend : divergence of the fluxes
127         DO jj= 2, jpjm1
128            DO ji = fs_2 , fs_jpim1   ! vector opt.
129               zdiv (ji,jj) = (  zflu(ji,jj) - zflu(ji-1,jj  )   &
130                  &            + zflv(ji,jj) - zflv(ji  ,jj-1)  ) / ( e1t (ji,jj) * e2t (ji,jj) )
131            END DO
132         END DO
133
134         ! save the first evaluation of the diffusive trend in zdiv0
135         IF( iter == 1 )   zdiv0(:,:) = zdiv(:,:)       
136
137         ! XXXX iterative evaluation?????
138         DO jj = 2, jpjm1
139            DO ji = fs_2 , fs_jpim1   ! vector opt.
140               zrlxint = (   ptab0(ji,jj)    &
141                  &       +  rdt_ice * (           zalfa   * ( zdiv(ji,jj) + zfact(ji,jj) * ptab(ji,jj) )   &
142                  &                      + ( 1.0 - zalfa ) *   zdiv0(ji,jj) )  )                             & 
143                  &    / ( 1.0 + zalfa * rdt_ice * zfact(ji,jj) )
144               zrlx(ji,jj) = ptab(ji,jj) + om * ( zrlxint - ptab(ji,jj) )
145            END DO
146         END DO
147
148         ! convergence test
149         zconv = 0.e0
150         DO jj = 2, jpjm1
151            DO ji = 2, jpim1
152               zconv = MAX( zconv, ABS( zrlx(ji,jj) - ptab(ji,jj) )  )
153            END DO
154         END DO
155         IF( lk_mpp )   CALL mpp_max( zconv )   ! max over the global domain
156
157         DO jj = 2, jpjm1
158            DO ji = 2 , jpim1
159               ptab(ji,jj) = zrlx(ji,jj)
160            END DO
161         END DO
162
163         ! lateral boundary condition on ptab
164         CALL lbc_lnk( ptab, 'T', 1. )
165         !                                         !==========================
166      END DO                                       ! end of sub-time step loop
167      !                                            !==========================
168
169      ptab(:,:) = ptab(:,:)
170
171      IF(l_ctl)   WRITE(numout,*) ' lim_hdf  : ', SUM( ptab(2:nictl,2:njctl)-ptab0(2:nictl,2:njctl) ), &
172      &                                        ' zconv= ', zconv, ' iter= ', iter
173
174   END SUBROUTINE lim_hdf
175
176#else
177   !!----------------------------------------------------------------------
178   !!   Default option          Dummy module           NO LIM sea-ice model
179   !!----------------------------------------------------------------------
180CONTAINS
181   SUBROUTINE lim_hdf         ! Empty routine
182   END SUBROUTINE lim_hdf
183#endif
184
185   !!======================================================================
186END MODULE limhdf
Note: See TracBrowser for help on using the repository browser.