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 trunk/NEMO/LIM_SRC_3 – NEMO

source: trunk/NEMO/LIM_SRC_3/limvar.F90 @ 867

Last change on this file since 867 was 834, checked in by ctlod, 16 years ago

Clean comments and useless lines, see ticket:#72

File size: 26.6 KB
Line 
1MODULE limvar
2#if defined key_lim3
3   !!----------------------------------------------------------------------
4   !!   'key_lim3'                                      LIM3 sea-ice model
5   !!----------------------------------------------------------------------
6   !!======================================================================
7   !!                       ***  MODULE limvar ***
8   !!                 Different sets of ice model variables
9   !!                   how to switch from one to another
10   !!
11   !!                 There are three sets of variables
12   !!                 VGLO : global variables of the model
13   !!                        - v_i (jpi,jpj,jpl)
14   !!                        - v_s (jpi,jpj,jpl)
15   !!                        - a_i (jpi,jpj,jpl)
16   !!                        - t_s (jpi,jpj,jpl)
17   !!                        - e_i (jpi,jpj,nlay_i,jpl)
18   !!                        - smv_i(jpi,jpj,jpl)
19   !!                        - oa_i (jpi,jpj,jpl)
20   !!                 VEQV : equivalent variables sometimes used in the model
21   !!                        - ht_i(jpi,jpj,jpl)
22   !!                        - ht_s(jpi,jpj,jpl)
23   !!                        - t_i (jpi,jpj,nlay_i,jpl)
24   !!                        ...
25   !!                 VAGG : aggregate variables, averaged/summed over all
26   !!                        thickness categories
27   !!                        - vt_i(jpi,jpj)
28   !!                        - vt_s(jpi,jpj)
29   !!                        - at_i(jpi,jpj)
30   !!                        - et_s(jpi,jpj)  !total snow heat content
31   !!                        - et_i(jpi,jpj)  !total ice thermal content
32   !!                        - smt_i(jpi,jpj) !mean ice salinity
33   !!                        - ot_i(jpi,jpj)  !average ice age
34   !!======================================================================
35
36   !!----------------------------------------------------------------------
37   !! * Modules used
38   USE dom_ice
39   USE par_oce          ! ocean parameters
40   USE phycst           ! physical constants (ocean directory)
41   USE ice_oce          ! ice variables
42   USE thd_ice
43   USE in_out_manager
44   USE ice
45   USE par_ice
46 
47   IMPLICIT NONE
48   PRIVATE
49
50   !! * Routine accessibility
51   PUBLIC lim_var_agg
52   PUBLIC lim_var_glo2eqv
53   PUBLIC lim_var_eqv2glo
54   PUBLIC lim_var_salprof
55   PUBLIC lim_var_bv
56   PUBLIC lim_var_salprof1d
57
58   !! * Module variables
59   REAL(wp)  ::           &  ! constant values
60      epsi20 = 1e-20   ,  &
61      epsi13 = 1e-13   ,  &
62      zzero  = 0.e0    ,  &
63      zone   = 1.e0
64
65   !!----------------------------------------------------------------------
66   !!   LIM3.0,  UCL-ASTR-LOCEAN-IPSL (2008)
67   !!   (c) UCL-ASTR and Martin Vancoppenolle
68   !!----------------------------------------------------------------------
69
70CONTAINS
71
72   SUBROUTINE lim_var_agg(n)
73        !!------------------------------------------------------------------
74        !!                ***  ROUTINE lim_var_agg  ***
75        !! ** Purpose :
76        !!        This routine aggregates ice-thickness-category variables to 
77        !!                                all-ice variables
78        !!        i.e. it turns VGLO into VAGG
79        !! ** Method  :
80        !!
81        !! ** Arguments :
82        !!           kideb , kiut : Starting and ending points on which the
83        !!                         the computation is applied
84        !!
85        !! ** Inputs / Ouputs : (global commons)
86        !! ** Arguments : n = 1, at_i vt_i only
87        !!                n = 2 everything
88        !!
89        !! ** External :
90        !!
91        !! ** References :
92        !!
93        !! ** History :
94        !!           (01-2006) Martin Vancoppenolle, UCL-ASTR
95        !!
96        !! note : you could add an argument when you need only at_i, vt_i
97        !!        and when you need everything
98        !!------------------------------------------------------------------
99        !! * Arguments
100
101        !! * Local variables
102        INTEGER ::   ji,       &   ! spatial dummy loop index
103                     jj,       &   ! spatial dummy loop index
104                     jk,       &   ! vertical layering dummy loop index
105                     jl            ! ice category dummy loop index
106 
107        REAL ::      zeps, epsi16, zinda, epsi06
108 
109        INTEGER, INTENT( in ) ::   n     ! describes what is needed
110       
111!!-- End of declarations
112!!----------------------------------------------------------------------------------------------
113       zeps = 1.0e-13
114       epsi16 = 1.0e-16
115       epsi06 = 1.0e-6
116
117       !------------------
118       ! Zero everything
119       !------------------
120
121       vt_i(:,:)  = 0.0
122       vt_s(:,:)  = 0.0
123       at_i(:,:)  = 0.0 
124       ato_i(:,:) = 1.0 
125
126       IF ( n .GT. 1 ) THEN
127          et_s(:,:)  = 0.0
128          ot_i(:,:)  = 0.0
129          smt_i(:,:) = 0.0
130          et_i(:,:)  = 0.0
131       ENDIF
132     
133       !--------------------
134       ! Compute variables
135       !--------------------
136
137       DO jl = 1, jpl
138          DO jj = 1, jpj
139             DO ji = 1, jpi
140
141                vt_i(ji,jj) = vt_i(ji,jj) + v_i(ji,jj,jl) ! ice volume
142                vt_s(ji,jj) = vt_s(ji,jj) + v_s(ji,jj,jl) ! snow volume
143                at_i(ji,jj) = at_i(ji,jj) + a_i(ji,jj,jl) ! ice concentration
144
145                zinda = MAX( zzero , SIGN( zone , at_i(ji,jj) - 0.10 ) ) 
146                icethi(ji,jj) = vt_i(ji,jj) / MAX(at_i(ji,jj),epsi16)*zinda 
147                   ! ice thickness
148             END DO
149          END DO
150       END DO
151
152       DO jj = 1, jpj
153          DO ji = 1, jpi
154             ato_i(ji,jj) = MAX(1.0 - at_i(ji,jj), 0.0)   ! open water fraction
155          END DO
156       END DO
157
158       IF ( n .GT. 1 ) THEN
159
160       DO jl = 1, jpl
161          DO jj = 1, jpj
162             DO ji = 1, jpi
163                et_s(ji,jj)  = et_s(ji,jj)  +     &       ! snow heat content
164                               e_s(ji,jj,1,jl)
165                zinda = MAX( zzero , SIGN( zone , vt_i(ji,jj) - 0.10 ) ) 
166                smt_i(ji,jj) = smt_i(ji,jj) +     &       ! ice salinity
167                               smv_i(ji,jj,jl) / MAX( vt_i(ji,jj) , zeps ) * &
168                               zinda
169                zinda = MAX( zzero , SIGN( zone , at_i(ji,jj) - 0.10 ) ) 
170                ot_i(ji,jj)  = ot_i(ji,jj)  +     &       ! ice age
171                               oa_i(ji,jj,jl)  / MAX( at_i(ji,jj) , zeps ) * &
172                               zinda
173             END DO
174          END DO
175       END DO
176               
177       DO jl = 1, jpl
178          DO jk = 1, nlay_i
179             DO jj = 1, jpj
180                DO ji = 1, jpi
181                   et_i(ji,jj) = et_i(ji,jj) + e_i(ji,jj,jk,jl) ! ice heat
182                                                                ! content
183                END DO
184             END DO
185          END DO
186       END DO
187
188       ENDIF ! n .GT. 1
189
190    END SUBROUTINE lim_var_agg
191
192!==============================================================================
193
194    SUBROUTINE lim_var_glo2eqv
195        !!------------------------------------------------------------------
196        !!                ***  ROUTINE lim_var_glo2eqv ***'
197        !! ** Purpose :
198        !!        This routine computes equivalent variables as function of   
199        !!                              global variables
200        !!        i.e. it turns VGLO into VEQV
201        !! ** Method  :
202        !!
203        !! ** Arguments :
204        !!           kideb , kiut : Starting and ending points on which the
205        !!                         the computation is applied
206        !!
207        !! ** Inputs / Ouputs :
208        !!
209        !! ** External :
210        !!
211        !! ** References :
212        !!
213        !! ** History :
214        !!           (01-2006) Martin Vancoppenolle, UCL-ASTR
215        !!
216        !!------------------------------------------------------------------
217
218        !! * Local variables
219        INTEGER ::   ji,       &   ! spatial dummy loop index
220                     jj,       &   ! spatial dummy loop index
221                     jk,       &   ! vertical layering dummy loop index
222                     jl            ! ice category dummy loop index
223
224        REAL :: zq_i, zaaa, zbbb, zccc, zdiscrim, &
225                ztmelts, zindb, zq_s, zfac1, zfac2
226
227        REAL :: zeps, epsi06
228
229        zeps    = 1.0e-10
230        epsi06  = 1.0e-06
231
232!!-- End of declarations
233!!------------------------------------------------------------------------------
234
235      !-------------------------------------------------------
236      ! Ice thickness, snow thickness, ice salinity, ice age
237      !-------------------------------------------------------
238      DO jl = 1, jpl
239         DO jj = 1, jpj
240            DO ji = 1, jpi
241               zindb          = 1.0-MAX(0.0,SIGN(1.0,- a_i(ji,jj,jl))) !0 if no ice and 1 if yes
242               ht_i(ji,jj,jl) = v_i(ji,jj,jl)   / MAX( a_i(ji,jj,jl) , zeps ) * zindb
243               ht_s(ji,jj,jl) = v_s(ji,jj,jl)   / MAX( a_i(ji,jj,jl) , zeps ) * zindb
244               o_i(ji,jj,jl)  = oa_i(ji,jj,jl)  / MAX( a_i(ji,jj,jl) , zeps ) * zindb
245            END DO
246         END DO
247      END DO
248
249      IF ( ( num_sal .EQ. 2 ) .OR. ( num_sal .EQ. 4 ) )THEN
250
251      DO jl = 1, jpl
252         DO jj = 1, jpj
253            DO ji = 1, jpi
254               zindb          = 1.0-MAX(0.0,SIGN(1.0,-a_i(ji,jj,jl))) !0 if no ice and 1 if yes
255               sm_i(ji,jj,jl) = smv_i(ji,jj,jl) / MAX(v_i(ji,jj,jl),zeps) * zindb
256            END DO
257         END DO
258      END DO
259
260      ENDIF
261
262      ! salinity profile
263      CALL lim_var_salprof
264
265      !-------------------
266      ! Ice temperatures
267      !-------------------
268      DO jl = 1, jpl
269        DO jk = 1, nlay_i
270          DO jj = 1, jpj
271            DO ji = 1, jpi
272              !Energy of melting q(S,T) [J.m-3]
273              zq_i       = e_i(ji,jj,jk,jl) / area(ji,jj) / &
274                           MAX( v_i(ji,jj,jl) , epsi06 ) * nlay_i 
275              ! zindb = 0 if no ice and 1 if yes
276              zindb      = 1.0 - MAX( 0.0 , SIGN( 1.0 , - v_i(ji,jj,jl) ) )
277              !convert units ! very important that this line is here
278              zq_i       = zq_i * unit_fac * zindb
279              !Ice layer melt temperature
280              ztmelts    =  -tmut*s_i(ji,jj,jk,jl) + rtt
281              !Conversion q(S,T) -> T (second order equation)
282              zaaa       =  cpic
283              zbbb       =  ( rcp - cpic ) * ( ztmelts - rtt ) + &
284                            zq_i / rhoic - lfus
285              zccc       =  lfus * (ztmelts-rtt)
286              zdiscrim   =  SQRT( MAX(zbbb*zbbb - 4.0*zaaa*zccc,0.0) )
287              t_i(ji,jj,jk,jl) = rtt + zindb *( - zbbb - zdiscrim ) / & 
288                                 ( 2.0 *zaaa )
289              t_i(ji,jj,jk,jl) = MIN( rtt, MAX(173.15, t_i(ji,jj,jk,jl) ) )
290            END DO
291          END DO
292        END DO
293      END DO
294
295      !--------------------
296      ! Snow temperatures
297      !--------------------
298      zfac1 = 1. / ( rhosn * cpic )
299      zfac2 = lfus / cpic 
300      DO jl = 1, jpl
301        DO jk = 1, nlay_s
302          DO jj = 1, jpj
303            DO ji = 1, jpi
304              !Energy of melting q(S,T) [J.m-3]
305              zq_s       = e_s(ji,jj,jk,jl) / area(ji,jj) / &
306                           MAX( v_s(ji,jj,jl) , epsi06 ) * nlay_s 
307              ! zindb = 0 if no ice and 1 if yes
308              zindb      = 1.0 - MAX( 0.0 , SIGN( 1.0 , - v_s(ji,jj,jl) ) )
309              !convert units ! very important that this line is here
310              zq_s       = zq_s * unit_fac * zindb
311              t_s(ji,jj,jk,jl) = rtt + zindb * ( - zfac1 * zq_s + zfac2 )
312              t_s(ji,jj,jk,jl) = MIN( rtt, MAX(173.15, t_s(ji,jj,jk,jl) ) )
313
314            END DO
315          END DO
316        END DO
317      END DO
318
319      !-------------------
320      ! Mean temperature
321      !-------------------
322      tm_i(:,:) = 0.0
323      DO jl = 1, jpl
324         DO jk = 1, nlay_i
325            DO jj = 1, jpj
326               DO ji = 1, jpi
327                  zindb          = 1.0-MAX(0.0,SIGN(1.0,-a_i(ji,jj,jl)))
328                  zindb          = zindb*1.0-MAX(0.0,SIGN(1.0,-v_i(ji,jj,jl)))
329                  tm_i(ji,jj) = tm_i(ji,jj) + t_i(ji,jj,jk,jl)*v_i(ji,jj,jl) / &
330                                REAL(nlay_i) / MAX( vt_i(ji,jj) , zeps )
331               END DO
332            END DO
333         END DO
334      END DO
335
336   END SUBROUTINE lim_var_glo2eqv
337
338!===============================================================================
339
340   SUBROUTINE lim_var_eqv2glo
341        !!------------------------------------------------------------------
342        !!                ***  ROUTINE lim_var_eqv2glo ***'
343        !! ** Purpose :
344        !!        This routine computes global     variables as function of   
345        !!                              equivalent variables
346        !!        i.e. it turns VEQV into VGLO
347        !! ** Method  :
348        !!
349        !! ** Arguments :
350        !!
351        !! ** Inputs / Ouputs : (global commons)
352        !!
353        !! ** External :
354        !!
355        !! ** References :
356        !!
357        !! ** History :
358        !!           (01-2006) Martin Vancoppenolle, UCL-ASTR
359        !!                     Take it easy man
360        !!                     Life is just a simple game, between
361        !!                     ups / and downs \ :@)
362        !!
363        !!------------------------------------------------------------------
364
365       v_i(:,:,:)   = ht_i(:,:,:) * a_i(:,:,:)
366       v_s(:,:,:)   = ht_s(:,:,:) * a_i(:,:,:)
367       smv_i(:,:,:) = sm_i(:,:,:) * v_i(:,:,:)
368       oa_i (:,:,:) = o_i (:,:,:) * a_i(:,:,:)
369
370    END SUBROUTINE lim_var_eqv2glo
371
372!===============================================================================
373
374    SUBROUTINE lim_var_salprof
375        !!------------------------------------------------------------------
376        !!                ***  ROUTINE lim_var_salprof ***'
377        !! ** Purpose :
378        !!        This routine computes salinity profile in function of
379        !!        bulk salinity     
380        !!
381        !! ** Method  : If bulk salinity greater than s_i_1,
382        !!              the profile is assumed to be constant (S_inf)
383        !!              If bulk salinity lower than s_i_0,
384        !!              the profile is linear with 0 at the surface (S_zero)
385        !!              If it is between s_i_0 and s_i_1, it is a
386        !!              alpha-weighted linear combination of s_inf and s_zero
387        !!
388        !! ** References : Vancoppenolle et al., 2007 (in preparation)
389        !!
390        !! ** History :
391        !!           (08-2006) Martin Vancoppenolle, UCL-ASTR
392        !!                     Take it easy man
393        !!                     Life is just a simple game, between ups
394        !!                     / and downs \ :@)
395        !!
396        !!------------------------------------------------------------------
397        !! * Arguments
398
399        !! * Local variables
400        INTEGER ::             &
401           ji            ,     & !: spatial dummy loop index
402           jj            ,     & !: spatial dummy loop index
403           jk            ,     & !: vertical layering dummy loop index
404           jl                    !: ice category dummy loop index
405
406        REAL(wp) ::            &
407           dummy_fac0    ,     & !: dummy factor used in computations
408           dummy_fac1    ,     & !: dummy factor used in computations
409           dummy_fac     ,     & !: dummy factor used in computations
410           zind0         ,     & !: switch, = 1 if sm_i lt s_i_0
411           zind01        ,     & !: switch, = 1 if sm_i between s_i_0 and s_i_1
412           zindbal       ,     & !: switch, = 1, if 2*sm_i gt sss_io
413           zargtemp              !: dummy factor
414
415        REAL(wp), DIMENSION(nlay_i) ::      &
416           zs_zero               !: linear salinity profile for salinities under
417                                 !: s_i_0
418
419        REAL(wp), DIMENSION(jpi,jpj,jpl) :: &
420           z_slope_s     ,     & !: slope of the salinity profile
421           zalpha                !: weight factor for s between s_i_0 and s_i_1
422
423!!-- End of declarations
424!!------------------------------------------------------------------------------
425
426      !---------------------------------------
427      ! Vertically constant, constant in time
428      !---------------------------------------
429
430      IF ( num_sal .EQ. 1 ) THEN
431
432         s_i(:,:,:,:) = bulk_sal
433
434      ENDIF
435
436      !-----------------------------------
437      ! Salinity profile, varying in time
438      !-----------------------------------
439
440      IF ( ( num_sal .EQ. 2 ) .OR. ( num_sal .EQ. 4 ) )THEN
441
442         DO jk = 1, nlay_i
443            s_i(:,:,jk,:)  = sm_i(:,:,:)
444         END DO ! jk
445
446         ! Slope of the linear profile zs_zero
447         !-------------------------------------
448         DO jl = 1, jpl
449            DO jj = 1, jpj
450               DO ji = 1, jpi
451                  z_slope_s(ji,jj,jl) = 2.0 * sm_i(ji,jj,jl) / MAX( 0.01      &
452                                      , ht_i(ji,jj,jl) )
453               END DO ! ji
454            END DO ! jj
455         END DO ! jl
456
457         ! Weighting factor between zs_zero and zs_inf
458         !---------------------------------------------
459         dummy_fac0 = 1. / ( ( s_i_0 - s_i_1 ) )
460         dummy_fac1 = s_i_1 / ( s_i_1 - s_i_0 )
461
462         zalpha(:,:,:) = 0.0
463
464         DO jl = 1, jpl
465            DO jj = 1, jpj
466               DO ji = 1, jpi
467                  ! zind0 = 1 if sm_i le s_i_0 and 0 otherwise
468                  zind0  = MAX( 0.0   , SIGN( 1.0  , s_i_0 - sm_i(ji,jj,jl) ) ) 
469                  ! zind01 = 1 if sm_i is between s_i_0 and s_i_1 and 0 othws
470                  zind01 = ( 1.0 - zind0 ) *                                  &
471                           MAX( 0.0   , SIGN( 1.0  , s_i_1 - sm_i(ji,jj,jl) ) ) 
472                  ! If 2.sm_i GE sss_io then zindbal = 1
473                  zindbal = MAX( 0.0 , SIGN( 1.0 , 2. * sm_i(ji,jj,jl) -      &
474                  sss_io(ji,jj) ) )
475                  zalpha(ji,jj,jl) = zind0  * 1.0                             &
476                                   + zind01 * ( sm_i(ji,jj,jl) * dummy_fac0 + &
477                                                dummy_fac1 )
478                  zalpha(ji,jj,jl) = zalpha(ji,jj,jl) * ( 1.0 - zindbal )
479               END DO
480            END DO
481         END DO
482
483         ! Computation of the profile
484         !----------------------------
485         dummy_fac = 1. / nlay_i
486
487         DO jl = 1, jpl
488            DO jk = 1, nlay_i
489               DO jj = 1, jpj
490                  DO ji = 1, jpi
491                     ! linear profile with 0 at the surface
492                     zs_zero(jk)      = z_slope_s(ji,jj,jl) * ( jk - 1./2. ) * &
493                                        ht_i(ji,jj,jl) * dummy_fac
494                     ! weighting the profile
495                     s_i(ji,jj,jk,jl) = zalpha(ji,jj,jl) * zs_zero(jk) +       &
496                                     ( 1.0 - zalpha(ji,jj,jl) ) * sm_i(ji,jj,jl)
497                  END DO ! ji
498               END DO ! jj
499            END DO ! jk
500         END DO ! jl
501
502      ENDIF ! num_sal
503
504      !-------------------------------------------------------
505      ! Vertically varying salinity profile, constant in time
506      !-------------------------------------------------------
507      ! Schwarzacher (1959) multiyear salinity profile (mean = 2.30)
508     
509      IF ( num_sal .EQ. 3 ) THEN
510
511         sm_i(:,:,:) = 2.30
512
513         DO jl = 1, jpl
514            DO jk = 1, nlay_i
515               DO jj = 1, jpj
516                  DO ji = 1, jpi
517                     zargtemp  = ( jk - 0.5 ) / nlay_i
518                     s_i(ji,jj,jk,jl) =  1.6 - 1.6 * COS( 3.14169265 * &
519                                         ( zargtemp**(0.407/           &
520                                         ( 0.573 + zargtemp ) ) ) )
521                  END DO ! ji
522               END DO ! jj
523            END DO ! jk
524         END DO ! jl
525
526      ENDIF ! num_sal
527
528   END SUBROUTINE lim_var_salprof
529
530!===============================================================================
531
532   SUBROUTINE lim_var_bv
533        !!------------------------------------------------------------------
534        !!                ***  ROUTINE lim_var_bv ***'
535        !! ** Purpose :
536        !!        This routine computes mean brine volume (%) in sea ice
537        !!
538        !! ** Method  : e = - 0.054 * S (ppt) / T (C)
539        !!
540        !! ** Arguments :
541        !!
542        !! ** Inputs / Ouputs : (global commons)
543        !!
544        !! ** External :
545        !!
546        !! ** References : Vancoppenolle et al., JGR, 2007
547        !!
548        !! ** History :
549        !!           (08-2006) Martin Vancoppenolle, UCL-ASTR
550        !!
551        !!------------------------------------------------------------------
552        !! * Arguments
553
554        !! * Local variables
555        INTEGER ::   ji,       &   ! spatial dummy loop index
556                     jj,       &   ! spatial dummy loop index
557                     jk,       &   ! vertical layering dummy loop index
558                     jl            ! ice category dummy loop index
559
560        REAL :: zbvi,          &   ! brine volume for a single ice category
561                zeps,          &   ! very small value
562                zindb              ! is there ice or not
563
564!!-- End of declarations
565!!------------------------------------------------------------------------------
566
567       zeps = 1.0e-13
568       bv_i(:,:) = 0.0
569       DO jl = 1, jpl
570          DO jk = 1, nlay_i
571             DO jj = 1, jpj
572                DO ji = 1, jpi
573                   zindb          = 1.0-MAX(0.0,SIGN(1.0,-a_i(ji,jj,jl))) !0 if no ice and 1 if yes
574                   zbvi = - zindb * tmut *s_i(ji,jj,jk,jl) /             & 
575                            MIN( t_i(ji,jj,jk,jl) - 273.15 , zeps )         &
576                            * v_i(ji,jj,jl) / REAL(nlay_i)
577                   bv_i(ji,jj) = bv_i(ji,jj) + zbvi  &
578                            / MAX( vt_i(ji,jj) , zeps )
579                END DO
580             END DO
581          END DO
582       END DO
583
584   END SUBROUTINE lim_var_bv 
585
586!===============================================================================
587
588   SUBROUTINE lim_var_salprof1d(kideb,kiut)
589      !!-------------------------------------------------------------------
590      !!                  ***  ROUTINE lim_thd_salprof1d  ***
591      !!
592      !! ** Purpose :   1d computation of the sea ice salinity profile
593      !!                Works with 1d vectors and is used by thermodynamic
594      !!                modules
595      !!
596      !! history :
597      !!   3.0  !  May  2007 M. Vancoppenolle  Original code
598      !!-------------------------------------------------------------------
599      INTEGER, INTENT(in) :: &
600         kideb, kiut             ! thickness category index
601
602      INTEGER ::             &
603         ji, jk,             &   ! geographic and layer index
604         zji, zjj
605
606      REAL(wp) ::            &
607         dummy_fac0,         &   ! dummy factors
608         dummy_fac1,         &
609         dummy_fac2,         &
610         zalpha    ,         &   ! weighting factor
611         zind0     ,         &   ! switches as in limvar
612         zind01    ,         &   ! switch
613         zindbal   ,         &   ! switch if in freshwater area
614         zargtemp
615   
616      REAL(wp), DIMENSION(jpij) ::            &
617         z_slope_s
618
619      REAL(wp), DIMENSION(jpij,jkmax) ::      &
620         zs_zero
621      !!-------------------------------------------------------------------
622         
623      !---------------------------------------
624      ! Vertically constant, constant in time
625      !---------------------------------------
626
627      IF ( num_sal .EQ. 1 ) THEN
628
629         s_i_b(:,:) = bulk_sal
630
631      ENDIF
632
633      !------------------------------------------------------
634      ! Vertically varying salinity profile, varying in time
635      !------------------------------------------------------
636
637      IF ( ( num_sal .EQ. 2 ) .OR. ( num_sal .EQ. 4 ) ) THEN
638
639         ! Slope of the linear profile zs_zero
640         !-------------------------------------
641         DO ji = kideb, kiut 
642               z_slope_s(ji) = 2.0 * sm_i_b(ji) / MAX( 0.01      &
643                                      , ht_i_b(ji) )
644         END DO ! ji
645
646         ! Weighting factor between zs_zero and zs_inf
647         !---------------------------------------------
648         dummy_fac0 = 1. / ( ( s_i_0 - s_i_1 ) )
649         dummy_fac1 = s_i_1 / ( s_i_1 - s_i_0 )
650         dummy_fac2 = 1. / nlay_i
651
652         DO jk = 1, nlay_i
653            DO ji = kideb, kiut
654               zji    =  MOD( npb(ji) - 1, jpi ) + 1
655               zjj    =  ( npb(ji) - 1 ) / jpi + 1
656               zalpha = 0.0
657               ! zind0 = 1 if sm_i le s_i_0 and 0 otherwise
658               zind0  = MAX( 0.0   , SIGN( 1.0  , s_i_0 - sm_i_b(ji) ) ) 
659               ! zind01 = 1 if sm_i is between s_i_0 and s_i_1 and 0 othws
660               zind01 = ( 1.0 - zind0 ) *                                  &
661                        MAX( 0.0   , SIGN( 1.0  , s_i_1 - sm_i_b(ji) ) ) 
662               ! if 2.sm_i GE sss_io then zindbal = 1
663               zindbal = MAX( 0.0 , SIGN( 1.0 , 2. * sm_i_b(ji) -      &
664               sss_io(zji,zjj) ) )
665
666               zalpha = zind0  * 1.0                                       &
667                      + zind01 * ( sm_i_b(ji) * dummy_fac0 +           &
668                                                dummy_fac1 )
669               zalpha = zalpha * ( 1.0 - zindbal )
670
671               zs_zero(ji,jk) = z_slope_s(ji) * ( jk - 1./2. ) * &
672                                ht_i_b(ji) * dummy_fac2
673               ! weighting the profile
674               s_i_b(ji,jk) = zalpha * zs_zero(ji,jk) +       &
675                           ( 1.0 - zalpha ) * sm_i_b(ji)
676            END DO ! ji
677         END DO ! jk
678
679      ENDIF ! num_sal
680
681      !-------------------------------------------------------
682      ! Vertically varying salinity profile, constant in time
683      !-------------------------------------------------------
684      ! Schwarzacher (1959) multiyear salinity profile (mean = 2.30)
685
686      IF ( num_sal .EQ. 3 ) THEN
687
688         sm_i_b(:) = 2.30
689
690         DO ji = kideb, kiut
691            DO jk = 1, nlay_i
692               zargtemp  = ( jk - 0.5 ) / nlay_i
693               s_i_b(ji,jk)  =  1.6 - 1.6*cos(3.14169265*(zargtemp**(0.407/ &
694                                (0.573+zargtemp))))
695            END DO ! jk
696         END DO ! ji
697
698      ENDIF ! num_sal
699
700   END SUBROUTINE lim_var_salprof1d
701
702!===============================================================================
703
704#else
705   !!======================================================================
706   !!                       ***  MODULE limvar  ***
707   !!                          no sea ice model
708   !!======================================================================
709CONTAINS
710   SUBROUTINE lim_var_agg          ! Empty routines
711   END SUBROUTINE lim_var_agg
712   SUBROUTINE lim_var_glo2eqv      ! Empty routines
713   END SUBROUTINE lim_var_glo2eqv
714   SUBROUTINE lim_var_eqv2glo      ! Empty routines
715   END SUBROUTINE lim_var_eqv2glo
716   SUBROUTINE lim_var_salprof      ! Empty routines
717   END SUBROUTINE lim_var_salprof
718   SUBROUTINE lim_var_bv           ! Emtpy routines
719   END SUBROUTINE lim_var_bv 
720   SUBROUTINE lim_var_salprof1d    ! Emtpy routines
721   END SUBROUTINE lim_var_salprof1d
722
723#endif
724END MODULE limvar
Note: See TracBrowser for help on using the repository browser.