lint: refuse an array element passed to a device routine that runs a seq loop (#1816)

Author: sbryngelson

docs/documentation/contributing.md

@@ -181,6 +181,7 @@ Both human reviewers and AI code reviewers reference this section.
 - **`stp` vs `wp` mixing:** In mixed-precision mode, `stp` (storage) may be half-precision while `wp` (working) is double. Conversions between them must be intentional, especially in MPI pack/unpack and RHS accumulation.
 - **No double-precision intrinsics:** `dsqrt`, `dexp`, `dlog`, `dble`, `dabs`, `real(8)`, `real(4)` are forbidden. Use generic intrinsics with `wp` kind.
 - **MPI type matching:** `mpi_p` must match `wp`; `mpi_io_p` must match `stp`. Mismatches corrupt communicated data.
+- **Scalars into device routines that loop:** a `GPU_ROUTINE` containing any `GPU_LOOP` (itself or through what it calls) must be called with scalars, never an array element (`q%%sf(j,k,l)`, `alpha(i)`). Copy the element to a local first and receive results into a local. Cray OpenACC 19 to 21 miscompiles the pair silently at every routine level, OpenMP offload does not ([#1815](https://github.com/MFlowCode/MFC/issues/1815)); the linter enforces it inside kernels and device routines.
 
 ### Memory and Allocation
 

toolchain/mfc/lint_source.py

@@ -554,6 +554,139 @@ def check_checker_input_constraints(repo_root: Path) -> list[str]:
     return errors
 
 
+# Intrinsics and MFC function prefixes that make `name(args)` a value, not an array element.
+_VALUE_CALL_NAMES = re.compile(r"^(?:f_\w+|real|int|nint|abs|max|min|sqrt|exp|log|sign|merge|mod|huge|tiny|epsilon|size|lbound|ubound|present|allocated|associated)$", re.IGNORECASE)
+# `a(i)`, `q(i)%sf(j, k, l)`, `s%vf(1)%sf(k, l, m)`: an element reference and nothing else.
+_ELEMENT_ARG = re.compile(r"^([A-Za-z_]\w*)(?:\([^()]*\))?(?:%\w+(?:\([^()]*\))?)*\([^()]*\)$")
+_PROCEDURE_DECL = re.compile(r"^(?:(?:impure|pure|elemental|recursive|module|non_recursive)\s+)*(?:subroutine|function)\s+(\w+)", re.IGNORECASE)
+_PROCEDURE_END = re.compile(r"^end\s+(?:subroutine|function)\b", re.IGNORECASE)
+
+
+def _split_top_level(text: str) -> list[str]:
+    """Split an argument list at the commas outside parentheses, array constructors and strings."""
+    parts, depth, quote, cur = [], 0, "", []
+    for ch in text:
+        if quote:
+            if ch == quote:
+                quote = ""
+        elif ch in "'\"":
+            quote = ch
+        elif ch in "([":
+            depth += 1
+        elif ch in ")]":
+            depth -= 1
+        elif ch == "," and depth == 0:
+            parts.append("".join(cur).strip())
+            cur = []
+            continue
+        cur.append(ch)
+    parts.append("".join(cur).strip())
+    return parts
+
+
+def _statements(lines: list[str]):
+    """Yield (first line number, statement) with Fortran continuation lines joined."""
+    buf, start = [], 0
+    for i, line in enumerate(lines, 1):
+        stripped = line.split("!", 1)[0].strip() if not line.strip().startswith("!") else ""
+        if not stripped:
+            continue
+        if not buf:
+            start = i
+        buf.append(stripped.lstrip("&").rstrip("&").strip())
+        if not stripped.endswith("&"):
+            yield start, " ".join(buf)
+            buf = []
+
+
+def _procedures(lines: list[str]):
+    """Yield (name, own line numbers, is device routine, has seq loop) per procedure.
+
+    Contained procedures nest; every line, directive and call belongs to the innermost one.
+    """
+    stack = []  # [name, own lines, device, looped]
+    for i, line in enumerate(lines, 1):
+        stripped = line.strip()
+        m = _PROCEDURE_DECL.match(stripped)
+        if m and not stripped.lower().startswith("end"):
+            stack.append([m.group(1), set(), False, False])
+        if not stack:
+            continue
+        stack[-1][1].add(i)
+        if "GPU_ROUTINE(" in stripped:
+            stack[-1][2] = True
+        elif "GPU_LOOP(" in stripped:
+            stack[-1][3] = True
+        elif _PROCEDURE_END.match(stripped):
+            name, own, device, looped = stack.pop()
+            yield name, own, device, looped
+
+
+def check_device_routine_element_args(repo_root: Path) -> list[str]:
+    """Flag an array element passed to a device routine that runs a seq loop.
+
+    CCE OpenACC (19.0.0 through 21.0.2, -O2; OpenMP offload is correct) miscompiles the pair: a
+    routine containing any `GPU_LOOP`, called with an array element as an actual argument, reads
+    the element as garbage and never writes it back. Either alone is fine, every `routine` level
+    is affected, and the loop counts when it sits in anything the routine calls. Copy the element
+    to a scalar before the call and receive results into a scalar. See
+    .claude/rules/common-pitfalls.md and sbryngelson/compiler-bugs cce/acc-routine-element-by-reference.
+    """
+    src_dir = repo_root / SRC_DIR
+    files = {src: src.read_text(encoding="utf-8").splitlines() for src in _fortran_fpp_files(src_dir)}
+
+    device, looped, bodies = set(), set(), {}
+    for lines in files.values():
+        for name, own, is_device, has_loop in _procedures(lines):
+            key = name.lower()
+            if is_device:
+                device.add(key)
+            if has_loop:
+                looped.add(key)
+            bodies[key] = [stmt for n, stmt in _statements(lines) if n in own]
+    if not device:
+        return []
+    # Any reference to a device routine by name, `call s_x(` or `y = f_x(`, is a call site.
+    site_re = re.compile(r"(?<![\w%])(" + "|".join(map(re.escape, sorted(device))) + r")\s*\(", re.IGNORECASE)
+    callees = {k: {m.group(1).lower() for stmt in v for m in site_re.finditer(stmt)} - {k} for k, v in bodies.items()}
+    # A device routine that calls a looped routine carries the loop once CCE inlines it.
+    tainted = looped & device
+    while True:
+        more = {r for r in device - tainted if callees.get(r, set()) & tainted}
+        if not more:
+            break
+        tainted |= more
+
+    errors: list[str] = []
+    for src, lines in files.items():
+        rel = src.relative_to(repo_root)
+        device_lines = set()
+        for name, own, is_device, _ in _procedures(lines):
+            if is_device:
+                device_lines |= own
+        kernel_depth = 0
+        for line_no, stmt in _statements(lines):
+            if "END_GPU_PARALLEL_LOOP" in stmt:
+                kernel_depth = max(0, kernel_depth - 1)
+            elif "GPU_PARALLEL_LOOP(" in stmt:
+                kernel_depth += 1
+            if kernel_depth == 0 and line_no not in device_lines:
+                continue  # host code passes elements freely
+            for m in site_re.finditer(stmt):
+                name = m.group(1).lower()
+                if name not in tainted:
+                    continue
+                depth, j = 1, m.end()
+                while j < len(stmt) and depth:
+                    depth += {"(": 1, ")": -1}.get(stmt[j], 0)
+                    j += 1
+                for arg in _split_top_level(stmt[m.end() : j - 1]):
+                    e = _ELEMENT_ARG.match(arg)
+                    if e and ":" not in arg and not _VALUE_CALL_NAMES.match(e.group(1)):
+                        errors.append(f"  {rel}:{line_no} `{arg}` into `{name}` (a device routine with a seq loop): pass a scalar, see common-pitfalls.md")
+    return errors
+
+
 def check_cluster_menu_slugs(repo_root: Path) -> list[str]:
     """Keep the ``./mfc.sh load`` cluster menu in sync with toolchain/modules.
 
@@ -614,6 +747,7 @@ def main():
     all_errors.extend(check_manual_registry_bcasts(repo_root))
     all_errors.extend(check_checker_input_constraints(repo_root))
     all_errors.extend(check_cluster_menu_slugs(repo_root))
+    all_errors.extend(check_device_routine_element_args(repo_root))
 
     if all_errors:
         print("Source lint failed:")

toolchain/mfc/test_lint_source.py

@@ -2,6 +2,7 @@
 
 from mfc.lint_source import (
     _extract_bcast_roots,
+    check_device_routine_element_args,
     check_double_precision,
     check_integer_wp,
     check_manual_registry_bcasts,
@@ -148,3 +149,123 @@ def test_manual_residue_is_clean(tmp_path):
     _write_proxy(tmp_path, "simulation", body)
 
     assert check_manual_registry_bcasts(tmp_path) == []
+
+
+_LOOPED = """    subroutine s_curve(rho, i, p)
+        $:GPU_ROUTINE(parallelism='[seq]')
+        real(wp), intent(in) :: rho
+        integer, intent(in) :: i
+        real(wp), intent(out) :: p
+        integer :: it
+        $:GPU_LOOP(parallelism='[seq]')
+        do it = 1, 8
+            p = p + rho
+        end do
+    end subroutine s_curve
+    subroutine s_wrap(rho, i, p)
+        $:GPU_ROUTINE(parallelism='[seq]')
+        real(wp), intent(in) :: rho
+        integer, intent(in) :: i
+        real(wp), intent(out) :: p
+        call s_curve(rho, i, p)
+    end subroutine s_wrap
+    subroutine s_plain(rho, i, p)
+        $:GPU_ROUTINE(parallelism='[seq]')
+        real(wp), intent(in) :: rho
+        integer, intent(in) :: i
+        real(wp), intent(out) :: p
+        p = rho
+    end subroutine s_plain
+"""
+
+
+def _KERNEL(body: str) -> str:
+    return "        $:GPU_PARALLEL_LOOP(collapse=3)\n        " + body + "\n        $:END_GPU_PARALLEL_LOOP()\n"
+
+
+def test_host_call_sites_are_not_flagged(tmp_path):
+    _write_src(tmp_path, "simulation/m_x.fpp", _LOOPED + "        call s_curve(q(1)%sf(j, k, l), 1, out(k, l, q))\n")
+    assert check_device_routine_element_args(tmp_path) == []
+
+
+def test_element_into_looped_device_routine_is_flagged(tmp_path):
+    _write_src(tmp_path, "simulation/m_x.fpp", _LOOPED + _KERNEL("call s_curve(q(1)%sf(j, k, l), 1, out(k, l, q))"))
+    errors = check_device_routine_element_args(tmp_path)
+    assert len(errors) == 2
+    assert "q(1)%sf(j, k, l)" in errors[0] and "out(k, l, q)" in errors[1]
+
+
+def test_element_reaches_the_loop_through_a_caller(tmp_path):
+    _write_src(tmp_path, "simulation/m_x.fpp", _LOOPED + _KERNEL("call s_wrap(pres, 1, blkmod(k, &\n            & l, q))"))
+    errors = check_device_routine_element_args(tmp_path)
+    assert len(errors) == 1 and "s_wrap" in errors[0]
+
+
+def test_scalars_expressions_and_loopless_routines_pass(tmp_path):
+    _write_src(
+        tmp_path,
+        "simulation/m_x.fpp",
+        _LOOPED
+        + _KERNEL("call s_curve(alpha_rho(i)/max(alpha(i), sgm_eps), i, p_i)\n        call s_plain(q(1)%sf(j, k, l), 1, out(k, l, q))\n        call s_curve(real(q(1)%sf(j, k, l), wp), 1, p_i)"),
+    )
+    assert check_device_routine_element_args(tmp_path) == []
+
+
+def test_loop_inside_a_device_function_counts_and_propagates(tmp_path):
+    src = """    function f_looped(x, i) result(y)
+        $:GPU_ROUTINE(function_name='f_looped', parallelism='[seq]')
+        real(wp), intent(in) :: x
+        integer, intent(in) :: i
+        real(wp) :: y
+        integer :: it
+        y = x
+        $:GPU_LOOP(parallelism='[seq]')
+        do it = 1, 8
+            y = y + 1._wp
+        end do
+    end function f_looped
+    subroutine s_via_function(x, i, y)
+        $:GPU_ROUTINE(parallelism='[seq]')
+        real(wp), intent(in) :: x
+        integer, intent(in) :: i
+        real(wp), intent(out) :: y
+        y = f_looped(x, i)
+    end subroutine s_via_function
+"""
+    calls = "out(k, l, q) = f_looped(q(1)%sf(k, l, q), 1)\n        call s_via_function(q(1)%sf(k, l, q), 1, tmp)\n        tmp = f_looped(p_scalar, 1)"
+    _write_src(tmp_path, "simulation/m_x.fpp", src + _KERNEL(calls))
+    errors = check_device_routine_element_args(tmp_path)
+    assert [e.split("`")[3] for e in errors] == ["f_looped", "s_via_function"]
+
+
+def test_constructor_commas_and_unprefixed_functions_and_contained_scoping(tmp_path):
+    src = """    function g_looped(x) result(y)
+        $:GPU_ROUTINE(function_name='g_looped', parallelism='[seq]')
+        real(wp), intent(in) :: x
+        real(wp) :: y
+        integer :: it
+        y = x
+        $:GPU_LOOP(parallelism='[seq]')
+        do it = 1, 8
+            y = y + 1._wp
+        end do
+    end function g_looped
+    subroutine s_outer(a, b)
+        real(wp), intent(in) :: a
+        real(wp), intent(out) :: b
+        b = a
+    contains
+        subroutine s_inner(x, y)
+            $:GPU_ROUTINE(parallelism='[seq]')
+            real(wp), intent(in) :: x
+            real(wp), intent(out) :: y
+            y = g_looped(x)
+        end subroutine s_inner
+    end subroutine s_outer
+"""
+    calls = "b = g_looped(q(1)%sf(k, l, q))\\n        c = g_looped(sum([v(1), v(2)]))\\n        call s_outer(q(1)%sf(k, l, q), tmp)"
+    _write_src(tmp_path, "simulation/m_x.fpp", src + _KERNEL(calls))
+    errors = check_device_routine_element_args(tmp_path)
+    # the unprefixed function is found by name; the array constructor is not split into a fake element;
+    # s_outer is not tainted by the loop that only its contained s_inner reaches (and is not a device routine)
+    assert [e.split("`")[1] for e in errors] == ["q(1)%sf(k, l, q)"]