source: NEMO/branches/2020/dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation/src/OCE/ISF/isfparmlt.F90 @ 13257

Last change on this file since 13257 was 12489, checked in by davestorkey, 9 months ago

Preparation for new timestepping scheme #2390.
Main changes:

  1. Initial euler timestep now handled in stp and not in TRA/DYN routines.
  2. Renaming of all timestep parameters. In summary, the namelist parameter is now rn_Dt and the current timestep is rDt (and rDt_ice, rDt_trc etc).
  3. Renaming of a few miscellaneous parameters, eg. atfp → rn_atfp (namelist parameter used everywhere) and rau0 → rho0.

This version gives bit-comparable results to the previous version of the trunk.

File size: 11.6 KB
Line 
1MODULE isfparmlt
2   !!======================================================================
3   !!                       ***  MODULE  isfparmlt  ***
4   !! Ice shelf parametrisation module :  update surface ocean boundary condition under ice
5   !!                   shelf using an ice shelf melt parametrisation
6   !!======================================================================
7   !! History :  4.0  ! original code
8   !!----------------------------------------------------------------------
9
10   USE isf_oce                  ! ice shelf
11   USE isftbl , ONLY: isf_tbl   ! ice shelf depth average
12
13   USE dom_oce                  ! ocean space and time domain
14   USE oce    , ONLY: ts        ! ocean dynamics and tracers
15   USE phycst , ONLY: rcp, rho0 ! physical constants
16   USE eosbn2 , ONLY: eos_fzp   ! equation of state
17
18   USE in_out_manager              ! I/O manager
19   USE iom        , ONLY: iom_put  ! I/O library
20   USE fldread    , ONLY: fld_read, FLD, FLD_N !
21   USE lib_fortran, ONLY: glob_sum !
22   USE lib_mpp    , ONLY: ctl_stop !
23
24   IMPLICIT NONE
25
26   PRIVATE
27
28   PUBLIC  isfpar_mlt 
29   
30   !!----------------------------------------------------------------------
31   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
32   !! $Id: sbcisf.F90 10536 2019-01-16 19:21:09Z mathiot $
33   !! Software governed by the CeCILL license (see ./LICENSE)
34   !!----------------------------------------------------------------------
35CONTAINS
36
37! -------------------------------------------------------------------------------------------------------
38! -------------------------------- PUBLIC SUBROUTINE ----------------------------------------------------
39! -------------------------------------------------------------------------------------------------------
40
41  SUBROUTINE isfpar_mlt( kt, Kmm, pqhc, pqoce, pqfwf )
42      !!---------------------------------------------------------------------
43      !!                  ***  ROUTINE isfpar_mlt  ***
44      !!
45      !! ** Purpose : Compute Salt and Heat fluxes related to ice_shelf
46      !!              melting and freezing
47      !!
48      !! ** Method  :  2 parameterizations are available according
49      !!                        1 : Specified melt flux
50      !!                        2 : Beckmann & Goose parameterization
51      !!
52      !!-------------------------- OUT -------------------------------------
53      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pqfwf, pqoce, pqhc  ! fresh water, ice-ocean heat and heat content fluxes
54      !!-------------------------- IN  -------------------------------------
55      INTEGER, INTENT(in) ::   kt   ! ocean time step
56      INTEGER, INTENT(in) ::   Kmm  ! ocean time level index
57      !!---------------------------------------------------------------------
58      !
59      ! Choose among the available ice shelf parametrisation
60      SELECT CASE ( cn_isfpar_mlt )
61      CASE ( 'spe' )    ! specified runoff in depth (Mathiot et al., 2017 in preparation)
62         CALL isfpar_mlt_spe(kt, Kmm, pqhc, pqoce, pqfwf)
63      CASE ( 'bg03' )    ! Beckmann and Goosse parametrisation
64         CALL isfpar_mlt_bg03(kt, Kmm, pqhc, pqoce, pqfwf)
65      CASE ( 'oasis' )
66         CALL isfpar_mlt_oasis( kt, Kmm, pqhc, pqoce, pqfwf)
67      CASE DEFAULT
68         CALL ctl_stop('STOP', 'unknown isf melt formulation : cn_isfpar (should not see this)')
69      END SELECT
70      !
71   END SUBROUTINE isfpar_mlt
72
73! -------------------------------------------------------------------------------------------------------
74! -------------------------------- PRIVATE SUBROUTINE ---------------------------------------------------
75! -------------------------------------------------------------------------------------------------------
76
77   SUBROUTINE isfpar_mlt_spe(kt, Kmm, pqhc, pqoce, pqfwf)
78      !!---------------------------------------------------------------------
79      !!                  ***  ROUTINE isfpar_mlt_spe  ***
80      !!
81      !! ** Purpose : prescribed ice shelf melting in case ice shelf cavities are closed.
82      !!              data read into a forcing files.
83      !!
84      !!-------------------------- OUT -------------------------------------
85      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pqhc, pqfwf, pqoce  ! fresh water and ice-ocean heat fluxes
86      !!-------------------------- IN  -------------------------------------
87      INTEGER,  INTENT(in) :: kt
88      INTEGER,  INTENT(in) :: Kmm    !  ocean time level index
89      !!--------------------------------------------------------------------
90      INTEGER :: jk
91      REAL(wp), DIMENSION(jpi,jpj,jpk)  :: ztfrz3d
92      REAL(wp), DIMENSION(jpi,jpj)      :: ztfrz
93      !!--------------------------------------------------------------------
94      !
95      ! 0. ------------Read specified runoff
96      CALL fld_read ( kt, 1, sf_isfpar_fwf   )
97      !
98      ! compute ptfrz
99      ! 1. ------------Mean freezing point
100      DO jk = 1,jpk
101         CALL eos_fzp(ts(:,:,jk,jp_sal,Kmm), ztfrz3d(:,:,jk), gdept(:,:,jk,Kmm))
102      END DO
103      CALL isf_tbl(Kmm, ztfrz3d, ztfrz, 'T', misfkt_par, rhisf_tbl_par, misfkb_par, rfrac_tbl_par )
104      !
105      pqfwf(:,:) = - sf_isfpar_fwf(1)%fnow(:,:,1)      ! fresh water flux from the isf (fwfisf <0 mean melting)
106      pqoce(:,:) =   pqfwf(:,:) * rLfusisf             ! ocean/ice shelf flux assume to be equal to latent heat flux
107      pqhc (:,:) =   pqfwf(:,:) * ztfrz(:,:) * rcp     ! heat content flux
108      !
109      CALL iom_put('isftfrz_par', ztfrz )
110      !
111   END SUBROUTINE isfpar_mlt_spe
112
113   SUBROUTINE isfpar_mlt_bg03(kt, Kmm, pqhc, pqoce, pqfwf)
114      !!---------------------------------------------------------------------
115      !!                  ***  ROUTINE isfpar_mlt_bg03  ***
116      !!
117      !! ** Purpose : compute an estimate of ice shelf melting and
118      !!              latent, ocean-ice and heat content heat fluxes
119      !!              in case cavities are closed based on the far fields T and S properties.
120      !!
121      !! ** Method  : The ice shelf melt is computed as proportional to the differences between the
122      !!              mean temperature and mean freezing point in front of the ice shelf averaged
123      !!              over the ice shelf min ice shelf draft and max ice shelf draft and the freezing point
124      !!
125      !! ** Reference : Beckmann and Goosse (2003), "A parameterization of ice shelf-ocean
126      !!                interaction for climate models", Ocean Modelling 5(2003) 157-170.
127      !!----------------------------------------------------------------------
128      !!-------------------------- OUT -------------------------------------
129      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pqhc, pqfwf, pqoce  ! fresh water and ice-ocean heat fluxes
130      !!-------------------------- IN  -------------------------------------
131      INTEGER,  INTENT(in) :: kt
132      INTEGER,  INTENT(in) :: Kmm    !  ocean time level index
133      !!--------------------------------------------------------------------
134      INTEGER :: jk
135      REAL(wp), DIMENSION(jpi,jpj,jpk) :: ztfrz3d        ! freezing point
136      REAL(wp), DIMENSION(jpi,jpj)     :: ztfrz          ! freezing point
137      REAL(wp), DIMENSION(jpi,jpj)     :: ztavg          ! temperature avg
138      !!----------------------------------------------------------------------
139      !
140      ! 0. ------------Mean freezing point
141      DO jk = 1,jpk
142         CALL eos_fzp(ts(:,:,jk,jp_sal,Kmm), ztfrz3d(:,:,jk), gdept(:,:,jk,Kmm))
143      END DO
144      CALL isf_tbl(Kmm, ztfrz3d, ztfrz, 'T', misfkt_par, rhisf_tbl_par, misfkb_par, rfrac_tbl_par )
145      !
146      ! 1. ------------Mean temperature
147      CALL isf_tbl(Kmm, ts(:,:,jk,jp_tem,Kmm), ztavg, 'T', misfkt_par, rhisf_tbl_par, misfkb_par, rfrac_tbl_par )
148      !
149      ! 2. ------------Net heat flux and fresh water flux due to the ice shelf
150      pqoce(:,:) =   rho0 * rcp * rn_gammat0 * risfLeff(:,:) * e1t(:,:) * ( ztavg(:,:) - ztfrz(:,:) ) * r1_e1e2t(:,:)
151      pqfwf(:,:) = - pqoce(:,:) / rLfusisf             ! derived from the latent heat flux
152      pqhc (:,:) =   pqfwf(:,:) * ztfrz(:,:) * rcp     ! heat content flux
153      !
154      ! 3. ------------BG03 output
155      ! output ttbl
156      CALL iom_put('ttbl_par', ztavg(:,:) * mskisf_par(:,:) )
157      !
158      ! output thermal driving
159      CALL iom_put('isfthermald_par',( ztfrz(:,:) - ztavg(:,:) ) * mskisf_par(:,:))
160      !
161      ! output freezing point used to define the thermal driving and heat content fluxes
162      CALL iom_put('isftfrz_par', ztfrz )
163      !
164   END SUBROUTINE isfpar_mlt_bg03
165
166   SUBROUTINE isfpar_mlt_oasis(kt, Kmm, pqhc , pqoce, pqfwf )
167      !!----------------------------------------------------------------------
168      !!                  ***  ROUTINE isfpar_mlt_oasis  ***
169      !!
170      !! ** Purpose    : scale the fwf read from input file by the total amount received by the sbccpl interface
171      !!
172      !! ** Purpose    : - read ice shelf melt from forcing file and scale it by the input file total amount => pattern
173      !!                 - compute total amount of fwf given by sbccpl (fwfisf_oasis)
174      !!                 - scale fwf and compute heat fluxes
175      !!
176      !!---------------------------------------------------------------------
177      !!-------------------------- OUT -------------------------------------
178      REAL(wp), DIMENSION(jpi,jpj), INTENT(  out) :: pqhc, pqoce, pqfwf  ! heat content, latent heat and fwf fluxes
179      !!-------------------------- IN  -------------------------------------
180      INTEGER                     , INTENT(in   ) :: kt                  ! current time step
181      INTEGER                     , INTENT(in   ) :: Kmm                 !  ocean time level index
182      !!--------------------------------------------------------------------
183      INTEGER                           :: jk                            ! loop index
184      REAL(wp)                          :: zfwf_fld, zfwf_oasis          ! total fwf in the forcing fields (pattern) and from the cpl interface (amount)
185      REAL(wp), DIMENSION(jpi,jpj)      :: ztfrz                         ! tbl freezing temperature
186      REAL(wp), DIMENSION(jpi,jpj)      :: zfwf                          ! 2d fwf map after scaling
187      REAL(wp), DIMENSION(jpi,jpj,jpk)  :: ztfrz3d
188      !!--------------------------------------------------------------------
189      !
190      ! 0. ------------Read specified runoff
191      CALL fld_read ( kt, 1, sf_isfpar_fwf   )
192      !
193      ! 1. ------------Mean freezing point (needed for heat content flux)
194      DO jk = 1,jpk
195         CALL eos_fzp(ts(:,:,jk,jp_sal,Kmm), ztfrz3d(:,:,jk), gdept(:,:,jk,Kmm))
196      END DO
197      CALL isf_tbl(Kmm, ztfrz3d, ztfrz, 'T', misfkt_par, rhisf_tbl_par, misfkb_par, rfrac_tbl_par )
198      !
199      ! 2. ------------Scale isf melt pattern with total amount from oasis
200      ! ice shelf 2d map
201      zfwf(:,:) = - sf_isfpar_fwf(1)%fnow(:,:,1)
202      !
203      ! compute glob sum from input file
204      ! (PM) should we consider delay sum as in fwb ? (it will offset by 1 time step if I understood well)
205      zfwf_fld = glob_sum('isfcav_mlt', e1e2t(:,:) * zfwf(:,:))
206      !
207      ! compute glob sum from atm->oce ice shelf fwf
208      ! (PM) should we consider delay sum as in fwb ?
209      zfwf_oasis = glob_sum('isfcav_mlt', e1e2t(:,:) * fwfisf_oasis(:,:))
210      !
211      ! scale fwf
212      zfwf(:,:) = zfwf(:,:) * zfwf_oasis / zfwf_fld
213      !
214      ! 3. -----------Define fwf and qoce
215      ! ocean heat flux is assume to be equal to the latent heat
216      pqfwf(:,:) =   zfwf(:,:)                         ! fwf                ( >0 out )
217      pqoce(:,:) = - pqfwf(:,:) * rLfusisf             ! ocean heat flux    ( >0 out ) (assumed to be the latent heat flux)
218      pqhc (:,:) =   pqfwf(:,:) * ztfrz(:,:) * rcp     ! heat content flux  ( >0 out )
219      !
220      CALL iom_put('isftfrz_par', ztfrz )
221      !
222   END SUBROUTINE isfpar_mlt_oasis
223
224END MODULE isfparmlt
Note: See TracBrowser for help on using the repository browser.