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.
limthd_da.F90 in branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3 – NEMO

source: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/limthd_da.F90 @ 8327

Last change on this file since 8327 was 8327, checked in by clem, 7 years ago

STEP4 (3): put all thermodynamics in 1D (limthd_da OK)

File size: 11.6 KB
Line 
1MODULE limthd_da
2   !!======================================================================
3   !!                       ***  MODULE limthd_da ***
4   !! LIM-3 sea-ice :  computation of lateral melting in the ice
5   !!======================================================================
6   !! History :   4.0   ! 2016-03 (C. Rousset) original code
7   !!---------------------------------------------------------------------
8#if defined key_lim3
9   !!----------------------------------------------------------------------
10   !!   'key_lim3'                                      LIM-3 sea-ice model
11   !!----------------------------------------------------------------------
12   !!   lim_thd_da   : sea ice lateral melting
13   !!----------------------------------------------------------------------
14   USE par_oce        ! ocean parameters
15   USE phycst         ! physical constants (ocean directory)
16   USE sbc_oce , ONLY : sst_m
17   USE ice            ! LIM variables
18   USE thd_ice        ! thermodynamic sea-ice variables
19   USE limtab         ! 1D <==> 2D transformation
20   !
21   USE lib_mpp        ! MPP library
22   USE wrk_nemo       ! work arrays
23   USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 
24
25   IMPLICIT NONE
26   PRIVATE
27
28   PUBLIC   lim_thd_da        ! called by limthd module
29
30   !!----------------------------------------------------------------------
31   !! NEMO/LIM3 4.0 , UCL - NEMO Consortium (2011)
32   !! $Id: limthd_da.F90 5123 2015-03-04 16:06:03Z clem $
33   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
34   !!----------------------------------------------------------------------
35CONTAINS
36
37   SUBROUTINE lim_thd_da
38      !!-------------------------------------------------------------------
39      !!                ***  ROUTINE lim_thd_da  ***   
40      !!   
41      !! ** Purpose :   computes sea ice lateral melting
42      !!
43      !! ** Method  :   dA/dt = - P * W   [s-1]
44      !!                   W = melting velocity [m.s-1]
45      !!                   P = perimeter of ice-ocean lateral interface normalized by grid cell area [m.m-2]
46      !!
47      !!                   W = m1 * (Tw -Tf)**m2                    --- originally from Josberger 1979 ---
48      !!                      (Tw - Tf) = elevation of water temp above freezing
49      !!                      m1 and m2 = (1.6e-6 , 1.36) best fit from field experiment near the coast of Prince Patrick Island (Perovich 1983) => static ice
50      !!                      m1 and m2 = (3.0e-6 , 1.36) best fit from MIZEX 84 experiment (Maykut and Perovich 1987) => moving ice
51      !!
52      !!                   P = N * pi * D                           --- from Rothrock and Thorndike 1984 ---
53      !!                      D = mean floe caliper diameter
54      !!                      N = number of floes = ice area / floe area(average) = A / (Cs * D**2)
55      !!                         A = ice concentration
56      !!                         Cs = deviation from a square (square:Cs=1 ; circle:Cs=pi/4 ; floe:Cs=0.66)
57      !!
58      !!                   D = Dmin * ( Astar / (Astar-A) )**beta   --- from Lupkes et al., 2012 (eq. 26-27) ---
59      !!                                                             
60      !!                      Astar = 1 / ( 1 - (Dmin/Dmax)**(1/beta) )
61      !!                      Dmin = minimum floe diameter (recommended to be 8m +- 20%)
62      !!                      Dmax = maximum floe diameter (recommended to be 300m, but it does not impact melting much except for Dmax<100m)
63      !!                      beta = 1.0 +-20% (recommended value)
64      !!                           = 0.3 best fit for western Fram Strait and Antarctica
65      !!                           = 1.4 best fit for eastern Fram Strait
66      !!
67      !! ** Tunable parameters  :   We propose to tune the lateral melting via 2 parameters
68      !!                               Dmin [6-10m]   => 6  vs 8m = +40% melting at the peak (A~0.5)
69      !!                                                 10 vs 8m = -20% melting
70      !!                               beta [0.8-1.2] => decrease = more melt and melt peaks toward higher concentration
71      !!                                                                  (A~0.5 for beta=1 ; A~0.8 for beta=0.2)
72      !!                                                 0.3 = best fit for western Fram Strait and Antarctica
73      !!                                                 1.4 = best fit for eastern Fram Strait
74      !!
75      !! ** Note   :   Former and more simple formulations for floe diameters can be found in Mai (1995),
76      !!               Birnbaum and Lupkes (2002), Lupkes and Birnbaum (2005). They are reviewed in Lupkes et al 2012
77      !!               A simpler implementation for CICE can be found in Bitz et al (2001) and Tsamados et al (2015)
78      !!
79      !! ** References
80      !!    Bitz, C. M., Holland, M. M., Weaver, A. J., & Eby, M. (2001).
81      !!              Simulating the ice‐thickness distribution in a coupled climate model.
82      !!              Journal of Geophysical Research: Oceans, 106(C2), 2441-2463.
83      !!    Josberger, E. G. (1979).
84      !!              Laminar and turbulent boundary layers adjacent to melting vertical ice walls in salt water
85      !!              (No. SCIENTIFIC-16). WASHINGTON UNIV SEATTLE DEPT OF ATMOSPHERIC SCIENCES.
86      !!    Lüpkes, C., Gryanik, V. M., Hartmann, J., & Andreas, E. L. (2012).
87      !!              A parametrization, based on sea ice morphology, of the neutral atmospheric drag coefficients
88      !!              for weather prediction and climate models.
89      !!              Journal of Geophysical Research: Atmospheres, 117(D13).
90      !!    Maykut, G. A., & Perovich, D. K. (1987).
91      !!              The role of shortwave radiation in the summer decay of a sea ice cover.
92      !!              Journal of Geophysical Research: Oceans, 92(C7), 7032-7044.
93      !!    Perovich, D. K. (1983).
94      !!              On the summer decay of a sea ice cover. (Doctoral dissertation, University of Washington).
95      !!    Rothrock, D. A., & Thorndike, A. S. (1984).
96      !!              Measuring the sea ice floe size distribution.
97      !!              Journal of Geophysical Research: Oceans, 89(C4), 6477-6486.
98      !!    Tsamados, M., Feltham, D., Petty, A., Schroeder, D., & Flocco, D. (2015).
99      !!              Processes controlling surface, bottom and lateral melt of Arctic sea ice in a state of the art sea ice model.
100      !!              Phil. Trans. R. Soc. A, 373(2052), 20140167.
101      !!---------------------------------------------------------------------
102      INTEGER             ::   ji, jj, jk, jl     ! dummy loop indices
103      INTEGER             ::   nidx
104      REAL(wp)            ::   zastar, zdfloe, zperi, zwlat, zda
105      REAL(wp), PARAMETER ::   zdmax = 300._wp
106      REAL(wp), PARAMETER ::   zcs   = 0.66_wp
107      REAL(wp), PARAMETER ::   zm1   = 3.e-6_wp
108      REAL(wp), PARAMETER ::   zm2   = 1.36_wp
109      !
110      REAL(wp), DIMENSION(jpij) ::   zda_tot
111      !!---------------------------------------------------------------------
112
113      ! select ice covered grid points
114      nidx = 0 ; idxice(:) = 0
115      DO jj = 2, jpjm1
116         DO ji = 2, jpim1
117            IF ( at_i(ji,jj) > epsi10 ) THEN
118               nidx         = nidx  + 1
119               idxice(nidx) = (jj - 1) * jpi + ji
120            ENDIF
121         END DO
122      END DO
123
124      !------------------------------------------------------------!
125      ! --- Calculate reduction of total sea ice concentration --- !
126      !------------------------------------------------------------!
127      zastar = 1._wp / ( 1._wp - (rn_dmin / zdmax)**(1._wp/rn_beta) )
128
129      CALL tab_2d_1d( nidx, at_i_1d(1:nidx), at_i , jpi, jpj, idxice(1:nidx) )
130      CALL tab_2d_1d( nidx, t_bo_1d(1:nidx), t_bo , jpi, jpj, idxice(1:nidx) )
131      CALL tab_2d_1d( nidx, sst_1d (1:nidx), sst_m, jpi, jpj, idxice(1:nidx) )
132      DO ji = 1, nidx   
133         zdfloe = rn_dmin * ( zastar / ( zastar - at_i_1d(ji) ) )**rn_beta         ! Mean floe caliper diameter [m]
134         zperi  = at_i_1d(ji) * rpi / ( zcs * zdfloe )                             ! Mean perimeter of the floe = N*pi*D = (A/cs*D^2)*pi*D [m.m-2]
135         zwlat  = zm1 * ( MAX( 0._wp, sst_1d(ji) - ( t_bo_1d(ji) - rt0 ) ) )**zm2  ! Melt speed rate [m/s]
136         
137         zda_tot(ji) = - MIN( zwlat * zperi * rdt_ice, at_i_1d(ji) )               ! sea ice concentration decrease
138      END DO
139     
140      !---------------------------------------------------------------------------------------------!
141      ! --- Distribute reduction among ice categories and calculate associated ice-ocean fluxes --- !
142      !---------------------------------------------------------------------------------------------!
143      DO jl = 1, jpl
144         CALL tab_2d_1d( nidx, a_i_1d    (1:nidx), a_i(:,:,jl) , jpi, jpj, idxice(1:nidx) )
145         CALL tab_2d_1d( nidx, ht_i_1d   (1:nidx), ht_i(:,:,jl), jpi, jpj, idxice(1:nidx) )
146         CALL tab_2d_1d( nidx, sm_i_1d   (1:nidx), sm_i(:,:,jl), jpi, jpj, idxice(1:nidx) )
147         CALL tab_2d_1d( nidx, sfx_lam_1d(1:nidx), sfx_lam     , jpi, jpj, idxice(1:nidx) )
148         CALL tab_2d_1d( nidx, hfx_thd_1d(1:nidx), hfx_thd     , jpi, jpj, idxice(1:nidx) )
149         CALL tab_2d_1d( nidx, wfx_lam_1d(1:nidx), wfx_lam     , jpi, jpj, idxice(1:nidx) )
150         DO jk = 1, nlay_i
151            CALL tab_2d_1d( nidx, e_i_1d(1:nidx,jk), e_i(:,:,jk,jl)  , jpi, jpj, idxice(1:nidx) )
152         END DO
153         DO jk = 1, nlay_s
154            CALL tab_2d_1d( nidx, e_s_1d(1:nidx,jk), e_s(:,:,jk,jl)  , jpi, jpj, idxice(1:nidx) )
155         END DO
156
157         DO ji = 1, nidx
158            ! decrease of concentration for the category jl
159            !    each category contributes to melting in proportion to its concentration
160            zda     = zda_tot(ji) * a_i_1d(ji) / at_i_1d(ji)
161           
162            ! Contribution to salt flux
163            sfx_lam_1d(ji) = sfx_lam_1d(ji) - rhoic *  ht_i_1d(ji) * zda * sm_i_1d(ji) * r1_rdtice
164           
165            ! Contribution to heat flux into the ocean [W.m-2], <0 
166            !clemX               hfx_thd_1d(ji) = hfx_thd_1d(ji) + zda * r1_rdtice * ( ht_i_1d(ji) * SUM( e_i_1d(ji,:) ) * r1_nlay_i  &
167            !                  &                                                     + ht_s_1d(ji) *      e_s_1d(ji,1)   * r1_nlay_s )
168            hfx_thd_1d(ji) = hfx_thd_1d(ji) + zda_tot(ji) / at_i_1d(ji) * r1_rdtice * ( SUM( e_i_1d(ji,:) ) + e_s_1d(ji,1) )
169           
170            ! Contribution to mass flux
171            wfx_lam_1d(ji) =  wfx_lam_1d(ji) - zda * r1_rdtice * ( rhoic * ht_i_1d(ji) + rhosn * ht_s_1d(ji) )
172           
173            ! new concentration
174            a_i_1d(ji) = a_i_1d(ji) + zda
175
176            ! ensure that ht_i = 0 where a_i = 0
177            IF( a_i_1d(ji) == 0._wp )   ht_i_1d(ji) = 0._wp 
178         END DO
179
180         CALL tab_1d_2d( nidx, a_i (:,:,jl), idxice, a_i_1d     (1:nidx), jpi, jpj )
181         CALL tab_1d_2d( nidx, ht_i(:,:,jl), idxice, ht_i_1d    (1:nidx), jpi, jpj )
182         CALL tab_1d_2d( nidx, sfx_lam     , idxice, sfx_lam_1d(1:nidx) , jpi, jpj )
183         CALL tab_1d_2d( nidx, hfx_thd     , idxice, hfx_thd_1d(1:nidx) , jpi, jpj )
184         CALL tab_1d_2d( nidx, wfx_lam     , idxice, wfx_lam_1d(1:nidx) , jpi, jpj )
185
186      END DO
187           
188      !
189   END SUBROUTINE lim_thd_da
190   
191#else
192   !!----------------------------------------------------------------------
193   !!   Default option         Dummy Module          No LIM-3 sea-ice model
194   !!----------------------------------------------------------------------
195CONTAINS
196   SUBROUTINE lim_thd_da          ! Empty routine
197   END SUBROUTINE lim_thd_da
198#endif
199   !!======================================================================
200END MODULE limthd_da
Note: See TracBrowser for help on using the repository browser.