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_2.F90 in branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2 – NEMO

source: branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limhdf_2.F90 @ 1855

Last change on this file since 1855 was 1855, checked in by gm, 14 years ago

ticket:#665 style change only, with the suppression of thd_ice_2 (merged in ice_2)
  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 6.7 KB
Line 
1MODULE limhdf_2
2   !!======================================================================
3   !!                    ***  MODULE limhdf_2   ***
4   !! LIM 2.0 ice model : horizontal diffusion of sea-ice quantities
5   !!======================================================================
6   !! History :  LIM  ! 2000-01 (UCL)  Original code
7   !!            1.0  ! 2001-05 (G. Madec, R. Hordoir) opa norm
8   !!            2.0  ! 2002-08 (C. Ethe)  F90, free form
9   !!----------------------------------------------------------------------
10#if defined key_lim2
11   !!----------------------------------------------------------------------
12   !!   'key_lim2'                                    LIM 2.0 sea-ice model
13   !!----------------------------------------------------------------------
14   !!   lim_hdf_2  : diffusion trend on sea-ice variable
15   !!----------------------------------------------------------------------
16   USE dom_oce          ! ocean domain
17   USE ice_2            ! LIM-2 variables
18   USE prtctl           ! Print control
19   USE lbclnk           !
20   USE lib_mpp          !
21   USE in_out_manager   ! I/O manager
22
23   IMPLICIT NONE
24   PRIVATE
25
26   PUBLIC   lim_hdf_2   ! called by lim_tra_2
27
28   LOGICAL  ::   linit = .TRUE.              ! indictor of initialisation
29   REAL(wp) ::   epsi04 = 1e-04              ! constant
30   REAL(wp), DIMENSION(jpi,jpj) ::   zfact   ! metric term
31
32   !! * Substitution
33#  include "vectopt_loop_substitute.h90"
34   !!----------------------------------------------------------------------
35   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)
36   !! $Id$
37   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
38   !!----------------------------------------------------------------------
39
40CONTAINS
41
42   SUBROUTINE lim_hdf_2( ptab )
43      !!-------------------------------------------------------------------
44      !!                  ***  ROUTINE lim_hdf_2  ***
45      !!
46      !! ** purpose :   Compute and add the diffusive trend on sea-ice variables
47      !!
48      !! ** method  :   Second order diffusive operator evaluated using a
49      !!              Cranck-Nicholson time Scheme.
50      !!
51      !! ** Action  :   update ptab with the diffusive contribution
52      !!-------------------------------------------------------------------
53      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj) ::   ptab   ! Field on which the diffusion is applied 
54      !!
55      INTEGER ::  ji, jj      ! dummy loop indices
56      INTEGER ::  its, iter   ! temporary integers
57      CHARACTER (len=55) :: charout
58      REAL(wp) ::   zalfa, zrlxint, zconv, zeps   ! temporary scalars
59      REAL(wp), DIMENSION(jpi,jpj) ::   zrlx, zflu, zflv, zdiv0, zdiv  ! temporary workspaces
60      REAL(wp), DIMENSION(jpi,jpj) ::   ztab0                ! ???
61      !!-------------------------------------------------------------------
62
63      ! Initialisation
64      ! ---------------   
65      ! Time integration parameters
66      zalfa = 0.5       ! =1.0/0.5/0.0 = implicit/Cranck-Nicholson/explicit
67      its   = 100       ! Maximum number of iteration
68      zeps  =  2. * epsi04
69
70      ! Arrays initialization
71      ztab0(:, : ) = ptab(:,:)
72      zdiv0(:, 1 ) = 0.e0
73      zdiv0(:,jpj) = 0.e0
74      IF( .NOT.lk_vopt_loop ) THEN
75         zflu (jpi,:) = 0.e0   
76         zflv (jpi,:) = 0.e0
77         zdiv0(1,  :) = 0.e0
78         zdiv0(jpi,:) = 0.e0
79      ENDIF
80
81      ! Metric coefficient (compute at the first call and saved in
82      IF( linit ) THEN
83         DO jj = 2, jpjm1 
84            DO ji = fs_2 , fs_jpim1   ! vector opt.
85               zfact(ji,jj) = ( e2u(ji,jj) + e2u(ji-1,jj) + e1v(ji,jj) + e1v(ji,jj-1) ) / ( e1t(ji,jj) * e2t(ji,jj) )
86            END DO
87         END DO
88         linit = .FALSE.
89      ENDIF
90
91
92      ! Sub-time step loop
93      zconv = 1.e0
94      iter  = 0
95
96      !                                                   !======================!
97      DO WHILE ( ( zconv > zeps ) .AND. (iter <= its) )   !  Sub-time step loop  !
98         !                                                !======================!
99         iter = iter + 1               ! incrementation of the sub-time step number
100
101         DO jj = 1, jpjm1              ! diffusive fluxes in U- and V- direction
102            DO ji = 1 , fs_jpim1   ! vector opt.
103               zflu(ji,jj) = pahu(ji,jj) * e2u(ji,jj) / e1u(ji,jj) * ( ptab(ji+1,jj) - ptab(ji,jj) )
104               zflv(ji,jj) = pahv(ji,jj) * e1v(ji,jj) / e2v(ji,jj) * ( ptab(ji,jj+1) - ptab(ji,jj) )
105            END DO
106         END DO
107         DO jj= 2, jpjm1               ! diffusive trend : divergence of the fluxes
108            DO ji = fs_2 , fs_jpim1   ! vector opt.
109               zdiv (ji,jj) = (  zflu(ji,jj) - zflu(ji-1,jj  )   &
110                  &            + zflv(ji,jj) - zflv(ji  ,jj-1)  ) / ( e1t (ji,jj) * e2t (ji,jj) )
111            END DO
112         END DO
113         !                             ! save the first evaluation of the diffusive trend in zdiv0
114         IF( iter == 1 )   zdiv0(:,:) = zdiv(:,:)       
115
116         DO jj = 2, jpjm1              ! XXXX iterative evaluation?????
117            DO ji = fs_2 , fs_jpim1   ! vector opt.
118               zrlxint = (   ztab0(ji,jj)    &
119                  &       +  rdt_ice * (           zalfa   * ( zdiv(ji,jj) + zfact(ji,jj) * ptab(ji,jj) )   &
120                  &                      + ( 1.0 - zalfa ) *   zdiv0(ji,jj) )  )                             & 
121                  &    / ( 1.0 + zalfa * rdt_ice * zfact(ji,jj) )
122               zrlx(ji,jj) = ptab(ji,jj) + om * ( zrlxint - ptab(ji,jj) )
123            END DO
124         END DO
125         CALL lbc_lnk( zrlx, 'T', 1. )
126
127         zconv = 0.e0                  ! convergence test
128         DO jj = 2, jpjm1
129            DO ji = 2, jpim1
130               zconv = MAX( zconv, ABS( zrlx(ji,jj) - ptab(ji,jj) )  )
131            END DO
132         END DO
133         IF( lk_mpp )   CALL mpp_max( zconv )   ! max over the global domain
134         !
135         ptab(:,:) = zrlx(:,:)         ! update value
136         !                                         !=============================!
137      END DO                                       !  end of sub-time step loop  !
138      !                                            !=============================!
139
140      IF(ln_ctl)   THEN
141         zrlx(:,:) = ptab(:,:) - ztab0(:,:)
142         WRITE(charout,FMT="(' lim_hdf  : zconv =',D23.16, ' iter =',I4,2X)") zconv, iter
143         CALL prt_ctl( tab2d_1=zrlx, clinfo1=charout )
144      ENDIF
145      !
146   END SUBROUTINE lim_hdf_2
147
148#else
149   !!----------------------------------------------------------------------
150   !!   Default option          Dummy module       NO LIM 2.0 sea-ice model
151   !!----------------------------------------------------------------------
152CONTAINS
153   SUBROUTINE lim_hdf_2       ! Empty routine
154   END SUBROUTINE lim_hdf_2
155#endif
156
157   !!======================================================================
158END MODULE limhdf_2
Note: See TracBrowser for help on using the repository browser.