source: trunk/NEMO/LIM_SRC/limistate.F90 @ 508

Last change on this file since 508 was 508, checked in by opalod, 15 years ago

nemo_v1_update_071:RB: add iom for restart and reorganization of restart

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 12.3 KB
Line 
1MODULE limistate
2   !!======================================================================
3   !!                     ***  MODULE  limistate  ***
4   !!              Initialisation of diagnostics ice variables
5   !!======================================================================
6   !! History :   2.0  !  01-04  (C. Ethe, G. Madec)  Original code
7   !!                  !  04-04  (S. Theetten) initialization from a file
8   !!                  !  06-07  (S. Masson)  IOM to read the restart
9   !!--------------------------------------------------------------------
10#if defined key_ice_lim
11   !!----------------------------------------------------------------------
12   !!   'key_ice_lim' :                                   LIM sea-ice model
13   !!----------------------------------------------------------------------
14   !!----------------------------------------------------------------------
15   !!   lim_istate      :  Initialisation of diagnostics ice variables
16   !!   lim_istate_init :  initialization of ice state and namelist read
17   !!----------------------------------------------------------------------
18   USE phycst
19   USE ocfzpt
20   USE oce             ! dynamics and tracers variables      !!gm used???
21   USE dom_oce                                                     !!gm used???
22   USE par_ice         ! ice parameters
23   USE ice_oce         ! ice variables
24   USE dom_ice
25   USE ice             ! ???
26   USE lbclnk
27   USE ice
28   USE iom
29   USE in_out_manager
30
31   IMPLICIT NONE
32   PRIVATE
33
34   PUBLIC lim_istate      ! routine called by lim_init.F90
35
36   REAL(wp) ::           &  !!! ** init namelist (namiceini) **
37      ttest  = 2.0  ,    &  ! threshold water temperature for initial sea ice
38      hninn  = 0.5  ,    &  ! initial snow thickness in the north
39      hginn  = 3.0  ,    &  ! initial ice thickness in the north
40      alinn  = 0.05 ,    &  ! initial leads area in the north
41      hnins  = 0.1  ,    &  ! initial snow thickness in the south
42      hgins  = 1.0  ,    &  ! initial ice thickness in the south
43      alins  = 0.1          ! initial leads area in the south
44
45   REAL(wp)  ::          &  ! constant values
46      zzero   = 0.e0  ,  &
47      zone    = 1.e0
48   !!----------------------------------------------------------------------
49   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2006)
50   !! $Header$
51   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
52   !!----------------------------------------------------------------------
53
54CONTAINS
55
56   SUBROUTINE lim_istate
57      !!-------------------------------------------------------------------
58      !!                    ***  ROUTINE lim_istate  ***
59      !!
60      !! ** Purpose :   defined the sea-ice initial state
61      !!
62      !! ** Method  :   restart from a state defined in a binary file
63      !!                or from arbitrary sea-ice conditions
64      !!--------------------------------------------------------------------
65      INTEGER  ::   ji, jj, jk                ! dummy loop indices
66      REAL(wp) ::   zidto, zs0, ztf, zbin     ! temporary scalar
67      REAL(wp), DIMENSION(jpi,jpj) ::   ztn   ! workspace
68      !--------------------------------------------------------------------
69
70       CALL lim_istate_init     !  reading the initials parameters of the ice
71
72      !-- Initialisation of sst,sss,u,v do i=1,jpi
73      u_io(:,:)  = 0.e0       ! ice velocity in x direction
74      v_io(:,:)  = 0.e0       ! ice velocity in y direction
75
76      IF( ln_limini ) THEN    !
77       
78         ! Initialisation at tn if no ice or sst_ini if ice
79         ! Idem for salinity
80
81      !--- Criterion for presence (zidto=1.) or absence (zidto=0.) of ice
82         DO jj = 1 , jpj
83            DO ji = 1 , jpi
84               
85               zidto = MAX(zzero, - SIGN(1.,frld(ji,jj) - 1.))
86               
87               sst_io(ji,jj) = ( nfice - 1 ) * (zidto * sst_ini(ji,jj)  + &   ! use the ocean initial values
88                    &          (1.0 - zidto ) * ( tn(ji,jj,1) + rt0 ))        ! tricky trick *(nfice-1) !
89               sss_io(ji,jj) = ( nfice - 1 ) * (zidto * sss_ini(ji,jj) + &
90                    &          (1.0 - zidto ) *  sn(ji,jj,1) )
91
92               ! to avoid the the melting of ice, several layers (mixed layer) should be
93               ! set to sst_ini (sss_ini) if there is ice
94               ! example for one layer
95               ! tn(ji,jj,1) = zidto * ( sst_ini(ji,jj) - rt0 )  + (1.0 - zidto ) *  tn(ji,jj,1)
96               ! sn(ji,jj,1) = zidto * sss_ini(ji,jj)  + (1.0 - zidto ) *  sn(ji,jj,1)
97               ! tb(ji,jj,1) = tn(ji,jj,1)
98               ! sb(ji,jj,1) = sn(ji,jj,1)
99            END DO
100         END DO
101         
102         
103         !  tfu: Melting point of sea water
104         tfu(:,:)  = ztf   
105         
106         tfu(:,:)  = ABS ( rt0 - 0.0575       * sss_ini(:,:)                               &
107              &                    + 1.710523e-03 * sss_ini(:,:) * SQRT( sss_ini(:,:) )    &
108              &                    - 2.154996e-04 * sss_ini(:,:) * sss_ini(:,:) )
109      ELSE                     !
110
111         
112         ! Initialisation at tn or -2 if ice
113         DO jj = 1, jpj
114            DO ji = 1, jpi
115               zbin = MAX( 0., SIGN( 1., fzptn(ji,jj) - tn(ji,jj,1) ) )
116               ztn(ji,jj) = ( (1.-zbin) * tn(ji,jj,1) - 2. * zbin + rt0 ) * tmask(ji,jj,1)
117            END DO
118         END DO
119         
120         u_io  (:,:) = 0.e0
121         v_io  (:,:) = 0.e0
122         sst_io(:,:) = ( nfice - 1 ) * ( tn(:,:,1) + rt0 )   ! use the ocean initial values
123         sss_io(:,:) = ( nfice - 1 ) *   sn(:,:,1)           ! tricky trick *(nfice-1) !
124         
125         ! reference salinity 34psu
126         zs0 = 34.e0
127         ztf = ABS ( rt0 - 0.0575       * zs0                           &
128              &                    + 1.710523e-03 * zs0 * SQRT( zs0 )   &
129              &                    - 2.154996e-04 * zs0 *zs0          )
130         
131         !  tfu: Melting point of sea water
132         tfu(:,:)  = ztf   
133         
134         DO jj = 1, jpj
135            DO ji = 1, jpi
136               !--- Criterion for presence (zidto=1) or absence (zidto=0) of ice
137               zidto  = tms(ji,jj) * ( 1.0 - MAX(zzero, SIGN( zone, ztn(ji,jj) - tfu(ji,jj) - ttest) ) )
138               
139               IF( fcor(ji,jj) >= 0.e0 ) THEN     !--  Northern hemisphere.
140                  hicif(ji,jj)   = zidto * hginn
141                  frld(ji,jj)    = zidto * alinn + ( 1.0 - zidto ) * 1.0
142                  hsnif(ji,jj)   = zidto * hninn
143               ELSE                               !---  Southern hemisphere.
144                  hicif(ji,jj)   = zidto * hgins
145                  frld(ji,jj)    = zidto * alins + ( 1.0 - zidto ) * 1.0
146                  hsnif(ji,jj)   = zidto * hnins
147               ENDIF
148            END DO
149         END DO
150         
151         sist  (:,:)   = tfu(:,:)
152         tbif  (:,:,1) = tfu(:,:)
153         tbif  (:,:,2) = tfu(:,:)
154         tbif  (:,:,3) = tfu(:,:)
155     
156      ENDIF
157      fsbbq (:,:)   = 0.e0
158      qstoif(:,:)   = 0.e0
159      u_ice (:,:)   = 0.e0
160      v_ice (:,:)   = 0.e0
161# if defined key_coupled
162      albege(:,:)   = 0.8 * tms(:,:)
163# endif
164
165      !---  Moments for advection.             
166
167      sxice (:,:)  = 0.e0   ;   sxsn (:,:)  = 0.e0   ;   sxa  (:,:)  = 0.e0
168      syice (:,:)  = 0.e0   ;   sysn (:,:)  = 0.e0   ;   sya  (:,:)  = 0.e0
169      sxxice(:,:)  = 0.e0   ;   sxxsn(:,:)  = 0.e0   ;   sxxa (:,:)  = 0.e0
170      syyice(:,:)  = 0.e0   ;   syysn(:,:)  = 0.e0   ;   syya (:,:)  = 0.e0
171      sxyice(:,:)  = 0.e0   ;   sxysn(:,:)  = 0.e0   ;   sxya (:,:)  = 0.e0
172
173      sxc0  (:,:)  = 0.e0   ;   sxc1 (:,:)  = 0.e0   ;   sxc2 (:,:)  = 0.e0
174      syc0  (:,:)  = 0.e0   ;   syc1 (:,:)  = 0.e0   ;   syc2 (:,:)  = 0.e0
175      sxxc0 (:,:)  = 0.e0   ;   sxxc1(:,:)  = 0.e0   ;   sxxc2(:,:)  = 0.e0
176      syyc0 (:,:)  = 0.e0   ;   syyc1(:,:)  = 0.e0   ;   syyc2(:,:)  = 0.e0
177      sxyc0 (:,:)  = 0.e0   ;   sxyc1(:,:)  = 0.e0   ;   sxyc2(:,:)  = 0.e0
178
179      sxst  (:,:)  = 0.e0
180      syst  (:,:)  = 0.e0
181      sxxst (:,:)  = 0.e0
182      syyst (:,:)  = 0.e0
183      sxyst (:,:)  = 0.e0
184
185      !-- lateral boundary conditions
186      CALL lbc_lnk( hicif, 'T', 1. )
187      CALL lbc_lnk( frld , 'T', 1. )
188
189      ! C A U T I O N  frld = 1 over land and lbc_lnk put zero along
190      ! *************  closed boundaries herefore we force to one over land
191      frld(:,:) = tms(:,:) * frld(:,:) + ( 1. - tms(:,:) )   
192
193      CALL lbc_lnk( hsnif, 'T', 1. )
194      CALL lbc_lnk( sist , 'T', 1. )
195      DO jk = 1, jplayersp1
196         CALL lbc_lnk(tbif(:,:,jk), 'T', 1. )
197      END DO
198      CALL lbc_lnk( fsbbq  , 'T', 1. )
199      CALL lbc_lnk( qstoif , 'T', 1. )
200      CALL lbc_lnk( sss_io , 'T', 1. )
201      !
202   END SUBROUTINE lim_istate
203
204   
205   SUBROUTINE lim_istate_init
206      !!-------------------------------------------------------------------
207      !!                   ***  ROUTINE lim_istate_init  ***
208      !!       
209      !! ** Purpose :   Definition of initial state of the ice
210      !!
211      !! ** Method  :   Read the namiceini namelist and check the parameter
212      !!                values called at the first timestep (nit000)
213      !!                or
214      !!                Read 7 variables from a previous restart file
215      !!                sst, sst, hicif, hsnif, frld, ts & tbif
216      !!
217      !! ** input   :   Namelist namiceini
218      !!-------------------------------------------------------------------
219      INTEGER :: inum_ice
220
221      NAMELIST/namiceini/ ln_limini, ln_limdmp, ttest, hninn, hginn, alinn, &
222         &                hnins, hgins, alins
223      !!-------------------------------------------------------------------
224     
225      ! Read Namelist namiceini
226      REWIND ( numnam_ice )
227      READ   ( numnam_ice , namiceini )
228     
229      IF(.NOT. ln_limini) THEN
230         IF(lwp) THEN
231            WRITE(numout,*)
232            WRITE(numout,*) 'lim_istate_init : ice parameters inititialisation '
233            WRITE(numout,*) '~~~~~~~~~~~~~~~'
234            WRITE(numout,*) '         threshold water temp. for initial sea-ice    ttest      = ', ttest
235            WRITE(numout,*) '         initial snow thickness in the north          hninn      = ', hninn
236            WRITE(numout,*) '         initial ice thickness in the north           hginn      = ', hginn 
237            WRITE(numout,*) '         initial leads area in the north              alinn      = ', alinn           
238            WRITE(numout,*) '         initial snow thickness in the south          hnins      = ', hnins 
239            WRITE(numout,*) '         initial ice thickness in the south           hgins      = ', hgins
240            WRITE(numout,*) '         initial leads area in the south              alins      = ', alins
241         ENDIF
242      ENDIF
243
244      IF( ln_limini ) THEN                      ! Ice initialization using input file
245
246         CALL iom_open( 'Ice_initialization.nc', inum_ice )
247
248         IF( inum_ice > 0 ) THEN
249            IF(lwp) THEN
250               WRITE(numout,*) ' '
251               WRITE(numout,*) 'lim_istate_init : ice state initialization with : Ice_initialization.nc'
252               WRITE(numout,*) '~~~~~~~~~~~~~~~'
253               WRITE(numout,*) '         Ice state initialization using input file    ln_limini  = ', ln_limini
254               WRITE(numout,*) '         Ice damping                                  ln_limdmp  = ', ln_limdmp
255               WRITE(numout,*) ' '
256            ENDIF
257           
258            CALL iom_get( inum_ice, jpdom_data, 'sst'  , sst_ini(:,:) )       
259            CALL iom_get( inum_ice, jpdom_data, 'sss'  , sss_ini(:,:) )       
260            CALL iom_get( inum_ice, jpdom_data, 'hicif', hicif  (:,:) )     
261            CALL iom_get( inum_ice, jpdom_data, 'hsnif', hsnif  (:,:) )     
262            CALL iom_get( inum_ice, jpdom_data, 'frld' , frld   (:,:) )     
263            CALL iom_get( inum_ice, jpdom_data, 'ts'   , sist   (:,:) )
264            CALL iom_get( inum_ice, jpdom_unknown, 'tbif', tbif(:,:,:),   &
265                 &        kstart = (/ mig(1),mjg(1),1 /), kcount = (/ nlci,nlcj,jplayersp1 /) )
266
267            CALL iom_close( inum_ice)
268           
269         ENDIF
270      ENDIF
271      !
272   END SUBROUTINE lim_istate_init
273
274#else
275   !!----------------------------------------------------------------------
276   !!   Default option :         Empty module          NO LIM sea-ice model
277   !!----------------------------------------------------------------------
278CONTAINS
279   SUBROUTINE lim_istate          ! Empty routine
280   END SUBROUTINE lim_istate
281#endif
282
283   !!======================================================================
284END MODULE limistate
Note: See TracBrowser for help on using the repository browser.