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

source: trunk/NEMOGCM/NEMO/LIM_SRC_3/limvar.F90 @ 2528

Last change on this file since 2528 was 2528, checked in by rblod, 13 years ago

Update NEMOGCM from branch nemo_v3_3_beta

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