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

Last change on this file since 4990 was 4990, checked in by timgraham, 7 years ago

Merged branches/2014/dev_MERGE_2014 back onto the trunk as follows:

In the working copy of branch ran:
svn merge svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/trunk@HEAD
1 conflict in LIM_SRC_3/limdiahsb.F90
Resolved by keeping the version from dev_MERGE_2014 branch
and commited at r4989

In working copy run:
svn switch svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/trunk
to switch working copy

Run:
svn merge —reintegrate svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/branches/2014/dev_MERGE_2014
to merge the branch into the trunk - no conflicts at this stage.

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