[11395] | 1 | MODULE isfcavmlt |
---|
| 2 | !!====================================================================== |
---|
| 3 | !! *** MODULE isfcav_mlt *** |
---|
| 4 | !! ice shelf module : update surface ocean boundary condition under ice |
---|
| 5 | !! shelf |
---|
| 6 | !!====================================================================== |
---|
| 7 | !! History : 4.0 ! 2019-09 (P. Mathiot) Original code |
---|
| 8 | !!---------------------------------------------------------------------- |
---|
| 9 | |
---|
| 10 | !!---------------------------------------------------------------------- |
---|
| 11 | !! isfcav_mlt : |
---|
| 12 | !!---------------------------------------------------------------------- |
---|
| 13 | USE oce ! ocean dynamics and tracers |
---|
[11403] | 14 | USE isf ! ice shelf public variables |
---|
[11395] | 15 | USE dom_oce ! ocean space and time domain |
---|
| 16 | USE phycst ! physical constants |
---|
| 17 | USE eosbn2 ! equation of state |
---|
| 18 | ! |
---|
| 19 | USE in_out_manager ! I/O manager |
---|
| 20 | USE iom ! I/O library |
---|
| 21 | USE fldread ! read input field at current time step |
---|
[11403] | 22 | USE lib_fortran |
---|
[11395] | 23 | |
---|
| 24 | IMPLICIT NONE |
---|
| 25 | PRIVATE |
---|
| 26 | |
---|
| 27 | PUBLIC isfcav_mlt |
---|
| 28 | |
---|
| 29 | !!---------------------------------------------------------------------- |
---|
| 30 | !! NEMO/OCE 4.0 , NEMO Consortium (2018) |
---|
| 31 | !! $Id: sbcisf.F90 10536 2019-01-16 19:21:09Z mathiot $ |
---|
| 32 | !! Software governed by the CeCILL license (see ./LICENSE) |
---|
| 33 | !!---------------------------------------------------------------------- |
---|
| 34 | CONTAINS |
---|
| 35 | |
---|
| 36 | ! ------------------------------------------------------------------------------------------------------- |
---|
| 37 | ! -------------------------------- PUBLIC SUBROUTINE ---------------------------------------------------- |
---|
| 38 | ! ------------------------------------------------------------------------------------------------------- |
---|
| 39 | |
---|
[11403] | 40 | SUBROUTINE isfcav_mlt(kt, pgt, pgs , pttbl, pstbl, & |
---|
| 41 | & pqhc, pqoce, pqfwf ) |
---|
[11395] | 42 | !!---------------------------------------------------------------------- |
---|
| 43 | !! |
---|
[11403] | 44 | !! ** Purpose : compute ice shelf fwf/heqt fluxes |
---|
| 45 | !! |
---|
[11395] | 46 | !!--------------------------------------------------------------------- |
---|
| 47 | !!-------------------------- OUT ------------------------------------- |
---|
[11403] | 48 | REAL(wp), DIMENSION(jpi,jpj), INTENT( out) :: pqhc, pqoce, pqfwf ! heat and fwf fluxes |
---|
[11395] | 49 | !!-------------------------- IN ------------------------------------- |
---|
| 50 | INTEGER, INTENT(in) :: kt |
---|
| 51 | REAL(wp), DIMENSION(jpi,jpj), INTENT(in ) :: pgt , pgs ! gamma t and gamma s |
---|
| 52 | REAL(wp), DIMENSION(jpi,jpj), INTENT(in ) :: pttbl, pstbl ! top boundary layer tracer |
---|
| 53 | !!--------------------------------------------------------------------- |
---|
| 54 | !!--------------------------------------------------------------------- |
---|
| 55 | ! |
---|
| 56 | ! compute latent heat and melt (2d) |
---|
| 57 | SELECT CASE ( cn_isfcav_mlt ) |
---|
| 58 | CASE ( 'spe' ) ! ice shelf melt specified (read input file, and heat fluxes derived from |
---|
[11425] | 59 | CALL isfcav_mlt_spe( kt, pstbl, & |
---|
[11403] | 60 | & pqhc, pqoce, pqfwf ) |
---|
[11395] | 61 | CASE ( '2eq' ) ! ISOMIP formulation (2 equations) for volume flux (Hunter et al., 2006) |
---|
[11425] | 62 | CALL isfcav_mlt_2eq( pgt, pttbl, pstbl, & |
---|
[11403] | 63 | & pqhc , pqoce, pqfwf ) |
---|
[11395] | 64 | CASE ( '3eq' ) ! ISOMIP+ formulation (3 equations) for volume flux (Asay-Davis et al., 2015) |
---|
[11425] | 65 | CALL isfcav_mlt_3eq( pgt, pgs , pttbl, pstbl, & |
---|
[11403] | 66 | & pqhc, pqoce, pqfwf ) |
---|
[11395] | 67 | CASE ( 'oasis' ) ! fwf pass trough oasis |
---|
[11425] | 68 | CALL isfcav_mlt_oasis( kt, pstbl, & |
---|
[11403] | 69 | & pqhc, pqoce, pqfwf ) |
---|
[11395] | 70 | CASE DEFAULT |
---|
| 71 | CALL ctl_stop('STOP', 'unknown isf melt formulation : cn_isfcav (should not see this)') |
---|
| 72 | END SELECT |
---|
| 73 | ! |
---|
| 74 | END SUBROUTINE isfcav_mlt |
---|
| 75 | |
---|
| 76 | ! ------------------------------------------------------------------------------------------------------- |
---|
| 77 | ! -------------------------------- PRIVATE SUBROUTINE --------------------------------------------------- |
---|
| 78 | ! ------------------------------------------------------------------------------------------------------- |
---|
| 79 | |
---|
[11425] | 80 | SUBROUTINE isfcav_mlt_spe(kt, pstbl, & ! <<== in |
---|
[11403] | 81 | & pqhc , pqoce, pqfwf ) ! ==>> out |
---|
[11395] | 82 | !!---------------------------------------------------------------------- |
---|
[11403] | 83 | !! ** Purpose : - read ice shelf melt from forcing file |
---|
| 84 | !! - compute ocea-ice heat flux (assuming it is equal to latent heat) |
---|
| 85 | !! - compute heat content flux |
---|
[11395] | 86 | !!--------------------------------------------------------------------- |
---|
| 87 | !!-------------------------- OUT ------------------------------------- |
---|
[11403] | 88 | REAL(wp), DIMENSION(jpi,jpj), INTENT( out) :: pqhc, pqoce, pqfwf ! heat content, latent heat and fwf fluxes |
---|
[11395] | 89 | !!-------------------------- IN ------------------------------------- |
---|
[11403] | 90 | INTEGER , INTENT(in ) :: kt ! current time step |
---|
| 91 | REAL(wp), DIMENSION(jpi,jpj), INTENT(in ) :: pstbl ! salinity in tbl |
---|
[11395] | 92 | !!-------------------------------------------------------------------- |
---|
[11403] | 93 | REAL(wp), DIMENSION(jpi,jpj) :: ztfrz ! tbl freezing temperature |
---|
[11395] | 94 | !!-------------------------------------------------------------------- |
---|
| 95 | ! |
---|
| 96 | ! Calculate freezing temperature |
---|
| 97 | CALL eos_fzp( pstbl(:,:), ztfrz(:,:), risfdep(:,:) ) |
---|
| 98 | ! |
---|
| 99 | ! read input file |
---|
| 100 | CALL fld_read ( kt, nn_fsbc, sf_isfcav_fwf ) |
---|
| 101 | ! |
---|
| 102 | ! define fwf and qoce |
---|
| 103 | ! ocean heat flux is assume to be equal to the latent heat |
---|
[11403] | 104 | pqfwf(:,:) = - sf_isfcav_fwf(1)%fnow(:,:,1) ! fwf ( >0 out) |
---|
| 105 | pqoce(:,:) = - pqfwf(:,:) * rLfusisf ! ocean heat flux ( >0 out) |
---|
| 106 | pqhc (:,:) = pqfwf(:,:) * ztfrz(:,:) * rcp ! heat content flux ( >0 out) |
---|
[11395] | 107 | ! |
---|
[11425] | 108 | END SUBROUTINE isfcav_mlt_spe |
---|
[11395] | 109 | |
---|
[11425] | 110 | SUBROUTINE isfcav_mlt_2eq(pgt, pttbl, pstbl, & ! <<== in |
---|
[11403] | 111 | & pqhc , pqoce, pqfwf ) ! ==>> out |
---|
[11395] | 112 | !!---------------------------------------------------------------------- |
---|
[11403] | 113 | !! ** Purpose : Compute ice shelf fwf/heqt fluxes using ISOMIP formulation (Hunter et al., 2006) |
---|
[11395] | 114 | !! |
---|
[11403] | 115 | !! ** Method : The ice shelf melt latent heat is defined as being equal to the ocean/ice heat flux. |
---|
| 116 | !! From this we can derived the fwf, ocean/ice heat flux and the heat content flux as being : |
---|
| 117 | !! qfwf = Gammat * Rau0 * Cp * ( Tw - Tfrz ) / Lf |
---|
| 118 | !! qhoce = qlat |
---|
| 119 | !! qhc = qfwf * Cp * Tfrz |
---|
| 120 | !! |
---|
| 121 | !! ** Reference : Hunter, J. R.: Specification for test models of ice shelf cavities, |
---|
| 122 | !! Tech. Rep. June, Antarctic Climate & Ecosystems Cooperative Research Centre, available at: |
---|
| 123 | !! http://staff.acecrc.org.au/~bkgalton/ISOMIP/test_cavities.pdf (last access: 21 July 2016), 2006. |
---|
[11395] | 124 | !!--------------------------------------------------------------------- |
---|
| 125 | !!-------------------------- OUT ------------------------------------- |
---|
[11403] | 126 | REAL(wp), DIMENSION(jpi,jpj), INTENT( out) :: pqhc, pqoce, pqfwf ! hean content, ocean-ice heat and fwf fluxes |
---|
[11395] | 127 | !!-------------------------- IN ------------------------------------- |
---|
| 128 | REAL(wp), DIMENSION(jpi,jpj), INTENT(in ) :: pgt ! temperature exchange coeficient |
---|
| 129 | REAL(wp), DIMENSION(jpi,jpj), INTENT(in ) :: pttbl, pstbl ! temperature and salinity in top boundary layer |
---|
| 130 | !!-------------------------------------------------------------------- |
---|
| 131 | REAL(wp), DIMENSION(jpi,jpj) :: ztfrz ! freezing temperature |
---|
[11403] | 132 | REAL(wp), DIMENSION(jpi,jpj) :: zthd ! thermal driving |
---|
[11395] | 133 | !!-------------------------------------------------------------------- |
---|
| 134 | ! |
---|
| 135 | ! Calculate freezing temperature |
---|
| 136 | CALL eos_fzp( pstbl(:,:), ztfrz(:,:), risfdep(:,:) ) |
---|
| 137 | ! |
---|
[11403] | 138 | ! thermal driving |
---|
[11488] | 139 | zthd (:,:) = ( pttbl(:,:) - ztfrz(:,:) ) * mskisf_cav(:,:) |
---|
[11403] | 140 | ! |
---|
[11395] | 141 | ! compute ocean-ice heat flux and then derive fwf assuming that ocean heat flux equal latent heat |
---|
[11403] | 142 | pqfwf(:,:) = - pgt(:,:) * rau0_rcp * zthd(:,:) * r1_Lfusisf ! fresh water flux ( > 0 out ) |
---|
| 143 | pqoce(:,:) = - pqfwf(:,:) * rLfusisf ! ocea-ice flux ( > 0 out ) |
---|
| 144 | pqhc (:,:) = pqfwf(:,:) * ztfrz(:,:) * rcp ! heat content flux ( > 0 out ) |
---|
[11395] | 145 | ! |
---|
| 146 | ! output thermal driving |
---|
[11403] | 147 | CALL iom_put('isfthermald_cav', zthd ) |
---|
[11395] | 148 | ! |
---|
[11425] | 149 | END SUBROUTINE isfcav_mlt_2eq |
---|
[11395] | 150 | |
---|
[11425] | 151 | SUBROUTINE isfcav_mlt_3eq(pgt, pgs , pttbl, pstbl, & |
---|
[11403] | 152 | & pqhc, pqoce, pqfwf ) |
---|
[11395] | 153 | !!---------------------------------------------------------------------- |
---|
[11403] | 154 | !! ** Purpose : Compute ice shelf fwf/heqt fluxes using the 3 equation formulation |
---|
[11395] | 155 | !! |
---|
[11403] | 156 | !! ** Method : The melt rate is determined considering the heat balance, the salt balance |
---|
| 157 | !! at the phase change interface and a linearisation of the equation of state. |
---|
| 158 | !! |
---|
| 159 | !! ** Reference : - Holland, D. M. and Jenkins, A., |
---|
| 160 | !! Modeling Thermodynamic Ice-Ocean Interactions at the Base of an Ice Shelf, |
---|
| 161 | !! J. Phys. Oceanogr., 29, 1999. |
---|
| 162 | !! - Asay-Davis, X. S., Cornford, S. L., Durand, G., Galton-Fenzi, B. K., Gladstone, |
---|
| 163 | !! R. M., Gudmundsson, G. H., Hattermann, T., Holland, D. M., Holland, D., Holland, |
---|
| 164 | !! P. R., Martin, D. F., Mathiot, P., Pattyn, F., and Seroussi, H.: |
---|
| 165 | !! Experimental design for three interrelated marine ice sheet and ocean model intercomparison projects: |
---|
| 166 | !! MISMIP v. 3 (MISMIP +), ISOMIP v. 2 (ISOMIP +) and MISOMIP v. 1 (MISOMIP1), |
---|
| 167 | !! Geosci. Model Dev., 9, 2471-2497, https://doi.org/10.5194/gmd-9-2471-2016, 2016. |
---|
[11395] | 168 | !!--------------------------------------------------------------------- |
---|
| 169 | !!-------------------------- OUT ------------------------------------- |
---|
[11403] | 170 | REAL(wp), DIMENSION(jpi,jpj), INTENT( out) :: pqhc, pqoce, pqfwf ! latent heat and fwf fluxes |
---|
[11395] | 171 | !!-------------------------- IN ------------------------------------- |
---|
[11403] | 172 | REAL(wp), DIMENSION(jpi,jpj), INTENT(in ) :: pgt , pgs ! heat/salt exchange coeficient |
---|
| 173 | REAL(wp), DIMENSION(jpi,jpj), INTENT(in ) :: pttbl, pstbl ! mean temperature and salinity in top boundary layer |
---|
[11395] | 174 | !!-------------------------------------------------------------------- |
---|
[11403] | 175 | REAL(wp) :: zeps1,zeps2,zeps3,zeps4,zeps6,zeps7 ! dummy local scalar for quadratic equation resolution |
---|
| 176 | REAL(wp) :: zaqe,zbqe,zcqe,zaqer,zdis,zsfrz,zcfac ! dummy local scalar for quadratic equation resolution |
---|
[11395] | 177 | REAL(wp) :: zeps = 1.e-20 |
---|
[11403] | 178 | REAL(wp), DIMENSION(jpi,jpj) :: ztfrz ! freezing point |
---|
| 179 | REAL(wp), DIMENSION(jpi,jpj) :: zqcon ! conductive flux through the ice shelf |
---|
| 180 | REAL(wp), DIMENSION(jpi,jpj) :: zthd ! thermal driving |
---|
[11395] | 181 | ! |
---|
| 182 | INTEGER :: ji, jj ! dummy loop indices |
---|
| 183 | !!-------------------------------------------------------------------- |
---|
| 184 | ! |
---|
| 185 | ! compute upward heat flux zhtflx and upward water flux zwflx |
---|
[11403] | 186 | ! Resolution of a 3d equation from equation 24, 25 and 26 (note conduction through the ice has been added to Eq 24) |
---|
[11395] | 187 | DO jj = 1, jpj |
---|
| 188 | DO ji = 1, jpi |
---|
| 189 | ! |
---|
| 190 | ! compute coeficient to solve the 2nd order equation |
---|
| 191 | zeps1 = rau0_rcp*pgt(ji,jj) |
---|
| 192 | zeps2 = rLfusisf*rau0*pgs(ji,jj) |
---|
| 193 | zeps3 = rhoisf*rcpisf*rkappa/MAX(risfdep(ji,jj),zeps) |
---|
| 194 | zeps4 = risf_lamb2+risf_lamb3*risfdep(ji,jj) |
---|
| 195 | zeps6 = zeps4-pttbl(ji,jj) |
---|
| 196 | zeps7 = zeps4-rtsurf |
---|
| 197 | ! |
---|
| 198 | ! solve the 2nd order equation to find zsfrz |
---|
| 199 | zaqe = risf_lamb1 * (zeps1 + zeps3) |
---|
[11423] | 200 | zaqer = 0.5_wp/MIN(zaqe,-zeps) |
---|
[11395] | 201 | zbqe = zeps1*zeps6+zeps3*zeps7-zeps2 |
---|
| 202 | zcqe = zeps2*pstbl(ji,jj) |
---|
| 203 | zdis = zbqe*zbqe-4.0_wp*zaqe*zcqe |
---|
| 204 | ! |
---|
| 205 | ! Presumably zdis can never be negative because gammas is very small compared to gammat |
---|
| 206 | zsfrz=(-zbqe-SQRT(zdis))*zaqer |
---|
| 207 | IF ( zsfrz < 0.0_wp ) zsfrz=(-zbqe+SQRT(zdis))*zaqer ! check this if this if is needed |
---|
| 208 | ! |
---|
| 209 | ! compute t freeze (eq. 25) |
---|
| 210 | ztfrz(ji,jj) = zeps4 + risf_lamb1 * zsfrz |
---|
| 211 | ! |
---|
[11403] | 212 | ! thermal driving |
---|
[11488] | 213 | zthd(ji,jj) = ( pttbl(ji,jj) - ztfrz(ji,jj) ) |
---|
[11403] | 214 | ! |
---|
[11395] | 215 | ! compute the upward water and heat flux (eq. 24 and eq. 26) |
---|
[11403] | 216 | pqfwf(ji,jj) = rau0 * pgs(ji,jj) * ( zsfrz - pstbl(ji,jj) ) / MAX(zsfrz,zeps) ! fresh water flux (> 0 out) |
---|
| 217 | pqoce(ji,jj) = rau0_rcp * pgt(ji,jj) * zthd (ji,jj) ! ocean-ice heat flux (> 0 out) |
---|
| 218 | pqhc (ji,jj) = rcp * pqfwf(ji,jj) * ztfrz(ji,jj) ! heat content flux (> 0 out) |
---|
[11395] | 219 | ! |
---|
[11488] | 220 | zqcon(ji,jj) = zeps3 * ( ztfrz(ji,jj) - rtsurf ) |
---|
[11403] | 221 | ! |
---|
[11395] | 222 | END DO |
---|
| 223 | END DO |
---|
| 224 | ! |
---|
| 225 | ! output conductive heat flux through the ice |
---|
[11488] | 226 | CALL iom_put('qconisf', zqcon(:,:) * mskisf_cav(:,:) ) |
---|
[11395] | 227 | ! |
---|
| 228 | ! output thermal driving |
---|
[11488] | 229 | CALL iom_put('isfthermald_cav', zthd(:,:) * mskisf_cav(:,:) ) |
---|
[11395] | 230 | ! |
---|
[11425] | 231 | END SUBROUTINE isfcav_mlt_3eq |
---|
[11395] | 232 | |
---|
[11425] | 233 | SUBROUTINE isfcav_mlt_oasis(kt, pstbl, & ! <<== in |
---|
[11403] | 234 | & pqhc , pqoce, pqfwf ) ! ==>> out |
---|
| 235 | !!---------------------------------------------------------------------- |
---|
| 236 | !! |
---|
| 237 | !! ** Purpose : scale the fwf read from input file by the total amount received by the sbccpl interface |
---|
| 238 | !! |
---|
| 239 | !! ** Purpose : - read ice shelf melt from forcing file => pattern |
---|
| 240 | !! - total amount of fwf is given by sbccpl (fwfisf_cpl) |
---|
| 241 | !! - scale fwf and compute heat fluxes |
---|
| 242 | !! |
---|
| 243 | !!--------------------------------------------------------------------- |
---|
| 244 | !!-------------------------- OUT ------------------------------------- |
---|
| 245 | REAL(wp), DIMENSION(jpi,jpj), INTENT( out) :: pqhc, pqoce, pqfwf ! heat content, latent heat and fwf fluxes |
---|
| 246 | !!-------------------------- IN ------------------------------------- |
---|
| 247 | INTEGER , INTENT(in ) :: kt ! current time step |
---|
| 248 | REAL(wp), DIMENSION(jpi,jpj), INTENT(in ) :: pstbl ! salinity in tbl |
---|
| 249 | !!-------------------------------------------------------------------- |
---|
[11423] | 250 | REAL(wp) :: zfwf_fld, zfwf_oasis ! total fwf in the forcing fields (pattern) and from the oasis interface (amount) |
---|
[11403] | 251 | REAL(wp), DIMENSION(jpi,jpj) :: ztfrz ! tbl freezing temperature |
---|
| 252 | REAL(wp), DIMENSION(jpi,jpj) :: zfwf ! 2d fwf map after scaling |
---|
| 253 | !!-------------------------------------------------------------------- |
---|
| 254 | ! |
---|
| 255 | ! Calculate freezing temperature |
---|
| 256 | CALL eos_fzp( pstbl(:,:), ztfrz(:,:), risfdep(:,:) ) |
---|
| 257 | ! |
---|
| 258 | ! read input file |
---|
| 259 | CALL fld_read ( kt, nn_fsbc, sf_isfcav_fwf ) |
---|
| 260 | ! |
---|
| 261 | ! ice shelf 2d map |
---|
| 262 | zfwf(:,:) = - sf_isfcav_fwf(1)%fnow(:,:,1) |
---|
| 263 | ! |
---|
| 264 | ! compute glob sum from input file |
---|
[11425] | 265 | ! (PM) should consider delay sum as in fwb (1 time step offset if I well understood) |
---|
| 266 | zfwf_fld = glob_sum('isfcav_mlt', e1e2t(:,:) * zfwf(:,:)) |
---|
[11403] | 267 | ! |
---|
| 268 | ! compute glob sum from atm->oce ice shelf fwf |
---|
[11425] | 269 | ! (PM) should consider delay sum as in fwb (1 time step offset if I well understood) |
---|
| 270 | zfwf_oasis = glob_sum('isfcav_mlt', e1e2t(:,:) * fwfisf_oasis(:,:)) |
---|
[11403] | 271 | ! |
---|
| 272 | ! scale fwf |
---|
[11423] | 273 | zfwf(:,:) = zfwf(:,:) * zfwf_oasis / zfwf_fld |
---|
[11403] | 274 | ! |
---|
| 275 | ! define fwf and qoce |
---|
| 276 | ! ocean heat flux is assume to be equal to the latent heat |
---|
| 277 | pqfwf(:,:) = zfwf(:,:) ! fwf ( >0 out) |
---|
| 278 | pqoce(:,:) = - pqfwf(:,:) * rLfusisf ! ocean heat flux ( >0 out) |
---|
| 279 | pqhc (:,:) = pqfwf(:,:) * ztfrz(:,:) * rcp ! heat content flux ( >0 out) |
---|
| 280 | ! |
---|
[11425] | 281 | END SUBROUTINE isfcav_mlt_oasis |
---|
[11403] | 282 | |
---|
[11395] | 283 | !SUBROUTINE isfmlt_3eq_frz_ktm1 |
---|
| 284 | ! compute tfrz based on sfrz value at kt-1 (need to be SAVED local array) |
---|
[11423] | 285 | ! => should reduce error due to linarisation |
---|
[11395] | 286 | ! compute qfwf (eq 24) |
---|
| 287 | ! compute zqoce, zqlat, zqcon, zqhc |
---|
| 288 | ! compute sfrz (eq 26) |
---|
| 289 | !END SUBROUTINE isfmlt_3eq_frz_ktm1 |
---|
| 290 | |
---|
| 291 | END MODULE isfcavmlt |
---|