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.
isfstp.F90 in NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF – NEMO

source: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfstp.F90 @ 11521

Last change on this file since 11521 was 11521, checked in by mathiot, 5 years ago

ENHANCE-02_ISF: fix issue with ice sheet coupling and conservation + other minor changes (ticket #2142)
  • Property svn:keywords set to Id
File size: 14.5 KB
Line 
1MODULE isfstp
2   !!======================================================================
3   !!                       ***  MODULE  isfstp  ***
4   !! Surface module :  compute iceshelf load, melt and heat flux
5   !!======================================================================
6   !! History :  3.2  !  2011-02  (C.Harris  ) Original code isf cav
7   !!            X.X  !  2006-02  (C. Wang   ) Original code bg03
8   !!            3.4  !  2013-03  (P. Mathiot) Merging + parametrization
9   !!            4.1  !  2019-09  (P. Mathiot) Split param/explicit ice shelf and re-organisation
10   !!----------------------------------------------------------------------
11
12   !!----------------------------------------------------------------------
13   !!   isfstp       : compute iceshelf melt and heat flux
14   !!----------------------------------------------------------------------
15   USE oce            ! ocean dynamics and tracers
16   USE dom_oce        ! ocean space and time domain
17   USE phycst         ! physical constants
18   USE eosbn2         ! equation of state
19   USE sbc_oce        ! surface boundary condition: ocean fields
20   USE zdfdrg         ! vertical physics: top/bottom drag coef.
21   !
22   USE in_out_manager ! I/O manager
23   USE iom            ! I/O library
24   USE fldread        ! read input field at current time step
25   USE lbclnk         !
26   USE lib_fortran    ! glob_sum
27   !
28   USE isfrst         ! iceshelf restart
29   USE isftbl         ! ice shelf boundary layer
30   USE isfpar         ! ice shelf parametrisation
31   USE isfcav         ! ice shelf cavity
32   USE isfload        ! ice shelf load
33   USE isfcpl         ! isf variables
34   USE isf            ! isf variables
35   USE isfutils
36
37   IMPLICIT NONE
38
39   PRIVATE
40
41   PUBLIC   isf_stp, isf_stp_init  ! routine called in sbcmod and divhor
42
43   !!----------------------------------------------------------------------
44   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
45   !! $Id$
46   !! Software governed by the CeCILL license (see ./LICENSE)
47   !!----------------------------------------------------------------------
48CONTAINS
49 
50  SUBROUTINE isf_stp( kt )
51      !!---------------------------------------------------------------------
52      !!                  ***  ROUTINE isf_stp  ***
53      !!
54      !! ** Purpose : compute total heat flux and total fwf due to ice shelf melt
55      !!
56      !! ** Method  : For each case (parametrisation or explicity cavity) :
57      !!              - define the before fields
58      !!              - compute top boundary layer properties
59      !!                (in case of parametrisation, this is the
60      !!                 depth range model array used to compute mean far fields properties)
61      !!              - compute fluxes
62      !!              - write restart variables
63      !!
64      !!----------------------------------------------------------------------
65      INTEGER, INTENT(in) ::   kt   ! ocean time step
66      !!---------------------------------------------------------------------
67      !
68      IF ( ln_isfcav_mlt ) THEN
69         !
70         ! before time step
71         IF ( kt /= nit000 ) THEN
72            risf_cav_tsc_b (:,:,:) = risf_cav_tsc (:,:,:)
73            fwfisf_cav_b(:,:)      = fwfisf_cav(:,:)
74         END IF
75         !
76         ! compute tbl lvl/h
77         CALL isf_tbl_lvl(ht_n, e3t_n, misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav)
78         !
79         ! compute ice shelf melt
80         CALL isf_cav( kt, risf_cav_tsc, fwfisf_cav)
81         !
82         ! write restart variables (risf_cav_tsc, fwfisf for now and before)
83         IF (lrst_oce) CALL isfrst_write(kt, 'cav', risf_cav_tsc, fwfisf_cav)
84         !
85      END IF
86      !
87      IF ( ln_isfpar_mlt ) THEN
88         !
89         ! before time step
90         IF ( kt /= nit000 ) THEN
91            risf_par_tsc_b(:,:,:) = risf_par_tsc(:,:,:)
92            fwfisf_par_b  (:,:)   = fwfisf_par  (:,:)
93         END IF
94         !
95         ! compute misfkb, rhisf_tbl, rfrac (deepest level, thickness, fraction of deepest cell affected by tbl)
96         CALL isf_tbl_lvl(ht_n, e3t_n, misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par)
97         !
98         ! compute ice shelf melt
99         CALL isf_par( kt, risf_par_tsc, fwfisf_par)
100         !
101         ! write restart variables (qoceisf, qhcisf, fwfisf for now and before)
102         IF (lrst_oce) CALL isfrst_write(kt, 'par', risf_par_tsc, fwfisf_par)
103         !
104      END IF
105
106      IF ( ll_isfcpl ) THEN
107
108         IF (lrst_oce) CALL isfcpl_rst_write(kt)
109
110      END IF
111      !
112   END SUBROUTINE isf_stp
113
114   SUBROUTINE isf_stp_init
115      !!---------------------------------------------------------------------
116      !!                  ***  ROUTINE isfstp_init  ***
117      !!
118      !! ** Purpose :   Initialisation of the ice shelf public variables
119      !!
120      !! ** Method  :   Read the namsbc namelist and set derived parameters
121      !!                Call init routines for all other SBC modules that have one
122      !!
123      !! ** Action  : - read namsbc parameters
124      !!              - allocate memory
125      !!              - call cav/param init routine
126      !!----------------------------------------------------------------------
127      INTEGER               :: inum, ierror
128      INTEGER               :: ios                  ! Local integer output status for namelist read
129      INTEGER               :: ikt, ikb
130      INTEGER               :: ji, jj
131      !!----------------------------------------------------------------------
132      NAMELIST/namisf/ ln_isf       ,                                                                               & 
133         &             ln_isfcav_mlt, cn_isfcav_mlt, cn_gammablk, rn_gammat0, rn_gammas0, rn_htbl, sn_isfcav_fwf,   &
134         &             ln_isfpar_mlt, cn_isfpar_mlt, sn_isfpar_fwf, sn_isfpar_zmin, sn_isfpar_zmax, sn_isfpar_Leff, &
135         &             ln_isfcpl    , nn_drown, ln_isfcpl_cons,                                                     &
136         &             cn_isfload   , cn_isfdir
137      !!----------------------------------------------------------------------
138      !
139      ! Allocate public array
140      CALL isf_alloc()
141      !
142      riceload(:,:)       = 0.0_wp
143      fwfisf_oasis(:,:)   = 0.0_wp
144      fwfisf_par(:,:)     = 0.0_wp    ; fwfisf_par_b(:,:)     = 0.0_wp
145      fwfisf_cav(:,:)     = 0.0_wp    ; fwfisf_cav_b(:,:)     = 0.0_wp
146      risf_cav_tsc(:,:,:) = 0.0_wp    ; risf_cav_tsc_b(:,:,:) = 0.0_wp
147      risf_par_tsc(:,:,:) = 0.0_wp    ; risf_par_tsc_b(:,:,:) = 0.0_wp
148      !
149      REWIND( numnam_ref )              ! Namelist namsbc_rnf in reference namelist : Runoffs
150      READ  ( numnam_ref, namisf, IOSTAT = ios, ERR = 901)
151901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namisf in reference namelist', lwp )
152      !
153      REWIND( numnam_cfg )              ! Namelist namsbc_rnf in configuration namelist : Runoffs
154      READ  ( numnam_cfg, namisf, IOSTAT = ios, ERR = 902 )
155902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namisf in configuration namelist', lwp )
156      IF(lwm) WRITE ( numond, namisf )
157      !
158      IF (lwp) THEN
159         WRITE(numout,*)
160         WRITE(numout,*) 'isf_init : ice shelf initialisation'
161         WRITE(numout,*) '~~~~~~~~~~~~'
162         WRITE(numout,*) '   Namelist namisf :'
163         !
164         WRITE(numout,*) '   ice shelf cavity (open or parametrised)  ln_isf = ', ln_isf
165         !
166         IF ( ln_isf ) THEN
167            WRITE(numout,*) '      melt inside the cavity                  ln_isfcav_mlt   = ', ln_isfcav_mlt
168            IF ( ln_isfcav ) THEN
169               WRITE(numout,*) '         melt formulation                        cn_isfcav_mlt   = ', TRIM(cn_isfcav_mlt)
170               WRITE(numout,*) '         thickness of the top boundary layer     rn_htbl     = ', rn_htbl
171               WRITE(numout,*) '         gamma formulation                       cn_gammablk = ', TRIM(cn_gammablk) 
172               IF ( TRIM(cn_gammablk) .NE. 'spe' ) THEN
173                  WRITE(numout,*) '            gammat coefficient                      rn_gammat0  = ', rn_gammat0 
174                  WRITE(numout,*) '            gammas coefficient                      rn_gammas0  = ', rn_gammas0 
175                  WRITE(numout,*) '            top drag coef. used (from namdrg_top)   rn_Cd0      = ', r_Cdmin_top 
176               END IF
177            END IF
178            WRITE(numout,*) ''
179            !
180            WRITE(numout,*) '      ice shelf melt parametrisation          ln_isfpar_mlt    = ', ln_isfpar_mlt
181            IF ( ln_isfpar_mlt ) THEN
182               WRITE(numout,*) '         isf parametrisation formulation         cn_isfpar_mlt   = ', TRIM(cn_isfpar_mlt)
183            END IF
184            WRITE(numout,*) ''
185            !
186         ELSE
187            IF ( ln_isfcav ) THEN
188               WRITE(numout,*) ''
189               WRITE(numout,*) '   W A R N I N G: ice shelf cavities are open BUT no melt will be computed or read from file !'
190               WRITE(numout,*) ''
191            END IF
192         END IF
193
194         WRITE(numout,*) '      Coupling to an ice sheet model          ln_isfcpl         = ', ln_isfcpl
195         IF ( ln_isfcpl ) THEN
196            WRITE(numout,*) '         conservation activated ln_isfcpl_cons           = ', ln_isfcpl_cons
197            WRITE(numout,*) '            number of call of the extrapolation loop = ', nn_drown
198         ENDIF
199         !
200         IF (ln_isfcav) WRITE(numout,*) '      Ice shelf load method                   cn_isfload        = ', TRIM(cn_isfload)
201         WRITE(numout,*) ''
202
203      END IF
204      !
205      !---------------------------------------------------------------------------------------------------------------------
206      ! initialisation ice shelf load
207      IF ( ln_isfcav ) THEN
208         !
209         ! compute ice shelf mask
210         mskisf_cav(:,:) = (1._wp - tmask(:,:,1)) * ssmask(:,:)
211         !
212         ! compute ice shelf load
213         CALL isf_load( risfload )
214         !
215      END IF
216      !
217      !---------------------------------------------------------------------------------------------------------------------
218      ! sanity check
219      ! melt in the cavity without cavity
220      IF ( ln_isfcav_mlt .AND. (.NOT. ln_isfcav) ) &
221         &   CALL ctl_stop('ice shelf melt in the cavity activated (ln_isfcav_mlt) but no cavity detected in domcfg (ln_isfcav), STOP' )
222      !
223      ! ice sheet coupling without cavity
224      IF ( ln_isfcpl .AND. (.NOT. ln_isfcav) ) &
225         &   CALL ctl_stop('coupling with an ice sheet model detected (ln_isfcpl) but no cavity detected in domcfg (ln_isfcav), STOP' )
226      !
227      IF ( ln_isfcpl .AND. ln_isfcpl_cons .AND. ln_linssh ) &
228         &   CALL ctl_stop( 'The coupling between NEMO and an ice sheet model with the conservation option does not work with the linssh option' )
229      !
230      IF ( l_isfoasis ) THEN
231         !
232         CALL ctl_stop( ' ln_ctl and ice shelf not tested' )
233         !
234         ! NEMO coupled to ATMO model with isf cavity need oasis method for melt computation
235         IF ( ln_isfcav_mlt .AND. TRIM(cn_isfcav_mlt) /= 'oasis' ) CALL ctl_stop( 'cn_isfcav_mlt = oasis is the only option availble if fwf send by oasis' )
236         IF ( ln_isfpar_mlt .AND. TRIM(cn_isfpar_mlt) /= 'oasis' ) CALL ctl_stop( 'cn_isfpar_mlt = oasis is the only option availble if fwf send by oasis' )
237         !
238         ! oasis melt computation not tested (coded but not tested)
239         IF ( ln_isfcav_mlt .OR. ln_isfpar_mlt ) THEN
240            IF ( TRIM(cn_isfcav_mlt) == 'oasis' ) CALL ctl_stop( 'cn_isfcav_mlt = oasis not tested' )
241            IF ( TRIM(cn_isfpar_mlt) == 'oasis' ) CALL ctl_stop( 'cn_isfpar_mlt = oasis not tested' )
242         END IF
243         !
244         ! oasis melt computation with cavity open and cavity parametrised (not coded)
245         IF ( ln_isfcav_mlt .AND. ln_isfpar_mlt ) THEN
246            IF ( TRIM(cn_isfpar_mlt) == 'oasis' .AND. TRIM(cn_isfcav_mlt) == 'oasis' ) CALL ctl_stop( 'cn_isfpar_mlt = oasis and cn_isfcav_mlt = oasis not coded' )
247         END IF
248      END IF
249      !
250      ! terminate routine now if no ice shelf melt formulation specify
251      IF ( ln_isf ) THEN
252         !
253         ! initialisation useful variable
254         r1_Lfusisf =  1._wp / rLfusisf
255         !
256         ! initialisation melt in the cavity
257         IF ( ln_isfcav_mlt ) THEN
258            !
259            ! initialisation  of cav variable
260            CALL isf_cav_init()
261            !
262            ! read cav variable from restart
263            IF ( ln_rstart ) CALL isfrst_read('cav', risf_cav_tsc, fwfisf_cav, risf_cav_tsc_b, fwfisf_cav_b)
264            !
265         END IF
266         !
267         !---------------------------------------------------------------------------------------------------------------------
268         ! initialisation parametrised melt
269         IF ( ln_isfpar_mlt ) THEN
270            !
271            ! initialisation  of par variable
272            CALL isf_par_init()
273            !
274            ! read par variable from restart
275            IF ( ln_rstart ) CALL isfrst_read('par', risf_par_tsc, fwfisf_par, risf_par_tsc_b, fwfisf_par_b)
276            !
277         END IF
278      END IF
279      !
280      !---------------------------------------------------------------------------------------------------------------------
281      ! initialisation ice sheet coupling
282      !
283      ll_isfcpl     = .FALSE.
284      ll_isfcpl_cons= .FALSE.
285      !
286      IF( ln_isfcpl ) THEN
287
288         ! prepare writing restart
289         IF( lwxios ) CALL iom_set_rstw_var_active('ssmask')
290         IF( lwxios ) CALL iom_set_rstw_var_active('tmask')
291         !IF( lwxios ) CALL iom_set_rstw_var_active('wmask')
292         !IF( lwxios ) CALL iom_set_rstw_var_active('gdepw_n')
293         IF( lwxios ) CALL iom_set_rstw_var_active('e3t_n')
294         IF( lwxios ) CALL iom_set_rstw_var_active('e3u_n')
295         IF( lwxios ) CALL iom_set_rstw_var_active('e3v_n')
296
297         IF( ln_rstart ) THEN
298            !
299            ll_isfcpl = .TRUE.
300            !
301            CALL isf_alloc_cpl()
302            !
303            ! extrapolation tracer properties
304            CALL isfcpl_tra()
305            !
306            ! correction of the horizontal divergence and associated temp. and salt content flux
307            CALL isfcpl_vol()
308            !
309            ! apply the 'conservation' method
310            IF ( ln_isfcpl_cons )  THEN
311               ll_isfcpl_cons = .TRUE.
312               CALL isfcpl_cons()
313            END IF
314            !
315            ! Need to include in the cpl cons the isfrst_cpl_div contribution
316            ! decide how to manage thickness level change in conservation
317            !
318            tsb    (:,:,:,:) = tsn (:,:,:,:)
319            sshb   (:,:)     = sshn(:,:)
320            !
321         END IF
322      END IF
323         
324  END SUBROUTINE isf_stp_init
325   !!======================================================================
326END MODULE isfstp
Note: See TracBrowser for help on using the repository browser.