Bug in the slope limitation in ISF case (v3.6_stable & trunk)

[gm]

It seems to me that, in ldfslp.F90, the tapering of the iso-neutral slope in vicinity of ocean surface is not properly done below the ice shelves.
This concerns both v3.6_stable and trunk version.
I give below a correction for the trunk version, that should be copy in the v3.6_stable, if ISF guys agree. I will be happy to further discuss with you the reason of this change.

in ldfslp.F90, lines 186-189

zslpml_hmlpu(ji,jj) = uslpml(ji,jj) / ( MAX(hmlpt(ji,jj), hmlpt(ji+1,jj  ), 5._wp)       &
   &                                   - 0.5_wp * ( risfdep(ji,jj) + risfdep(ji+1,jj  ) ) )
zslpml_hmlpv(ji,jj) = vslpml(ji,jj) / ( MAX(hmlpt(ji,jj), hmlpt(ji  ,jj+1), 5._wp)       &
   &                                   - 0.5_wp * ( risfdep(ji,jj) + risfdep(ji  ,jj+1) ) )

should be replaced by:

zslp_mldu(ji,jj) = uslpml(ji,jj) / (  MAX( hmlpt  (ji,jj) , hmlpt  (ji+1,jj), 5._wp)  &
   &                                - MAX( risfdep(ji,jj) , risfdep(ji+1,jj)       )  )
zslp_mldv(ji,jj) = vslpml(ji,jj) / (  MAX( hmlpt  (ji,jj) , hmlpt  (ji,jj+1), 5._wp)  &
   &                                - MAX( risfdep(ji,jj) , risfdep(ji,jj+1)       )  )

as the MLD at u- and v-points should be the max of the 2 surrounding t-points

and the lines 217-220:

zdepu = 0.5_wp * ( ( gdept_n(ji,jj,jk) + gdept_n(ji+1,jj  ,jk) )                            &
   &              - ( risfdep(ji,jj)    + risfdep(ji+1,jj)    ) - e3u_n(ji,jj,miku(ji,jj)) )
zdepv = 0.5_wp * ( ( gdept_n(ji,jj,jk) + gdept_n(ji,jj+1,jk) )                              &
   &              - ( risfdep(ji,jj)    + risfdep(ji,jj+1)    ) - e3v_n(ji,jj,mikv(ji,jj)) )

should be replaced by:

zdepu = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji+1,jj,jk) ) - MAX( risfdep(ji,jj) , risfdep(ji+1,jj) )
zdepv = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji,jj+1,jk) ) - MAX( risfdep(ji,jj) , risfdep(ji,jj+1) )

[mathiot]

After discussion, we agreed with the modifications concerning zslpml_hmlpu. About zdepu, we agreed about the MAX(risfdep(ji,jj),risfdep(ji+1,jj)) but we agreed to keep the -0.5*e3u(ji,jj,miku(ji,jj)) to set slope to 0 at the surface level.

I corrected this for the trunk and the 3.6_STABLE respectively in revision 6352 and 6353.

To implement this on the 3.6_STABLE, I changed the definition of hmlp (and hmlp and hmld for consistency) to the depth of the bottom of the mixed instead of the thickness of the mixed layer (no change if no ice shelf cavity because ice shelf draft = 0). It is the same choice made in the trunk. However for output, mldr10_1 and mldkz5 output are respectively defined as the differences between the mixed layer/turbocline depth and the ice shelf draft (again no change if no ice shelf cavity)

This change fixes, by the way, a bug in dyndmp, trcdmp and tradmp with nn_zdmp = 2 (ie no damping into the mixed layer).

A test on compatibility between ln_isfcav and ln_traldf_grif is also added in traldf.F90 for the 3.6 and the trunk because it has not been tested.

SETTE tests OK and solver.stat + tracer.stat are the same on MetO machine for reference configurations. No change in ISOMIP because ISOMIP is not using iso neutral diffusion. Changes are expected in configuration with ice shelf cavity (ln_isfcav=true) and iso neutral diffusion.