Feat: N-15 step 2 parallel i64 reductions (#277)

docs/SPEC.md

@@ -892,7 +892,19 @@ All in `PM/commonTest`, so they run wherever the interpreter loads.
     `min(a, b)` compares with `<` and keeps the first element in iteration order on ties (including
     `+0.0` versus `-0.0`), NaN behaviour depends on position because comparisons with NaN are false,
     and chunked combination must reproduce that chain bit for bit; if an implementation cannot, they
-    are left to a later step. **Float `+` and `*` reductions are excluded in v1**: chunked evaluation
+    are left to a later step. Step 2 implements `+` and `*` only: the compiled subset has no bitwise
+    operator (the front end skips a function using `&`, `|` or `^` on ints, there is no IR operation),
+    so `& | ^` are left until the subset has them.
+    A reduction loop is `acc = acc op e` statements only (no element-wise stores in the same loop, so
+    a rerun needs no restore), each `acc` a function-local i64 read nowhere else in the loop and
+    assigned before it, `e` built from i64 constants, locals, `b[i]` of a list or local int array
+    and i64 `+ - *`. Each chunk computes its partial from the identity with the sequential code's
+    overflow checks; an overflow in any chunk, in an element expression, or in the combination
+    discards the region (no accumulator written) and the sequential loop runs from its start state,
+    which deopts or promotes exactly as it does today, so the result is Python's unbounded int. The
+    exact integer result of a region that does not overflow equals the sequential one; a region may
+    complete where the sequential order would have overflowed in an intermediate prefix and the
+    interpreter's redo would have given the same value (the deopt counter is not part of the identity). **Float `+` and `*` reductions are excluded in v1**: chunked evaluation
     is a different expression tree and not bit-identical. This is not permanent: a fastmath opt-in
     (#153) awaits the owner's decision, and N-15 follows it. Partials combine in chunk order.
   - **Observability.** N-14 holds (a traced call runs the interpreted function, never a region);
@@ -914,7 +926,7 @@ All in `PM/commonTest`, so they run wherever the interpreter loads.
     chunk, NaN/inf/-0.0), identical bytes for every thread count, failure at first/middle/last
     iteration, negative qualification cases, disabled-path test, reductions, race and GC tests,
     TSan/ASan, N-14, signals, fork and shutdown, and measurement tests for the threshold.
-  `Status: step 1 implemented (#277): element-wise F64 stores a[i] = <F64 + - * and negation of constants, locals and b[i]> into list copy-in or local arrays; no reductions, no failure path (eligible loops cannot fail or deopt); the rest planned (design)`.
+  `Status: steps 1 and 2 implemented (#277): element-wise F64 stores a[i] = <F64 + - * and negation of constants, locals and b[i]> into list copy-in or local arrays (no failure path: eligible loops cannot fail or deopt); i64 reductions acc += e and acc *= e with the overflow rerun above (stores and reductions never in one loop); & | ^ and min/max reductions, float reductions, the poll inside a region and the rest planned (design)`.
 
 ---
 

python-multiplatform-ksp/src/main/python/typedpython/cgen.py

@@ -170,6 +170,7 @@
 """
 from __future__ import annotations
 
+import re
 from pathlib import Path
 
 from typedpython import ir
@@ -2145,11 +2146,32 @@ def _emit_region(self, s: ir.ForRange, ts: str, tn: str, mag: int) -> None:
                             "function, bad step or no ParInfo)")
         scalars: set[str] = set()
         arrays: set[str] = set()
-        for st in s.body:
-            if not isinstance(st, ir.StoreIndex):
-                raise CGenError(f"{f.name}: parallel loop body has a {type(st).__name__}")
-            arrays.add(st.array)
-            self._par_collect(st.value, scalars, arrays)
+        reds: list[tuple[str, str, ir.Expr]] = []         # (accumulator, op name, the other operand)
+        if info.reduce:
+            if info.stored or len(info.reduce) != len(s.body):
+                raise CGenError(f"{f.name}: parallel reduction loop with stores or a missing reduction")
+            for st, (acc, opn) in zip(sorted(s.body, key=lambda x: getattr(x, "target", "")), info.reduce):
+                v = st.value if isinstance(st, ir.Assign) else None
+                if (not isinstance(v, ir.BinOp) or st.target != acc or v.op.name != opn
+                        or opn not in ("ADD", "MUL") or v.type is not Type.I64
+                        or self.types.get(acc) is not Type.I64 or acc == s.var):
+                    raise CGenError(f"{f.name}: parallel reduction statement is not `acc = acc op e`")
+                if isinstance(v.left, ir.Local) and v.left.name == acc:
+                    other = v.right
+                elif isinstance(v.right, ir.Local) and v.right.name == acc:
+                    other = v.left
+                else:
+                    raise CGenError(f"{f.name}: parallel reduction of {acc} does not read {acc}")
+                reds.append((acc, opn, other))
+                self._par_collect(other, scalars, arrays)
+            if any(x in scalars for x, _, _ in reds):
+                raise CGenError(f"{f.name}: a parallel reduction reads its accumulator")
+        else:
+            for st in s.body:
+                if not isinstance(st, ir.StoreIndex):
+                    raise CGenError(f"{f.name}: parallel loop body has a {type(st).__name__}")
+                arrays.add(st.array)
+                self._par_collect(st.value, scalars, arrays)
         scalars.discard(s.var)
         for a in arrays:
             if a not in self.arrays and a not in self.local_arrays:
@@ -2180,6 +2202,8 @@ def _emit_region(self, s: ir.ForRange, ts: str, tn: str, mag: int) -> None:
             fields.append(f"    {ct} {ident};")
             decls.append(f"    {ct} {ident} = c->{ident};")
             sets.append(f"tp_pr.{ident} = {ident};")
+        if reds:
+            fields += ["    unsigned char *tp_pst;"] + [f"    int64_t *tp_part{j};" for j in range(len(reds))]
         saved = (self.lines, self.temps, self.depth, self.brk)
         self.lines, self.temps, self.depth, self.brk = [], [], 2, []
         try:
@@ -2188,30 +2212,73 @@ def _emit_region(self, s: ir.ForRange, ts: str, tn: str, mag: int) -> None:
         finally:
             self.lines, self.temps, self.depth, self.brk = saved
         text = "\n".join(body_lines)
+        if reds:
+            # the one failure path a reduction body has: the sequential overflow check, which here
+            # marks the chunk bad (the worker's own tp_rc / tp_exit); nothing else may remain
+            text = re.sub(r"^\s*tp_s = tp_(add|sub|mul)_i64\(.*\);\n?", "", text, flags=re.M)
+            text = re.sub(r"^\s*if \(tp_s != 0\) \{ tp_rc = tp_s; goto tp_exit; \}\n?", "", text, flags=re.M)
         for banned in ("goto", "tp_rc", "tp_s ", "return", "tp_st", "tp_dict", "PyObject", "tp_release"):
             if banned in text:
                 raise CGenError(f"{f.name}: parallel loop body emitted `{banned}`")
         sign_var = _ident("l", s.var)
         out = [f"typedef struct {{", *fields, f"}} {ctx_t};", "",
-               f"static void {body_f}(const void *tp_c, uint64_t tp_k0, uint64_t tp_k1)", "{",
+               f"static void {body_f}(const void *tp_c, unsigned tp_ci, uint64_t tp_k0, uint64_t tp_k1)", "{",
                f"    const {ctx_t} *c = (const {ctx_t} *)tp_c;",
                "    uint64_t tp_k;", *decls]
+        if reds:
+            out += ["    int tp_rc = 0;", "    int tp_s = 0;"]
+            out += [f"    int64_t {_ident('l', acc)} = {'1' if opn == 'MUL' else '0'};" for acc, opn, _ in reds]
+            out.append("    (void)tp_rc; (void)tp_s;")
         for ct, n in body_temps:
             out.append(f"    {ct}{'' if ct.endswith('*') else ' '}{n} = 0;")
         out.append("    (void)c;" + "".join(f" (void){_ident('l', x)};" for x in sorted(scalars)))
         out.append(f"    for (tp_k = tp_k0; tp_k < tp_k1; tp_k++) {{")
         out.append(f"        {sign_var} = (int64_t)((uint64_t)c->tp_ts + tp_k * UINT64_C({mag}));")
         out += body_lines
-        out += ["    }", "}", ""]
+        out.append("    }")
+        if reds:
+            out += [f"    c->tp_part{j}[tp_ci] = {_ident('l', acc)};" for j, (acc, _, _) in enumerate(reds)]
+            out += ["    c->tp_pst[tp_ci] = TP_PAR_OK;", "    return;", "tp_exit:",
+                    "    c->tp_pst[tp_ci] = TP_PAR_BAD;"]
+        out += ["}", ""]
         self.par_defs.append("\n".join(out))
         w = max(1, info.weight)
+        pok = self.tmp("int") if reds else None
+        if reds:
+            self.emit(f"{pok} = 0;")
         self.emit(f"if ({tn} >= (uint64_t)((TP_PAR_MIN_WORK + {w - 1}) / {w}) && tp_par_threads() > 1) {{")
         self.depth += 1
         self.emit(f"{ctx_t} tp_pr;")
+        if reds:
+            self.emit("unsigned char tp_pst[TP_PAR_HARD_CAP] = {0};")
+            for j in range(len(reds)):
+                self.emit(f"int64_t tp_pp{j}[TP_PAR_HARD_CAP];")
         self.emit(f"tp_pr.tp_ts = {ts};")
         for line in sets:
             self.emit(line)
+        if reds:
+            self.emit("tp_pr.tp_pst = tp_pst;")
+            for j in range(len(reds)):
+                self.emit(f"tp_pr.tp_part{j} = tp_pp{j};")
         self.emit(f"(void)tp_par_run({body_f}, &tp_pr, {tn});")
+        if reds:
+            # Combine in chunk order with the sequential checks. All or none: an overflow anywhere
+            # (a chunk, or here) leaves every accumulator untouched and the sequential loop runs.
+            self.emit(f"{pok} = 1;")
+            for j, (acc, opn, _) in enumerate(reds):
+                self.emit(f"int64_t tp_rv{j} = {_ident('l', acc)};")
+                self.emit(f"{pok} = {pok} && tp_par_combine(TP_PAR_OP_{opn}, &tp_rv{j}, tp_pp{j}, tp_pst) == 0;")
+            self.emit(f"if ({pok}) {{")
+            for j, (acc, _, _) in enumerate(reds):
+                self.emit(f"    {_ident('l', acc)} = tp_rv{j};")
+            self.emit(f"    {_ident('l', s.var)} = (int64_t)((uint64_t){ts} + ({tn} - 1u) * UINT64_C({mag}));")
+            self.emit("} else {")
+            self.emit("    tp_par_rerun_note();")
+            self.emit("}")
+            self.depth -= 1
+            self.emit("}")
+            self.emit(f"if (!{pok}) {{")
+            return
         self.emit(f"{_ident('l', s.var)} = (int64_t)((uint64_t){ts} + ({tn} - 1u) * UINT64_C({mag}));")
         self.depth -= 1
         self.emit("} else {")
@@ -2756,14 +2823,14 @@ def _state_functions(mg: _ModGen) -> str:
 
 def _par_stats_function() -> str:
     """`__typedpython_parallel__()`: (regions that used workers, workers started, failed thread
-    creations) since import, so a test can prove the parallel path ran (SPEC N-15)."""
+    creations, reduction regions rerun sequentially after an overflow) since import, so a test can prove the parallel path ran (SPEC N-15)."""
     return "\n".join([
         "static PyObject *tp_par_stats_fn(PyObject *tp_module, PyObject *tp_unused)",
         "{",
-        "    int64_t tp_v[3];",
+        "    int64_t tp_v[4];",
         "    (void)tp_module; (void)tp_unused;",
         "    tp_par_stats(tp_v);",
-        "    return Py_BuildValue(\"(LLL)\", (long long)tp_v[0], (long long)tp_v[1], (long long)tp_v[2]);",
+        "    return Py_BuildValue(\"(LLLL)\", (long long)tp_v[0], (long long)tp_v[1], (long long)tp_v[2], (long long)tp_v[3]);",
         "}",
         "",
         "static PyMethodDef tp_par_stats_def = {\"__typedpython_parallel__\", tp_par_stats_fn, METH_NOARGS, NULL};",

python-multiplatform-ksp/src/main/python/typedpython/ir.py

@@ -653,7 +653,7 @@ class While(Stmt):
 @dataclass(frozen=True)
 class ParInfo:
     """Why a `ForRange` may run its iterations on several threads (SPEC N-15, step 1): the body is
-    only `a[var] = <F64 expression>` stores, `stored` are the arrays stored into (each at the loop
+    only `a[var] = <F64 expression>` stores (or, step 2, only i64 `acc op= e` reductions: see `reduce`), `stored` are the arrays stored into (each at the loop
     variable, proven in bounds, so distinct iterations write distinct elements) and `read` the arrays
     read, all at the loop variable and proven in bounds. The expression nodes are Const, scalar
     Local (not assigned in the body), ToFloat of an I64, F64 + - * and negation, so no node can fail,
@@ -664,6 +664,11 @@ class ParInfo:
     stored: tuple[str, ...]
     read: tuple[str, ...]
     weight: int = 1
+    # Step 2 (#277): the i64 reductions `acc = acc op e` of a loop whose body is only such statements,
+    # as sorted (acc, BinOpKind name) pairs, op ADD or MUL. Then `stored` is empty. Each chunk computes
+    # its partial from the op's identity with the sequential code's overflow checks; an overflow in any
+    # chunk or in the combination discards the region and the sequential loop runs from its start.
+    reduce: tuple[tuple[str, str], ...] = ()
 
 
 @dataclass(frozen=True)

python-multiplatform-ksp/src/main/python/typedpython/runtime/API.md

@@ -492,7 +492,7 @@ for the epoch's increment. `tp_runtime.h` needs the GCC/Clang `__atomic` builtin
 
 ## Parallel regions (SPEC N-15, issues #151, #277)
 A `ForRange` the verifier marked `par` (`ir.ParInfo`; marked function, `a[i] = F64 expression` stores
-only) runs its iterations on per-region C threads above a work threshold. Workers touch no Python
+or i64 `+ *` reductions) runs its iterations on per-region C threads above a work threshold. Workers touch no Python
 object and no Python API, are not Python threads and never outlive the call.
     #define TP_PAR_MIN_WORK 65536          /* PLACEHOLDER, to be measured on #151: weighted iterations (trip count * ParInfo.weight) */
     #define TP_PAR_MAX_THREADS 8           /* default cap; TYPEDPYTHON_THREADS overrides (1 disables), hard cap 64 */
@@ -501,21 +501,33 @@ object and no Python API, are not Python threads and never outlive the call.
     #define TP_PAR_ALIGN 8u                /* chunk starts are multiples of 8 iterations (one cache line of doubles) */
     #define TP_PAR_HAVE_THREADS            /* defined with pthreads; -DTP_NO_THREADS compiles the threads out */
     #define TP_PAR_TEST_FAIL_CREATE        /* test hook: every pthread_create reports failure */
-    typedef void (*tp_par_fn)(const void *ctx, uint64_t lo, uint64_t hi);
+    typedef void (*tp_par_fn)(const void *ctx, unsigned ci, uint64_t lo, uint64_t hi);   /* ci: chunk index, 0.. */
     int  tp_par_threads(void);                                   /* the cap, read once from TYPEDPYTHON_THREADS, else min(cores, 8) */
     int  tp_par_run(tp_par_fn fn, const void *ctx, uint64_t n);  /* runs [0, n) in contiguous chunks, joined before return; returns workers started */
-    void tp_par_stats(int64_t out[3]);                           /* regions with workers, workers started, failed creations */
+    void tp_par_stats(int64_t out[4]);                           /* regions with workers, workers started, failed creations, reduction reruns */
+    #define TP_PAR_NONE 0 / TP_PAR_OK 1 / TP_PAR_BAD 2    /* per-chunk state of a reduction chunk (zero-initialised array of TP_PAR_HARD_CAP) */
+    #define TP_PAR_OP_ADD 0 / TP_PAR_OP_MUL 1
+    int  tp_par_combine(int op, int64_t *acc, const int64_t *part, const unsigned char *st);  /* acc op= part[0], part[1], ... in chunk order, checked; 1 if a chunk was BAD or it overflowed (acc untouched) */
+    void tp_par_rerun_note(void);                                /* counts a discarded reduction region */
 `tp_par_run` always completes the whole range. The calling thread runs the first chunk; a chunk whose
 `pthread_create` failed runs on the calling thread after it (fewer workers, down to none: the sequential
 loop, no exception, not a semantic event). Workers start with all signals blocked, so Python signal
 handlers run on Python threads only. The chunk function must be safe for any chunk on any thread;
 generated code passes a context struct of raw array pointers and scalar copies. A module with a
 `par` loop also exports `__typedpython_parallel__()` returning `(regions, workers_started,
-failed_creations)`, so tests can prove the parallel path ran. Free-threaded state: the four counters and
+failed_creations, reduction_reruns)`, so tests can prove the parallel path ran. Free-threaded state: the four counters and
 the thread cap are atomics (relaxed); a region's chunk bounds, jobs and thread handles live on the
 caller's stack. The caller does not poll the eval breaker inside a region (array functions do not poll
 anyway, #152).
 
+Reduction regions (step 2). The chunk function starts each accumulator at the op's identity (0 for `+`,
+1 for `*`), runs the sequential statements `acc = acc op e` with the sequential overflow checks, and
+writes `part[ci]` and `st[ci] = TP_PAR_OK`; a failed check sets `st[ci] = TP_PAR_BAD` and stops the chunk.
+After `tp_par_run` the caller combines every accumulator with `tp_par_combine`; if any accumulator fails
+(a BAD chunk or an overflow in the combination) no accumulator and no loop variable is written, the rerun
+counter is incremented and the sequential loop runs from the loop's start state, which deopts or not
+exactly as today. (`& | ^` and min/max are not reductions yet: the compiled subset has no bitwise operator.)
+
 ## Development-mode import gate (Python API, not C; SPEC N-7b)
     typedpython.importgate.install(source_roots)   # meta path finder; checks project modules before they run
     typedpython.importgate.uninstall()