bug in ggl90/IDEMIX code?

[mjlosch]

In pkg/ggl90/ggl90_calc.F, recip_hFacI (inverse cell thickness of w-cells) is used a few times to scale variables and coefficients, but only when useIDEMIX=.TRUE.. Two places, where I think this is wrong are:

MITgcm/pkg/ggl90/ggl90_calc.F

Lines 696 to 697 in 1575dc0

	verticalShear(i,j) = verticalShear(i,j)
	& * recip_hFacI(i,j,k)*recip_hFacI(i,j,k)

because verticalShear is also used to compute the Ri-number, where sigmaR has been computed without taking hFacs into account (grad_sigma.F). So in my view it would be consistent to not scale verticalShear with recip_hFacI**2
and

MITgcm/pkg/ggl90/ggl90_calc.F

Lines 831 to 842 in 1575dc0

	#ifdef ALLOW_GGL90_IDEMIX
	IF ( useIDEMIX ) THEN
	DO k=2,Nr
	DO j=jMin,jMax
	DO i=iMin,iMax
	a3d(i,j,k) = a3d(i,j,k)*recip_hFacI(i,j,k)
	c3d(i,j,k) = c3d(i,j,k)*recip_hFacI(i,j,k)
	ENDDO
	ENDDO
	ENDDO
	ENDIF
	#endif /* ALLOW_GGL90_IDEMIX */

The implicit solver code was initially more or less a copy of impldiff.F, so I think that the GGL90 code is consistent with how we do things in general (i.e. only use recip_drF*recip_hFacC but not recip_hFacI with recip_drC).

As a consequence verticalShear can be too large near the bottom leading to too much forcing, and similarly for the coefficients a3d and c3d. Any opinions are welcome.

[jm-c]

@mjlosch Did the changes in PR #714 (also dealing with some recip_hFacI) address this issue above ?

[mjlosch]

Yes, but in #714, we took a different stance and agreed that using recip_hFacI in the coefficients a3d and c3d is correct and now we do it for all cases and in #714 we did not address the first part with the verticalShear, and I am not longer certain, that the current code is wrong.