`ComputeStrain_S/T` high-side `y` upwind gates use the wrong sign for `tau12` and `tau23`

[asalmgren]

Location

• Source/Diffusion/ERF_ComputeStrain_S.cpp:252-257
• Source/Diffusion/ERF_ComputeStrain_S.cpp:299-304
• Source/Diffusion/ERF_ComputeStrain_T.cpp:365-372
• Source/Diffusion/ERF_ComputeStrain_T.cpp:441-448
• Reference implementations:
  • Source/Diffusion/ERF_ComputeStrain_N.cpp:227-232
  • Source/Diffusion/ERF_ComputeStrain_N.cpp:267-272
  • Source/Diffusion/ERF_ComputeStrain_EB.cpp:215-220
  • Source/Diffusion/ERF_ComputeStrain_EB.cpp:251-256

Problem

The high-y ext_dir_upwind gates in the stretched and terrain strain kernels test for v >= 0 before using the
one-sided boundary stencil:

if (!need_to_test || v(i,dom_hi.y+1,k) >= zero) {
    tau12(i,j,k) = ...

and similarly for tau23.

Why this is a bug

At the high-y boundary, inflow corresponds to v <= 0, not v >= 0.

The normal-grid and EB kernels use v(i,dom_hi.y+1,k) <= zero for the matching high-side y stencils. The current
S/T logic therefore applies the one-sided inflow stencil on outflow and the fallback first-order stencil on inflow.

That reverses the intended ext_dir_upwind behavior for both tau12 and tau23 on the high-y boundary in
stretched and terrain-fitted runs.

Suggested patch

Use the same high-side inflow test as the N and EB kernels.

--- a/Source/Diffusion/ERF_ComputeStrain_S.cpp
+++ b/Source/Diffusion/ERF_ComputeStrain_S.cpp
@@
-            if (!need_to_test || v(i,dom_hi.y+1,k) >= zero) {
+            if (!need_to_test || v(i,dom_hi.y+1,k) <= zero) {
                 tau12(i,j,k) = myhalf * ( -(-(Real(8.)/three) * u(i,j,k) + three * u(i,j-1,k) - (one/three) * u(i,j-2,k))*dxInv[1]*mfy +
                                      + (v(i, j, k) - v(i-1, j, k))*dxInv[0]*mfx);
@@
-            if (!need_to_test || v(i,dom_hi.y+1,k) >= zero) {
+            if (!need_to_test || v(i,dom_hi.y+1,k) <= zero) {
                 tau23(i,j,k) = myhalf * ( dv_dz
                                      - (-(Real(8.)/three) * w(i,j  ,k) + three * w(i,j-1,k) - (one/three) * w(i,j-2,k))*dxInv[1] * mfy );

--- a/Source/Diffusion/ERF_ComputeStrain_T.cpp
+++ b/Source/Diffusion/ERF_ComputeStrain_T.cpp
@@
-            if (!need_to_test || v(i,dom_hi.y+1,k) >= zero) {
+            if (!need_to_test || v(i,dom_hi.y+1,k) <= zero) {
                 tau12(i,j,k) = myhalf * ( -(-(Real(8.)/three) * u(i,j,k) + three * u(i,j-1,k) - (one/three) * u(i,j-2,k))*dxInv[1]*mfy +
                                      + (v(i, j, k) - v(i-1, j, k))*dxInv[0]*mfx
@@
-            if (!need_to_test || v(i,dom_hi.y+1,k) >= zero) {
+            if (!need_to_test || v(i,dom_hi.y+1,k) <= zero) {
                 tau23(i,j,k) = myhalf * ( dv_dz
                                      - ( (-(Real(8.)/three) * w(i,j  ,k) + three * w(i,j-1,k) - (one/three) * w(i,j-2,k))*dxInv[1]
                                        - (met_h_eta)*GradWz ) * mfy );