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.
limvar.F90 in branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3 – NEMO

source: branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limvar.F90 @ 4332

Last change on this file since 4332 was 4332, checked in by clem, 10 years ago

update LIM3 to fix remaining bugs. Now working in global and regional config.
  • Property svn:keywords set to Id
File size: 22.9 KB
Line 
1MODULE limvar
2   !!======================================================================
3   !!                       ***  MODULE limvar ***
4   !!                 Different sets of ice model variables
5   !!                   how to switch from one to another
6   !!
7   !!                 There are three sets of variables
8   !!                 VGLO : global variables of the model
9   !!                        - v_i (jpi,jpj,jpl)
10   !!                        - v_s (jpi,jpj,jpl)
11   !!                        - a_i (jpi,jpj,jpl)
12   !!                        - t_s (jpi,jpj,jpl)
13   !!                        - e_i (jpi,jpj,nlay_i,jpl)
14   !!                        - smv_i(jpi,jpj,jpl)
15   !!                        - oa_i (jpi,jpj,jpl)
16   !!                 VEQV : equivalent variables sometimes used in the model
17   !!                        - ht_i(jpi,jpj,jpl)
18   !!                        - ht_s(jpi,jpj,jpl)
19   !!                        - t_i (jpi,jpj,nlay_i,jpl)
20   !!                        ...
21   !!                 VAGG : aggregate variables, averaged/summed over all
22   !!                        thickness categories
23   !!                        - vt_i(jpi,jpj)
24   !!                        - vt_s(jpi,jpj)
25   !!                        - at_i(jpi,jpj)
26   !!                        - et_s(jpi,jpj)  !total snow heat content
27   !!                        - et_i(jpi,jpj)  !total ice thermal content
28   !!                        - smt_i(jpi,jpj) !mean ice salinity
29   !!                        - ot_i(jpi,jpj)  !average ice age
30   !!======================================================================
31   !! History :   -   ! 2006-01 (M. Vancoppenolle) Original code
32   !!            4.0  ! 2011-02 (G. Madec) dynamical allocation
33   !!----------------------------------------------------------------------
34#if defined key_lim3
35   !!----------------------------------------------------------------------
36   !!   'key_lim3'                                      LIM3 sea-ice model
37   !!----------------------------------------------------------------------
38   !!   lim_var_agg       :
39   !!   lim_var_glo2eqv   :
40   !!   lim_var_eqv2glo   :
41   !!   lim_var_salprof   :
42   !!   lim_var_salprof1d :
43   !!   lim_var_bv        :
44   !!----------------------------------------------------------------------
45   USE par_oce        ! ocean parameters
46   USE phycst         ! physical constants (ocean directory)
47   USE sbc_oce        ! Surface boundary condition: ocean fields
48   USE ice            ! ice variables
49   USE par_ice        ! ice parameters
50   USE thd_ice        ! ice variables (thermodynamics)
51   USE dom_ice        ! ice domain
52   USE in_out_manager ! I/O manager
53   USE lib_mpp        ! MPP library
54   USE wrk_nemo       ! work arrays
55   USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 
56
57   IMPLICIT NONE
58   PRIVATE
59
60   PUBLIC   lim_var_agg          !
61   PUBLIC   lim_var_glo2eqv      !
62   PUBLIC   lim_var_eqv2glo      !
63   PUBLIC   lim_var_salprof      !
64   PUBLIC   lim_var_icetm        !
65   PUBLIC   lim_var_bv           !
66   PUBLIC   lim_var_salprof1d    !
67
68   REAL(wp) ::   epsi10 = 1.e-10_wp   !    -       -
69   REAL(wp) ::   zzero = 0.e0        !    -       -
70   REAL(wp) ::   zone  = 1.e0        !    -       -
71
72   !!----------------------------------------------------------------------
73   !! NEMO/LIM3 4.0 , UCL - NEMO Consortium (2011)
74   !! $Id$
75   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
76   !!----------------------------------------------------------------------
77CONTAINS
78
79   SUBROUTINE lim_var_agg( kn )
80      !!------------------------------------------------------------------
81      !!                ***  ROUTINE lim_var_agg  ***
82      !!
83      !! ** Purpose :   aggregates ice-thickness-category variables to all-ice variables
84      !!              i.e. it turns VGLO into VAGG
85      !! ** Method  :
86      !!
87      !! ** Arguments : n = 1, at_i vt_i only
88      !!                n = 2 everything
89      !!
90      !! note : you could add an argument when you need only at_i, vt_i
91      !!        and when you need everything
92      !!------------------------------------------------------------------
93      INTEGER, INTENT( in ) ::   kn     ! =1 at_i & vt only ; = what is needed
94      !
95      INTEGER  ::   ji, jj, jk, jl   ! dummy loop indices
96      REAL(wp) ::   zinda, zindb
97      !!------------------------------------------------------------------
98
99      !--------------------
100      ! Compute variables
101      !--------------------
102      vt_i (:,:) = 0._wp
103      vt_s (:,:) = 0._wp
104      at_i (:,:) = 0._wp
105      ato_i(:,:) = 1._wp
106      !
107      DO jl = 1, jpl
108         DO jj = 1, jpj
109            DO ji = 1, jpi
110               !
111               vt_i(ji,jj) = vt_i(ji,jj) + v_i(ji,jj,jl) ! ice volume
112               vt_s(ji,jj) = vt_s(ji,jj) + v_s(ji,jj,jl) ! snow volume
113               at_i(ji,jj) = at_i(ji,jj) + a_i(ji,jj,jl) ! ice concentration
114               !
115               zinda = MAX( zzero , SIGN( zone , at_i(ji,jj) - epsi10 ) ) 
116               icethi(ji,jj) = vt_i(ji,jj) / MAX( at_i(ji,jj) , epsi10 ) * zinda  ! ice thickness
117            END DO
118         END DO
119      END DO
120
121      DO jj = 1, jpj
122         DO ji = 1, jpi
123            ato_i(ji,jj) = MAX( 1._wp - at_i(ji,jj), 0._wp )   ! open water fraction
124         END DO
125      END DO
126
127      IF( kn > 1 ) THEN
128         et_s (:,:) = 0._wp
129         ot_i (:,:) = 0._wp
130         smt_i(:,:) = 0._wp
131         et_i (:,:) = 0._wp
132         !
133         DO jl = 1, jpl
134            DO jj = 1, jpj
135               DO ji = 1, jpi
136                  zinda = MAX( zzero , SIGN( zone , vt_i(ji,jj) - epsi10 ) ) 
137                  zindb = MAX( zzero , SIGN( zone , at_i(ji,jj) - epsi10 ) ) 
138                  et_s(ji,jj)  = et_s(ji,jj)  + e_s(ji,jj,1,jl)                                       ! snow heat content
139                  smt_i(ji,jj) = smt_i(ji,jj) + smv_i(ji,jj,jl) / MAX( vt_i(ji,jj) , epsi10 ) * zinda   ! ice salinity
140                  ot_i(ji,jj)  = ot_i(ji,jj)  + oa_i(ji,jj,jl)  / MAX( at_i(ji,jj) , epsi10 ) * zindb   ! ice age
141               END DO
142            END DO
143         END DO
144         !
145         DO jl = 1, jpl
146            DO jk = 1, nlay_i
147               et_i(:,:) = et_i(:,:) + e_i(:,:,jk,jl)       ! ice heat content
148            END DO
149         END DO
150         !
151      ENDIF
152      !
153   END SUBROUTINE lim_var_agg
154
155
156   SUBROUTINE lim_var_glo2eqv
157      !!------------------------------------------------------------------
158      !!                ***  ROUTINE lim_var_glo2eqv ***
159      !!
160      !! ** Purpose :   computes equivalent variables as function of global variables
161      !!              i.e. it turns VGLO into VEQV
162      !!------------------------------------------------------------------
163      INTEGER  ::   ji, jj, jk, jl   ! dummy loop indices
164      REAL(wp) ::   zq_i, zaaa, zbbb, zccc, zdiscrim     ! local scalars
165      REAL(wp) ::   ztmelts, zindb, zq_s, zfac1, zfac2   !   -      -
166      !!------------------------------------------------------------------
167
168      !-------------------------------------------------------
169      ! Ice thickness, snow thickness, ice salinity, ice age
170      !-------------------------------------------------------
171      DO jl = 1, jpl
172         DO jj = 1, jpj
173            DO ji = 1, jpi
174               zindb = 1._wp - MAX( 0._wp , SIGN( 1._wp,- a_i(ji,jj,jl) + epsi10 ) )   !0 if no ice and 1 if yes
175               ht_i(ji,jj,jl) = v_i (ji,jj,jl) / MAX( a_i(ji,jj,jl) , epsi10 ) * zindb
176               ht_s(ji,jj,jl) = v_s (ji,jj,jl) / MAX( a_i(ji,jj,jl) , epsi10 ) * zindb
177               o_i(ji,jj,jl)  = oa_i(ji,jj,jl) / MAX( a_i(ji,jj,jl) , epsi10 ) * zindb
178            END DO
179         END DO
180      END DO
181
182      IF(  num_sal == 2  )THEN
183         DO jl = 1, jpl
184            DO jj = 1, jpj
185               DO ji = 1, jpi
186                  zindb = 1._wp - MAX( 0._wp , SIGN( 1._wp,- a_i(ji,jj,jl) + epsi10 ) )   !0 if no ice and 1 if yes
187                  sm_i(ji,jj,jl) = smv_i(ji,jj,jl) / MAX( v_i(ji,jj,jl) , epsi10 ) * zindb
188               END DO
189            END DO
190         END DO
191      ENDIF
192
193      CALL lim_var_salprof      ! salinity profile
194
195      !-------------------
196      ! Ice temperatures
197      !-------------------
198!CDIR NOVERRCHK
199      DO jl = 1, jpl
200!CDIR NOVERRCHK
201         DO jk = 1, nlay_i
202!CDIR NOVERRCHK
203            DO jj = 1, jpj
204!CDIR NOVERRCHK
205               DO ji = 1, jpi
206                  !                                                              ! Energy of melting q(S,T) [J.m-3]
207                  zq_i    = e_i(ji,jj,jk,jl) / area(ji,jj) / MAX( v_i(ji,jj,jl) , epsi10 ) * REAL(nlay_i,wp) 
208                  zindb   = 1.0 - MAX( 0.0 , SIGN( 1.0 , - v_i(ji,jj,jl) + epsi10 ) )     ! zindb = 0 if no ice and 1 if yes
209                  zq_i    = zq_i * unit_fac * zindb                              !convert units
210                  ztmelts = -tmut * s_i(ji,jj,jk,jl) + rtt                       ! Ice layer melt temperature
211                  !
212                  zaaa       =  cpic                  ! Conversion q(S,T) -> T (second order equation)
213                  zbbb       =  ( rcp - cpic ) * ( ztmelts - rtt ) + zq_i / rhoic - lfus
214                  zccc       =  lfus * (ztmelts-rtt)
215                  zdiscrim   =  SQRT( MAX(zbbb*zbbb - 4._wp*zaaa*zccc , 0._wp) )
216                  t_i(ji,jj,jk,jl) = rtt + zindb *( - zbbb - zdiscrim ) / ( 2.0 *zaaa )
217                  t_i(ji,jj,jk,jl) = MIN( rtt, MAX( 173.15_wp, t_i(ji,jj,jk,jl) ) )       ! 100-rtt < t_i < rtt
218               END DO
219            END DO
220         END DO
221      END DO
222
223      !--------------------
224      ! Snow temperatures
225      !--------------------
226      zfac1 = 1._wp / ( rhosn * cpic )
227      zfac2 = lfus / cpic 
228      DO jl = 1, jpl
229         DO jk = 1, nlay_s
230            DO jj = 1, jpj
231               DO ji = 1, jpi
232                  !Energy of melting q(S,T) [J.m-3]
233                  zq_s  = e_s(ji,jj,jk,jl) / ( area(ji,jj) * MAX( v_s(ji,jj,jl) , epsi10 ) ) * REAL(nlay_s,wp)
234                  zindb = 1._wp - MAX( 0._wp , SIGN( 1._wp , - v_s(ji,jj,jl) + epsi10 ) )     ! zindb = 0 if no ice and 1 if yes
235                  zq_s  = zq_s * unit_fac * zindb                                    ! convert units
236                  !
237                  t_s(ji,jj,jk,jl) = rtt + zindb * ( - zfac1 * zq_s + zfac2 )
238                  t_s(ji,jj,jk,jl) = MIN( rtt, MAX( 173.15, t_s(ji,jj,jk,jl) ) )     ! 100-rtt < t_i < rtt
239               END DO
240            END DO
241         END DO
242      END DO
243
244      !-------------------
245      ! Mean temperature
246      !-------------------
247      tm_i(:,:) = 0._wp
248      DO jl = 1, jpl
249         DO jk = 1, nlay_i
250            DO jj = 1, jpj
251               DO ji = 1, jpi
252                  zindb = (  1._wp - MAX( 0._wp , SIGN( 1._wp , - vt_i(ji,jj) + epsi10 ) )  )
253                  tm_i(ji,jj) = tm_i(ji,jj) + zindb * t_i(ji,jj,jk,jl) * v_i(ji,jj,jl)   &
254                     &                      / (  REAL(nlay_i,wp) * MAX( vt_i(ji,jj) , epsi10 )  )
255               END DO
256            END DO
257         END DO
258      END DO
259      !
260   END SUBROUTINE lim_var_glo2eqv
261
262
263   SUBROUTINE lim_var_eqv2glo
264      !!------------------------------------------------------------------
265      !!                ***  ROUTINE lim_var_eqv2glo ***
266      !!
267      !! ** Purpose :   computes global variables as function of equivalent variables
268      !!                i.e. it turns VEQV into VGLO
269      !! ** Method  :
270      !!
271      !! ** History :  (01-2006) Martin Vancoppenolle, UCL-ASTR
272      !!------------------------------------------------------------------
273      !
274      v_i(:,:,:)   = ht_i(:,:,:) * a_i(:,:,:)
275      v_s(:,:,:)   = ht_s(:,:,:) * a_i(:,:,:)
276      smv_i(:,:,:) = sm_i(:,:,:) * v_i(:,:,:)
277      oa_i (:,:,:) = o_i (:,:,:) * a_i(:,:,:)
278      !
279   END SUBROUTINE lim_var_eqv2glo
280
281
282   SUBROUTINE lim_var_salprof
283      !!------------------------------------------------------------------
284      !!                ***  ROUTINE lim_var_salprof ***
285      !!
286      !! ** Purpose :   computes salinity profile in function of bulk salinity     
287      !!
288      !! ** Method  : If bulk salinity greater than s_i_1,
289      !!              the profile is assumed to be constant (S_inf)
290      !!              If bulk salinity lower than s_i_0,
291      !!              the profile is linear with 0 at the surface (S_zero)
292      !!              If it is between s_i_0 and s_i_1, it is a
293      !!              alpha-weighted linear combination of s_inf and s_zero
294      !!
295      !! ** References : Vancoppenolle et al., 2007 (in preparation)
296      !!------------------------------------------------------------------
297      INTEGER  ::   ji, jj, jk, jl   ! dummy loop index
298      REAL(wp) ::   dummy_fac0, dummy_fac1, dummy_fac, zsal      ! local scalar
299      REAL(wp) ::   zind0, zind01, zindbal, zargtemp , zs_zero   !   -      -
300      REAL(wp), POINTER, DIMENSION(:,:,:) ::   z_slope_s, zalpha   ! 3D pointer
301      !!------------------------------------------------------------------
302
303      CALL wrk_alloc( jpi, jpj, jpl, z_slope_s, zalpha )
304
305      !---------------------------------------
306      ! Vertically constant, constant in time
307      !---------------------------------------
308      IF(  num_sal == 1  )   s_i(:,:,:,:) = bulk_sal
309
310      !-----------------------------------
311      ! Salinity profile, varying in time
312      !-----------------------------------
313      IF(  num_sal == 2  ) THEN
314         !
315         DO jk = 1, nlay_i
316            s_i(:,:,jk,:)  = sm_i(:,:,:)
317         END DO
318         !
319         DO jl = 1, jpl                               ! Slope of the linear profile
320            DO jj = 1, jpj
321               DO ji = 1, jpi
322                  z_slope_s(ji,jj,jl) = 2._wp * sm_i(ji,jj,jl) / MAX( 0.01 , ht_i(ji,jj,jl) )
323               END DO
324            END DO
325         END DO
326         !
327         dummy_fac0 = 1._wp / ( s_i_0 - s_i_1 )       ! Weighting factor between zs_zero and zs_inf
328         dummy_fac1 = s_i_1 / ( s_i_1 - s_i_0 )
329         !
330         zalpha(:,:,:) = 0._wp
331         DO jl = 1, jpl
332            DO jj = 1, jpj
333               DO ji = 1, jpi
334                  ! zind0 = 1 if sm_i le s_i_0 and 0 otherwise
335                  zind0  = MAX( 0._wp   , SIGN( 1._wp  , s_i_0 - sm_i(ji,jj,jl) ) ) 
336                  ! zind01 = 1 if sm_i is between s_i_0 and s_i_1 and 0 othws
337                  zind01 = ( 1._wp - zind0 ) * MAX( 0._wp   , SIGN( 1._wp  , s_i_1 - sm_i(ji,jj,jl) ) ) 
338                  ! If 2.sm_i GE sss_m then zindbal = 1
339                  ! this is to force a constant salinity profile in the Baltic Sea
340                  zindbal = MAX( 0._wp , SIGN( 1._wp , 2._wp * sm_i(ji,jj,jl) - sss_m(ji,jj) ) )
341                  zalpha(ji,jj,jl) = zind0  + zind01 * ( sm_i(ji,jj,jl) * dummy_fac0 + dummy_fac1 )
342                  zalpha(ji,jj,jl) = zalpha(ji,jj,jl) * ( 1._wp - zindbal )
343               END DO
344            END DO
345         END DO
346
347         dummy_fac = 1._wp / REAL( nlay_i )                   ! Computation of the profile
348         DO jl = 1, jpl
349            DO jk = 1, nlay_i
350               DO jj = 1, jpj
351                  DO ji = 1, jpi
352                     !                                      ! linear profile with 0 at the surface
353                     zs_zero = z_slope_s(ji,jj,jl) * ( REAL(jk,wp) - 0.5_wp ) * ht_i(ji,jj,jl) * dummy_fac
354                     !                                      ! weighting the profile
355                     s_i(ji,jj,jk,jl) = zalpha(ji,jj,jl) * zs_zero + ( 1._wp - zalpha(ji,jj,jl) ) * sm_i(ji,jj,jl)
356                  END DO ! ji
357               END DO ! jj
358            END DO ! jk
359         END DO ! jl
360         !
361      ENDIF ! num_sal
362
363      !-------------------------------------------------------
364      ! Vertically varying salinity profile, constant in time
365      !-------------------------------------------------------
366
367      IF(  num_sal == 3  ) THEN      ! Schwarzacher (1959) multiyear salinity profile (mean = 2.30)
368         !
369         sm_i(:,:,:) = 2.30_wp
370         !
371         DO jl = 1, jpl
372!CDIR NOVERRCHK
373            DO jk = 1, nlay_i
374               zargtemp  = ( REAL(jk,wp) - 0.5_wp ) / REAL(nlay_i,wp)
375               zsal =  1.6_wp * (  1._wp - COS( rpi * zargtemp**(0.407_wp/(0.573_wp+zargtemp)) )  )
376               s_i(:,:,jk,jl) =  zsal
377            END DO
378         END DO
379         !
380      ENDIF ! num_sal
381      !
382      CALL wrk_dealloc( jpi, jpj, jpl, z_slope_s, zalpha )
383      !
384   END SUBROUTINE lim_var_salprof
385
386
387   SUBROUTINE lim_var_icetm
388      !!------------------------------------------------------------------
389      !!                ***  ROUTINE lim_var_icetm ***
390      !!
391      !! ** Purpose :   computes mean sea ice temperature
392      !!------------------------------------------------------------------
393      INTEGER  ::   ji, jj, jk, jl   ! dummy loop indices
394      REAL(wp) ::   zindb   !   -      -
395      !!------------------------------------------------------------------
396
397      ! Mean sea ice temperature
398      tm_i(:,:) = 0._wp
399      DO jl = 1, jpl
400         DO jk = 1, nlay_i
401            DO jj = 1, jpj
402               DO ji = 1, jpi
403                  zindb = (  1._wp - MAX( 0._wp , SIGN( 1._wp , - vt_i(ji,jj) + epsi10 ) )  )
404                  tm_i(ji,jj) = tm_i(ji,jj) + zindb * t_i(ji,jj,jk,jl) * v_i(ji,jj,jl)   &
405                     &                      / (  REAL(nlay_i,wp) * MAX( vt_i(ji,jj) , epsi10 )  )
406               END DO
407            END DO
408         END DO
409      END DO
410
411   END SUBROUTINE lim_var_icetm
412
413
414   SUBROUTINE lim_var_bv
415      !!------------------------------------------------------------------
416      !!                ***  ROUTINE lim_var_bv ***
417      !!
418      !! ** Purpose :   computes mean brine volume (%) in sea ice
419      !!
420      !! ** Method  : e = - 0.054 * S (ppt) / T (C)
421      !!
422      !! References : Vancoppenolle et al., JGR, 2007
423      !!------------------------------------------------------------------
424      INTEGER  ::   ji, jj, jk, jl   ! dummy loop indices
425      REAL(wp) ::   zbvi, zinda, zindb      ! local scalars
426      !!------------------------------------------------------------------
427      !
428      bv_i(:,:) = 0._wp
429      DO jl = 1, jpl
430         DO jk = 1, nlay_i
431            DO jj = 1, jpj
432               DO ji = 1, jpi
433                  zinda = (  1._wp - MAX( 0._wp , SIGN( 1._wp , (t_i(ji,jj,jk,jl) - rtt) + epsi10 ) )  )
434                  zindb = (  1._wp - MAX( 0._wp , SIGN( 1._wp , - vt_i(ji,jj) + epsi10 ) )  )
435                  zbvi  = - zinda * tmut * s_i(ji,jj,jk,jl) / MIN( t_i(ji,jj,jk,jl) - rtt, - epsi10 )   &
436                     &                   * v_i(ji,jj,jl)    / REAL(nlay_i,wp)
437                  bv_i(ji,jj) = bv_i(ji,jj) + zindb * zbvi  / MAX( vt_i(ji,jj) , epsi10 )
438               END DO
439            END DO
440         END DO
441      END DO
442      !
443   END SUBROUTINE lim_var_bv
444
445
446   SUBROUTINE lim_var_salprof1d( kideb, kiut )
447      !!-------------------------------------------------------------------
448      !!                  ***  ROUTINE lim_thd_salprof1d  ***
449      !!
450      !! ** Purpose :   1d computation of the sea ice salinity profile
451      !!                Works with 1d vectors and is used by thermodynamic modules
452      !!-------------------------------------------------------------------
453      INTEGER, INTENT(in) ::   kideb, kiut   ! thickness category index
454      !
455      INTEGER  ::   ji, jk    ! dummy loop indices
456      INTEGER  ::   ii, ij  ! local integers
457      REAL(wp) ::   dummy_fac0, dummy_fac1, dummy_fac2, zargtemp, zsal   ! local scalars
458      REAL(wp) ::   zalpha, zind0, zind01, zindbal, zs_zero              !   -      -
459      !
460      REAL(wp), POINTER, DIMENSION(:) ::   z_slope_s
461      !!---------------------------------------------------------------------
462
463      CALL wrk_alloc( jpij, z_slope_s )
464
465      !---------------------------------------
466      ! Vertically constant, constant in time
467      !---------------------------------------
468      IF( num_sal == 1 )   s_i_b(:,:) = bulk_sal
469
470      !------------------------------------------------------
471      ! Vertically varying salinity profile, varying in time
472      !------------------------------------------------------
473
474      IF(  num_sal == 2  ) THEN
475         !
476         DO ji = kideb, kiut          ! Slope of the linear profile zs_zero
477            z_slope_s(ji) = 2._wp * sm_i_b(ji) / MAX( 0.01 , ht_i_b(ji) )
478         END DO
479
480         ! Weighting factor between zs_zero and zs_inf
481         !---------------------------------------------
482         dummy_fac0 = 1._wp / ( s_i_0 - s_i_1 )
483         dummy_fac1 = s_i_1 / ( s_i_1 - s_i_0 )
484         dummy_fac2 = 1._wp / REAL(nlay_i,wp)
485
486!CDIR NOVERRCHK
487         DO jk = 1, nlay_i
488!CDIR NOVERRCHK
489            DO ji = kideb, kiut
490               ii =  MOD( npb(ji) - 1 , jpi ) + 1
491               ij =     ( npb(ji) - 1 ) / jpi + 1
492               ! zind0 = 1 if sm_i le s_i_0 and 0 otherwise
493               zind0  = MAX( 0._wp , SIGN( 1._wp  , s_i_0 - sm_i_b(ji) ) ) 
494               ! zind01 = 1 if sm_i is between s_i_0 and s_i_1 and 0 othws
495               zind01 = ( 1._wp - zind0 ) * MAX( 0._wp , SIGN( 1._wp , s_i_1 - sm_i_b(ji) ) ) 
496               ! if 2.sm_i GE sss_m then zindbal = 1
497               ! this is to force a constant salinity profile in the Baltic Sea
498               zindbal = MAX( 0._wp , SIGN( 1._wp , 2._wp * sm_i_b(ji) - sss_m(ii,ij) ) )
499               !
500               zalpha = (  zind0 + zind01 * ( sm_i_b(ji) * dummy_fac0 + dummy_fac1 )  ) * ( 1.0 - zindbal )
501               !
502               zs_zero = z_slope_s(ji) * ( REAL(jk,wp) - 0.5_wp ) * ht_i_b(ji) * dummy_fac2
503               ! weighting the profile
504               s_i_b(ji,jk) = zalpha * zs_zero + ( 1._wp - zalpha ) * sm_i_b(ji)
505            END DO ! ji
506         END DO ! jk
507
508      ENDIF ! num_sal
509
510      !-------------------------------------------------------
511      ! Vertically varying salinity profile, constant in time
512      !-------------------------------------------------------
513
514      IF( num_sal == 3 ) THEN      ! Schwarzacher (1959) multiyear salinity profile (mean = 2.30)
515         !
516         sm_i_b(:) = 2.30_wp
517         !
518!CDIR NOVERRCHK
519         DO jk = 1, nlay_i
520            zargtemp  = ( REAL(jk,wp) - 0.5_wp ) / REAL(nlay_i,wp)
521            zsal =  1.6_wp * (  1._wp - COS( rpi * zargtemp**(0.407_wp/(0.573_wp+zargtemp)) )  )
522            DO ji = kideb, kiut
523               s_i_b(ji,jk) = zsal
524            END DO
525         END DO
526         !
527      ENDIF
528      !
529      CALL wrk_dealloc( jpij, z_slope_s )
530      !
531   END SUBROUTINE lim_var_salprof1d
532
533#else
534   !!----------------------------------------------------------------------
535   !!   Default option         Dummy module          NO  LIM3 sea-ice model
536   !!----------------------------------------------------------------------
537CONTAINS
538   SUBROUTINE lim_var_agg          ! Empty routines
539   END SUBROUTINE lim_var_agg
540   SUBROUTINE lim_var_glo2eqv      ! Empty routines
541   END SUBROUTINE lim_var_glo2eqv
542   SUBROUTINE lim_var_eqv2glo      ! Empty routines
543   END SUBROUTINE lim_var_eqv2glo
544   SUBROUTINE lim_var_salprof      ! Empty routines
545   END SUBROUTINE lim_var_salprof
546   SUBROUTINE lim_var_bv           ! Emtpy routines
547   END SUBROUTINE lim_var_bv
548   SUBROUTINE lim_var_salprof1d    ! Emtpy routines
549   END SUBROUTINE lim_var_salprof1d
550#endif
551
552   !!======================================================================
553END MODULE limvar
Note: See TracBrowser for help on using the repository browser.