tra_adv_qck may give incorrect results on the northfold

[hadcv]

NOTE: This bug fix changes results when using nn_hls = 1

Context

tra_adv_qck gives different run.stat results for nn_hls = 1 compared to nn_hls = 2.

Analysis

In ​tra_adv_qck_j, the calculation of zfu on the northfold does not appear to be consistent with the rest of the domain.

231            DO jk = 1, jpkm1
232               !
233               !--- Computation of the ustream and downstream value of the tracer and the mask
234               DO_2D( 0, 0, nn_hls-1, nn_hls-1 )
235                  ! Upstream in the x-direction for the tracer
236                  zfc(ji,jj,jk) = pt(ji,jj-1,jk,jn,Kbb)
237                  ! Downstream in the x-direction for the tracer
238                  zfd(ji,jj,jk) = pt(ji,jj+1,jk,jn,Kbb)
239               END_2D
240            END DO
241            IF (nn_hls.EQ.1) CALL lbc_lnk( 'traadv_qck', zfc(:,:,:), 'T', 1.0_wp , zfd(:,:,:), 'T', 1.0_wp )   ! Lateral boundary conditions
242   
243            !
244            ! Horizontal advective fluxes
245            ! ---------------------------
246            !
247            DO_3D( 0, 0, nn_hls-1, 0, 1, jpkm1 )
248               zdir = 0.5 + SIGN( 0.5_wp, pV(ji,jj,jk) )   ! if pU > 0 : zdir = 1 otherwise zdir = 0
249               zfu(ji,jj,jk) = zdir * zfc(ji,jj,jk ) + ( 1. - zdir ) * zfd(ji,jj+1,jk)  ! FU in the x-direction for T
250            END_3D

Substituting zfc and zfd into zfu gives:

zfu(ji,jj,jk) = zdir * pt(ji,jj-1,jk,jn,Kbb) + ( 1. - zdir ) * pt(ji,jj+2,jk,jn,Kbb)

The issue is that the lbc_lnk in the nn_hls = 1 case will mirror zfd, so that:

zfd(:,Nje0,  :) = pt(:,Nje0+1,:)
zfd(:,Nje0+1,:) = pt(:,Nje0,  :)

and therefore the calculation of zfu on the northfold effectively uses a different stencil to the rest of the domain:

zfu(:,Nje0-2,:) = zdir * pt(:,Nje0-3,:,jn,Kbb) + ( 1. - zdir ) * pt(:,Nje0,  :,jn,Kbb)
zfu(:,Nje0-1,:) = zdir * pt(:,Nje0-2,:,jn,Kbb) + ( 1. - zdir ) * pt(:,Nje0+1,:,jn,Kbb)
zfu(:,Nje0,  :) = zdir * pt(:,Nje0-1,:,jn,Kbb) + ( 1. - zdir ) * pt(:,Nje0,  :,jn,Kbb)

In the nn_hls = 2 case, zfd is calculated explicitly and the northfold stencil for zfu is consistent with the rest of the domain, therefore giving different results to the nn_hls = 1 case:

zfd(:,Nje0+1,:) = pt(:,Nje0+2,:)

zfu(:,Nje0-2,:) = zdir * pt(:,Nje0-3,:,jn,Kbb) + ( 1. - zdir ) * pt(:,Nje0,  :,jn,Kbb)
zfu(:,Nje0-1,:) = zdir * pt(:,Nje0-2,:,jn,Kbb) + ( 1. - zdir ) * pt(:,Nje0+1,:,jn,Kbb)
zfu(:,Nje0,  :) = zdir * pt(:,Nje0-1,:,jn,Kbb) + ( 1. - zdir ) * pt(:,Nje0+2,:,jn,Kbb)

Fix

The correct values for zfd on the northfold can be recovered from zfc after the lbc_lnk:

+         IF( nn_hls == 1 .AND. l_IdoNFold .AND. ntej == Nje0 ) THEN
+            DO jk = 1, jpkm1
+               WHERE( tmask_i(ntsi:ntei,ntej:jpj) == 0._wp ) zfd(ntsi:ntei,ntej:jpj,jk) = zfc(ntsi:ntei,ntej:jpj,jk)
+            END DO
+         ENDIF
         !
         ! Horizontal advective fluxes
         ! ---------------------------
         !
         DO_3D( 0, 0, nn_hls-1, 0, 1, jpkm1 )
            zdir = 0.5 + SIGN( 0.5_wp, pV(ji,jj,jk) )   ! if pU > 0 : zdir = 1 otherwise zdir = 0
            zfu(ji,jj,jk) = zdir * zfc(ji,jj,jk ) + ( 1. - zdir ) * zfd(ji,jj+1,jk)  ! FU in the x-direction for T
         END_3D

This makes the results consistent for both nn_hls, but it would be useful to get an independent confirmation that this is indeed a bug.

[smasson]

For T points I would say that

jj-1 -> jj+1
jj+1 -> jj-1

jj begin the position of the same point on each side of the folding, its jj value may not be the same on each side of the folding.
So in the 2 lines above the jj on each side of the -> may not be the same.
See the attached pdf and look at the colors of the points "above" and "below" a T point.

[smasson]

after our discussion and some more drawings on a piece of paper, I fully agree with your diagnostic and your correction that is btw really smart! It works also for F-folding

[hadcv]

In SETTE, this fix causes compilation errors in AGRIF_DEMO when using loop fusion:

/home/d00/hadcv/cylc-run/u-bs939_tiling_sette/work/19760101T0000Z/nemo_sette-lf_hls_2/nemo/cfgs/AGRIF_DEMO_ST/BLD/ppsrc/nemo/traadv_qck_lf.f90(445): error #6633: The type of the actual argument differs from the type of the dummy argument.   [IARRAY0]
     &array2, Agrif_tabvars(10)%array2, Agrif_tabvars_i(16)%iarray0, Agr&
------------------------------------------------------------^
/home/d00/hadcv/cylc-run/u-bs939_tiling_sette/work/19760101T0000Z/nemo_sette-lf_hls_2/nemo/cfgs/AGRIF_DEMO_ST/BLD/ppsrc/nemo/traadv_qck_lf.f90(446): error #6633: The type of the actual argument differs from the type of the dummy argument.   [IARRAY0]
     &if_tabvars_i(25)%iarray0, Agrif_tabvars_i(23)%iarray0, Agrif_tabva&
----------------------------------------------------^
/home/d00/hadcv/cylc-run/u-bs939_tiling_sette/work/19760101T0000Z/nemo_sette-lf_hls_2/nemo/cfgs/AGRIF_DEMO_ST/BLD/ppsrc/nemo/traadv_qck_lf.f90(440): error #6631: A non-optional actual argument must be present when invoking a procedure with an explicit interface.   [JPK]
        call Sub_Loop_tra_adv_qck_j(kt,cdtype,p2dt,pV,Kbb,Kmm,pt,kjpt,Kr&
-------------^
/home/d00/hadcv/cylc-run/u-bs939_tiling_sette/work/19760101T0000Z/nemo_sette-lf_hls_2/nemo/cfgs/AGRIF_DEMO_ST/BLD/ppsrc/nemo/traadv_qck_lf.f90(440): error #6631: A non-optional actual argument must be present when invoking a procedure with an explicit interface.   [JPJ]
        call Sub_Loop_tra_adv_qck_j(kt,cdtype,p2dt,pV,Kbb,Kmm,pt,kjpt,Kr&
-------------^
/home/d00/hadcv/cylc-run/u-bs939_tiling_sette/work/19760101T0000Z/nemo_sette-lf_hls_2/nemo/cfgs/AGRIF_DEMO_ST/BLD/ppsrc/nemo/traadv_qck_lf.f90(440): error #6631: A non-optional actual argument must be present when invoking a procedure with an explicit interface.   [JPI]
        call Sub_Loop_tra_adv_qck_j(kt,cdtype,p2dt,pV,Kbb,Kmm,pt,kjpt,Kr&
-------------^
/home/d00/hadcv/cylc-run/u-bs939_tiling_sette/work/19760101T0000Z/nemo_sette-lf_hls_2/nemo/cfgs/AGRIF_DEMO_ST/BLD/ppsrc/nemo/traadv_qck_lf.f90(445): error #7976: An allocatable dummy argument may only be argument associated with an allocatable actual argument.   [IARRAY0]
     &array2, Agrif_tabvars(10)%array2, Agrif_tabvars_i(16)%iarray0, Agr&
------------------------------------------------------------^
/home/d00/hadcv/cylc-run/u-bs939_tiling_sette/work/19760101T0000Z/nemo_sette-lf_hls_2/nemo/cfgs/AGRIF_DEMO_ST/BLD/ppsrc/nemo/traadv_qck_lf.f90(445): error #6634: The shape matching rules of actual arguments and dummy arguments have been violated.   [IARRAY0]
     &array2, Agrif_tabvars(10)%array2, Agrif_tabvars_i(16)%iarray0, Agr&
------------------------------------------------------------^
compilation aborted for /home/d00/hadcv/cylc-run/u-bs939_tiling_sette/work/19760101T0000Z/nemo_sette-lf_hls_2/nemo/cfgs/AGRIF_DEMO_ST/BLD/ppsrc/nemo/traadv_qck_lf.f90 (code 1)
fcm_internal compile failed (256)
gmake: *** [traadv_qck_lf.o] Error 1

conv seems to be getting confused. In the ppsrc file for traadv_qck_lf.F90, the tra_adv_qck* functions are being taken from traadv_qck.F90.

Renaming the functions to e.g. tra_adv_qck_i_lf prevents this from happening; they are taken from traadv_qck_lf.F90 correctly.

[hadcv]

In 14978:

Error: Failed to load processor CommitTicketReference

No macro or processor named 'CommitTicketReference' found