Bug in traadv_fct.F90

[acc]

A minor bug exists in the interpolation scheme for the FCT 4th order vertical scheme. The bug was introduced when the code was adapted to cope with the under ice shelf cavities but affects all runs using this option. It manifests as spurious temperatures outside the initial range early in runs of the OVERFLOW test case. Even with the fix shown below, spurious values still occur with this test but much later in the run (> 3000 time steps) and these are in response to excessive vertical velocities which stop the run after 3329 time steps.

Note the pre-ISF enabled version is still in the code (as interp_4th_cpt_org) but is not used. It can be removed subject to successful testing of the new version in an ISF case.

Context

The wrong index is used when setting the bottom boundary condition. The use of jk does not stray out of bounds because jk has a residual valid value from an earlier loop but the boundary condition will be wrong for anything other than flat-bottomed cases.

Analysis

The fix is simple:

src/OCE/TRA/traadv_fct.F90

@@ -555 +555 @@
             zwd (ji,jj,ikt) = 1._wp          ! top
             zwi (ji,jj,ikt) = 0._wp
             zws (ji,jj,ikt) = 0._wp
-            zwrm(ji,jj,ikt) = 0.5_wp * ( pt_in(ji,jj,jk-1) + pt_in(ji,jj,jk) )
+            zwrm(ji,jj,ikt) = 0.5_wp * ( pt_in(ji,jj,ikt-1) + pt_in(ji,jj,ikt) )
             !
             zwd (ji,jj,ikb) = 1._wp          ! bottom
             zwi (ji,jj,ikb) = 0._wp
             zws (ji,jj,ikb) = 0._wp
-            zwrm(ji,jj,ikb) = 0.5_wp * ( pt_in(ji,jj,jk-1) + pt_in(ji,jj,jk) )            
+            zwrm(ji,jj,ikb) = 0.5_wp * ( pt_in(ji,jj,ikb-1) + pt_in(ji,jj,ikb) )            
          END DO
       END DO

Recommendation

[acc]

In 9901:

Error: Failed to load processor CommitTicketReference

No macro or processor named 'CommitTicketReference' found

[acc]

Final fix also includes code to set level 2 zwd values in the ice cavity case which could be left unset for thick ice shelves and pose a risk of generating Nans. Extra code is simply:

      IF ( ln_isfcav ) THEN            ! set level two values which may not be set in ISF case
         zwd(:,:,2) = 1._wp  ;  zwi(:,:,2) = 0._wp  ;  zws(:,:,2) = 0._wp  ;  zwrm(:,:,2) = 0._wp
      END IF