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.
isfpar.F90 in NEMO/trunk/src/OCE/ISF – NEMO

source: NEMO/trunk/src/OCE/ISF/isfpar.F90 @ 15055

Last change on this file since 15055 was 15004, checked in by mathiot, 3 years ago

ticket #2960: commit fix to the trunk (WARNING: output convention of isf fluxes changed from oce->isf to isf->oce), no impact on the input file needed for some options
File size: 8.7 KB
Line 
1MODULE isfpar
2   !!======================================================================
3   !!                       ***  MODULE  isfpar  ***
4   !! ice shelf module :  update ocean boundary condition under ice
5   !!                   shelf
6   !!======================================================================
7   !! History :  3.2  !  2011-02  (C.Harris  ) Original code isf cav
8   !!            X.X  !  2006-02  (C. Wang   ) Original code bg03
9   !!            3.4  !  2013-03  (P. Mathiot) Merging + parametrization
10   !!            4.1  !  2019-09  (P. Mathiot) Restructuration
11   !!            4.2  !  2021-05  (C. Ethe   ) Test and fix oasis case
12   !!----------------------------------------------------------------------
13
14   !!----------------------------------------------------------------------
15   !!   isfpar       : compute ice shelf melt using a prametrisation of ice shelf cavities
16   !!----------------------------------------------------------------------
17   USE isf_oce        ! ice shelf
18   !
19   USE isfrst   , ONLY: isfrst_write, isfrst_read ! ice shelf restart read/write subroutine
20   USE isftbl   , ONLY: isf_tbl_ktop, isf_tbl_lvl ! ice shelf top boundary layer properties subroutine
21   USE isfparmlt, ONLY: isfpar_mlt                ! ice shelf melt formulation subroutine
22   USE isfdiags , ONLY: isf_diags_flx             ! ice shelf diags subroutine
23   USE isfutils , ONLY: debug, read_2dcstdta      ! ice shelf debug subroutine
24   !
25   USE dom_oce  , ONLY: bathy          ! ocean space and time domain
26   USE par_oce  , ONLY: jpi,jpj        ! ocean space and time domain
27   USE phycst   , ONLY: r1_rho0_rcp    ! physical constants
28   !
29   USE in_out_manager ! I/O manager
30   USE iom            ! I/O library
31   USE fldread        ! read input field at current time step
32   USE lbclnk         ! lbc_lnk
33
34   IMPLICIT NONE
35   PRIVATE
36
37   PUBLIC   isf_par, isf_par_init
38
39   !!----------------------------------------------------------------------
40   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
41   !! $Id: sbcisf.F90 10536 2019-01-16 19:21:09Z mathiot $
42   !! Software governed by the CeCILL license (see ./LICENSE)
43   !!----------------------------------------------------------------------
44CONTAINS
45 
46   SUBROUTINE isf_par( kt, Kmm, ptsc, pqfwf )
47      !!---------------------------------------------------------------------
48      !!                     ***  ROUTINE isf_par ***     
49      !!
50      !! ** Purpose : compute the heat and fresh water due to ice shelf melting/freezing using a parametrisation
51      !!
52      !! ** Comment : in isf_par and all its call tree,
53      !!              'tbl' means parametrisation layer (ie how the far field temperature/salinity is computed)
54      !!              instead of in a proper top boundary layer as at the ice shelf ocean interface
55      !!              as the action to compute the properties of the tbl or the parametrisation layer are the same,
56      !!              (ie average T/S over a specific depth (can be across multiple levels))
57      !!              the name tbl was kept.
58      !!
59      !! ** Convention : all fluxes are from isf to oce
60      !!
61      !!---------------------------------------------------------------------
62      !!-------------------------- OUT --------------------------------------
63      REAL(wp), DIMENSION(jpi,jpj)     , INTENT(inout) :: pqfwf
64      REAL(wp), DIMENSION(jpi,jpj,jpts), INTENT(inout) :: ptsc
65      !!-------------------------- IN  --------------------------------------
66      INTEGER, INTENT(in) ::   kt                                     ! ocean time step
67      INTEGER, INTENT(in) ::   Kmm                                    ! ocean time level index
68      !!---------------------------------------------------------------------
69      REAL(wp), DIMENSION(jpi,jpj) :: zqoce, zqhc, zqlat, zqh
70      !!---------------------------------------------------------------------
71      !
72      ! compute heat content, latent heat and melt fluxes (2d)
73      CALL isfpar_mlt( kt, Kmm, zqhc, zqoce, pqfwf  )
74      !
75      ! compute heat and water flux (from isf to oce)
76      pqfwf(:,:) = pqfwf(:,:) * mskisf_par(:,:)
77      zqoce(:,:) = zqoce(:,:) * mskisf_par(:,:)
78      zqhc (:,:) = zqhc(:,:)  * mskisf_par(:,:)
79      !
80      ! compute latent heat flux (from isf to oce)
81      zqlat(:,:) = - pqfwf(:,:) * rLfusisf    ! 2d latent heat flux (W/m2)
82      !
83      ! total heat flux (from isf to oce)
84      zqh(:,:) = ( zqhc (:,:) + zqoce(:,:) )
85      !
86      ! lbclnk on melt and heat fluxes
87      CALL lbc_lnk( 'isfmlt', zqh, 'T', 1.0_wp, pqfwf, 'T', 1.0_wp)
88      !
89      ! output fluxes
90      CALL isf_diags_flx( Kmm, misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par, 'par', pqfwf, zqoce, zqlat, zqhc)
91      !
92      ! set temperature content
93      ptsc(:,:,jp_tem) = zqh(:,:) * r1_rho0_rcp
94      !
95      ! write restart variables (qoceisf, qhcisf, fwfisf for now and before)
96      IF (lrst_oce) CALL isfrst_write(kt, 'par', ptsc, pqfwf)
97      !
98      IF ( ln_isfdebug ) THEN
99         IF(lwp) WRITE(numout,*)
100         CALL debug('isf_par: ptsc T',ptsc(:,:,1))
101         CALL debug('isf_par: ptsc S',ptsc(:,:,2))
102         CALL debug('isf_par: pqfwf fwf',pqfwf(:,:))
103         IF(lwp) WRITE(numout,*)
104      END IF
105      !
106   END SUBROUTINE isf_par
107
108   SUBROUTINE isf_par_init
109      !!---------------------------------------------------------------------
110      !!                  ***  ROUTINE isf_par_init  ***
111      !!
112      !! ** Purpose : initialisation of the variable needed for the parametrisation of ice shelf melt
113      !!
114      !!----------------------------------------------------------------------
115      INTEGER               :: ierr
116      REAL(wp), DIMENSION(jpi,jpj) :: ztblmax, ztblmin
117      !!----------------------------------------------------------------------
118      !
119      ! allocation
120      CALL isf_alloc_par()
121      !
122      ! initialisation
123      misfkt_par(:,:)     = 1         ; misfkb_par(:,:)       = 1         
124      rhisf_tbl_par(:,:)  = 1e-20     ; rfrac_tbl_par(:,:)    = 0.0_wp
125      !
126      ! define isf tbl tickness, top and bottom indice
127      CALL read_2dcstdta(TRIM(sn_isfpar_zmax%clname), TRIM(sn_isfpar_zmax%clvar), ztblmax)
128      CALL read_2dcstdta(TRIM(sn_isfpar_zmin%clname), TRIM(sn_isfpar_zmin%clvar), ztblmin)
129      !
130      ! mask ice shelf parametrisation location
131      ztblmax(:,:) = ztblmax(:,:) * ssmask(:,:)
132      ztblmin(:,:) = ztblmin(:,:) * ssmask(:,:)
133      !
134      ! if param used under an ice shelf overwrite ztblmin by the ice shelf draft
135      WHERE ( risfdep > 0._wp .AND. ztblmin > 0._wp )
136         ztblmin(:,:) = risfdep(:,:)
137      END WHERE
138      !
139      ! ensure ztblmax <= bathy
140      WHERE ( ztblmax(:,:) > bathy(:,:) )
141         ztblmax(:,:) = bathy(:,:)
142      END WHERE
143      !
144      ! compute ktop and update ztblmin to gdepw_0(misfkt_par)
145      CALL isf_tbl_ktop(ztblmin, misfkt_par) !   out: misfkt_par
146      !                                      ! inout: ztblmin
147      !
148      ! initial tbl thickness
149      rhisf0_tbl_par(:,:) = ztblmax(:,:) - ztblmin(:,:)
150      !
151      ! define iceshelf parametrisation mask
152      mskisf_par = 0
153      WHERE ( rhisf0_tbl_par(:,:) > 0._wp )
154         mskisf_par(:,:) = 1._wp
155      END WHERE
156      !
157      ! read par variable from restart
158      IF ( ln_rstart ) CALL isfrst_read('par', risf_par_tsc, fwfisf_par, risf_par_tsc_b, fwfisf_par_b)
159      !
160      SELECT CASE ( TRIM(cn_isfpar_mlt) )
161         !
162      CASE ( 'spe' )
163         !
164         ALLOCATE( sf_isfpar_fwf(1), STAT=ierr )
165         ALLOCATE( sf_isfpar_fwf(1)%fnow(jpi,jpj,1), sf_isfpar_fwf(1)%fdta(jpi,jpj,1,2) )
166         CALL fld_fill( sf_isfpar_fwf, (/ sn_isfpar_fwf /), cn_isfdir, 'isf_par_init', 'read fresh water flux isf data', 'namisf' )
167         !
168         IF(lwp) WRITE(numout,*)
169         IF(lwp) WRITE(numout,*) '      ==>>>   ice melt read from forcing field (cn_isfmlt_par = spe)'
170         !
171      CASE ( 'bg03' )
172         !
173         IF(lwp) WRITE(numout,*)
174         IF(lwp) WRITE(numout,*) '      ==>>>   bg03 parametrisation (cn_isfmlt_par = bg03)'
175         !
176         ! read effective length
177         CALL read_2dcstdta(TRIM(sn_isfpar_Leff%clname), TRIM(sn_isfpar_Leff%clvar), risfLeff)
178         risfLeff = risfLeff*1000.0_wp           !: convertion in m
179         !
180      CASE ( 'oasis' )
181         !
182         ALLOCATE( sf_isfpar_fwf(1), STAT=ierr )
183         ALLOCATE( sf_isfpar_fwf(1)%fnow(jpi,jpj,1), sf_isfpar_fwf(1)%fdta(jpi,jpj,1,2) )
184         CALL fld_fill( sf_isfpar_fwf, (/ sn_isfpar_fwf /), cn_isfdir, 'isf_par_init', 'read fresh water flux isf data', 'namisf' )
185         !
186         IF(lwp) WRITE(numout,*)
187         IF(lwp) WRITE(numout,*) '      ==>>>    isf melt provided by OASIS (cn_isfmlt_par = oasis)'
188         !
189      CASE DEFAULT
190         CALL ctl_stop( 'sbc_isf_init: wrong value of nn_isf' )
191      END SELECT
192      !
193   END SUBROUTINE isf_par_init
194
195END MODULE isfpar
Note: See TracBrowser for help on using the repository browser.