misplaced parentheses in sfc_src usage in tracer_vertdiff, tracer_vertdiff_Eulerian

[klindsay28]

The variable sfc_src in subroutine tracer_vertdiff is declared as
https://github.com/NOAA-GFDL/MOM6/blob/c432872bebda98b0acaff7d2378dc5a8523d3377/src/tracer/MOM_tracer_diabatic.F90#L54-L56

It is subsequently used in the block
https://github.com/NOAA-GFDL/MOM6/blob/c432872bebda98b0acaff7d2378dc5a8523d3377/src/tracer/MOM_tracer_diabatic.F90#L150-L156
and similar blocks further down.

Its usage would be dimensionally consistent if its units were CU, but it is documented to have units CU H. So somethings seems amiss. It looks like the parentheses on these usage lines changed in commit 982c373. The git diff for this block is

@@ -156,7 +156,7 @@ subroutine tracer_vertdiff(h_old, ea, eb, dt, tr, G, GV, &
         b1(i) = 1.0 / (b_denom_1 + eb(i,j,1))
         d1(i) = b_denom_1 * b1(i)
         h_tr = h_work(i,j,1) + h_neglect
-        tr(i,j,1) = b1(i)*(h_tr*tr(i,j,1) + sfc_src(i,j))
+        tr(i,j,1) = (b1(i)*h_tr)*tr(i,j,1) + sfc_src(i,j)
       endif ; enddo
       do k=2,nz-1 ; do i=is,ie ; if (G%mask2dT(i,j) > 0.5) then
         c1(i,k) = eb(i,j,k-1) * b1(i)

@ashao, the author of this commit, agrees that this change in parentheses is problematic.

This change in parentheses occurred in 2 blocks in tracer_vertdiff and appears to have been duplicated into tracer_vertdiff_Eulerian with commit 7598390.

[klindsay28]

Some followup...

The effect of the change of parentheses in 982c373 was nullified, when the argument evap_CFL_limit was present, because sfc_src was set to zero after a call to applyTracerBoundaryFluxesInOut in the following block that occurs before the line with the parentheses change:
https://github.com/NOAA-GFDL/MOM6/blob/982c3733be0f9e335cbb6fa2fbc6173fa3cfb36b/src/tracer/MOM_tracer_diabatic.F90#L111-L119
However, the this call to applyTracerBoundaryFluxesInOut was migrated to tracer module's column_physics subroutines in a subsequent commit, f38bd7b. It doesn't look like the effect of zeroing out sfc_src was preserved when the applyTracerBoundaryFluxesInOut call was moved in this latter commit.

Related to the setting of sfc_src to zero, there is a comment with applyTracerBoundaryFluxesInOut that states

!! NOTE: Please note that in this routine sfc_flux gets set to zero to ensure that the surface
!! flux of the tracer does not get applied again during a subsequent call to tracer_vertdif

This comment seems to be out of date since commit dfd22a2, when the argument sfc_flux of applyTracerBoundaryFluxesInOut was replaced with sfc_src, and the the setting of sfc_flux to zero was migrated to zeroing sfc_src after the call to applyTracerBoundaryFluxesInOut.

[gustavo-marques]

We have a fix for this and would like to know if the bug needs to be preserved for the purposeses of recovering old answers. Answers will change for experiments with USE_OCMIP2_CFC = True (e.g., this configuration in MOM6-examples)

[Hallberg-NOAA]

I think that in this case, we need not reproduce the old answers, for several reasons:

1. These answers are very demonstrably wrong, in that if they grossly violate conservation of the tracer amounts and clearly fail the dimensional consistency tests. In other words, there are effectively no plausible answers to reproduce.
2. Because these corrections only apply to passive tracers, they will be easy enough to trace the changes, should anyone ever need to reproduce the old CFC calculations.
3. To the best of my knowledge, these CFCs are not documented in any publications. The closely related ESM2G configurations predated the introduction of these bugs, and all of CFCs in the CM4/ESM4 family of models used the CFCs from the generic tracer code, and not this package.

[ashao]

@Hallberg-NOAA talked about this yesterday and we both agree on the above points. There's some work that I will do to address some of the more esoteric changes that might need to be done for the pseudosalt tracer. Could you ping me as a reviewer on your PR to dev/ncar?

[gustavo-marques]

@Hallberg-NOAA and @ashao, thank you. Yes, we will ping you to review this PR.

[Hallberg-NOAA]

This correction to this problem has been merged into the main branch of MOM6, and into the dev/gfdl branch as a part of MOM6 PR #1491.