dom_vvl_2 with vvl and s-coordinates

[jchanut]

This bug may concern people using both key_vvl and sigma coordinate. The problem is in dom_vvl_2 routine (module domvvl.F90) that sets "before" vertical scale factors at velocity points from time filtered vertical scale factors at T-points. The surface weighted update between adjacent T-points conflict with bathymetry definition at velocity points (unweighted calculation in that case). As a result, the vertical integral of fse3u_b/fse3v_b does not match hu_0+sshu_b / hv_0+sshv_b.
This problem can be detected in AMM12 configuration, but in the meridional direction only (the grid is regular otherwise). The (time invariant) maximum error on the actual bathymetry used at before time step is less than 1 meter in that case. Certainly not a big deal, except if like me, you try to find out why the total depth you get is not what you should have...

As a remedy I suggest whether:

• to initialize scale factors at velocity points (in domzgr.F90) thanks to surface weighted averages.

• or change what is done in dom_vvl_2 in the form: fse3u_b(:,:,:) = fse3u_0(:,:,:)*(1+sshu_b(:,:)*muu(:,:,:))

where sshu_b holds the right sea level.

Jérôme

[jchanut]

To fix this, I finally suggest to work on thickness "anomalies" in dom_vvl_2 such that lines below:

DO jk = 1, jpkm1       ! set the before scale factors at u- & v-points
  DO jj = 2, jpjm1
    DO ji = fs_2, fs_jpim1
      zvt = fse3t_b(ji,jj,jk) * e1e2t(ji,jj)
      pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1e2t(ji+1,jj) ) / ( e1u(ji,jj) * e2u(ji,jj) )
      pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e1e2t(ji,jj+1) ) / ( e1v(ji,jj) * e2v(ji,jj) )
    END DO
  END DO
END DO

could be replaced by:

DO jk = 1, jpkm1           ! set the before scale factors at u- & v-points
  DO jj = 2, jpjm1
    DO ji = fs_2, fs_jpim1
      zvt    = ( fse3t_b(ji  ,jj  ,jk) - fse3t_0(ji  ,jj  ,jk) ) * e1e2t(ji  ,jj  )
      zvtip1 = ( fse3t_b(ji+1,jj  ,jk) - fse3t_0(ji+1,jj  ,jk) ) * e1e2t(ji+1,jj  )
      zvtjp1 = ( fse3t_b(ji  ,jj+1,jk) - fse3t_0(ji  ,jj+1,jk) ) * e1e2t(ji  ,jj+1)
      pe3u_b(ji,jj,jk) = fse3u_0(ji,jj,jk) + 0.5_wp * ( zvt + zvtip1 ) / ( e1u(ji,jj) * e2u(ji,jj) )
      pe3v_b(ji,jj,jk) = fse3v_0(ji,jj,jk) + 0.5_wp * ( zvt + zvtjp1 ) / ( e1v(ji,jj) * e2v(ji,jj) )
    END DO
  END DO
END DO

[detraceyb]

This revision also included a fix to a non-minor error on line 204.