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

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

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