simd(s390x): vector max/min on z14+ (VFMAXDB Java mode)

README.md

@@ -68,7 +68,7 @@ loops are [go-asmgen](https://github.com/go-asmgen)-generated SIMD kernels: sum,
 add/sub/mul/div, sqrt and the dot product (SSE2 or AVX2/FMA on amd64, NEON on
 arm64, VSX on ppc64le, LASX on loong64, RVV on riscv64 and the vector facility
 on s390x, the last three when the CPU has them), max/min too on
-amd64 and arm64, and a **panel-packed, cache-blocked GEMM** with an SIMD-FMA
+amd64, arm64 and s390x (z14 and later), and a **panel-packed, cache-blocked GEMM** with an SIMD-FMA
 micro-kernel (NEON 4×8; AVX2/FMA 6×8 chosen at run time, SSE2 fallback; VSX
 8×8; LASX 8×8; RVV 6×8; s390x 4×8), the OpenBLAS/BLIS structure. The 32-bit
 targets run the same pure-Go code those kernels are tested against. `Exp` and `Log` are ports of Arm's optimized-routines

SECURITY.md

@@ -57,6 +57,14 @@ standard library's known vulnerabilities; `govulncheck ./...` reports none.
 
 ## Audit log
 
+**2026-10-09 (v0.9.1).** New s390x `maxVXE`/`minVXE` kernels, gated on
+AT_HWCAP VXRS_EXT (z14 and later), each behind its wrapper: from eight
+elements up, with no pointer arithmetic beyond `len(a)`. On the real z17 they
+pass `TestMaxMinSIMD`, the specials-in-every-lane test, the fence tests and the
+parallel stress test. A copy using mode 0 (IEEE maxNum, which returns the
+number when the other operand is a NaN) fails there, so the NaN rule is
+checked, not assumed.
+
 **2026-10-09 (v0.9.0).**
 
 - **The Go floor is 1.27.2.** At 1.27.1, `govulncheck` listed standard-library

docs/perf.md

@@ -880,8 +880,16 @@ stay dynamic, so its slower efficiency cores are handled as before.
   two-lane registers plus a scalar tail, every instruction a Go mnemonic.
   Since v0.9.0 the GEMM micro-kernel is a 4×8 tile: sixteen two-lane
   accumulators, the B row in one `VLM`, each A value replicated with `VLREPG`
-  and fused in with `VFMADB`. max/min stay scalar: the vector max with NumPy's
-  NaN rule needs the z14, beyond the z13 Go supports. The operand order of
+  and fused in with `VFMADB`. Since v0.9.1, max/min are vector too on a z14
+  or later: `VFMAXDB`/`VFMINDB` in mode 1, Java's Math.max/min, which gives NaN
+  for a NaN operand and orders −0 below +0, NumPy's rule. They come with the
+  vector-enhancements facility 1 (AT_HWCAP VXRS_EXT), beyond the z13 that Go
+  also supports, where the scalar reducers remain. On the z17, one core,
+  `benchstat` over 8 interleaved runs, all at p = 0.000: `Max` takes 92% to 96%
+  less time (1 Ki: 2.58 µs → 112 ns; 4 Mi: 10.7 ms → 0.81 ms), and `MaxAxis(1)`
+  94% less. `Sum`, untouched, measured within ±0.6%: the control. A sabotaged copy in
+  mode 0 (IEEE maxNum, which drops a NaN) failed `TestMaxMinSIMD` on the z17.
+  The operand order of
   `VFSDB`/`VFDDB`/`VFMADB` was **measured**, not read: each run once on known
   inputs (`VFSDB V1, V2, V3` is V3 = V2 − V1). AT_HWCAP is read in the
   machine's byte order, which on big-endian s390x is not the little-endian the

docs/plan-ndarray.md

@@ -127,7 +127,8 @@ suite (plus the per-arch CI jobs already wired in `.github/workflows/ci.yml`).
     reports LASX), and since v0.7.0 riscv64 has RVV kernels (sum, dot, sqrt,
     elementwise, and since v0.7.1 a 6×8 GEMM micro-kernel; used when the CPU
     reports V), and since v0.8.0 s390x has vector-facility kernels (sum, dot,
-    sqrt, elementwise, and since v0.9.0 a 4×8 GEMM micro-kernel; measured on a
+    sqrt, elementwise, since v0.9.0 a 4×8 GEMM micro-kernel and since v0.9.1
+    max/min on z14 and later; measured on a
     real IBM z17). The
     32-bit targets keep the
     validated scalar oracles plus a scalar 4×4 GEMM micro-kernel over the

internal/kernels/asmgen/s390x/gen.go

@@ -37,6 +37,8 @@ func main() {
 		f.Add(binKernel(k.name, k.vec, k.sca))
 	}
 	f.Add(gemmKernel())
+	f.Add(extremumKernel("maxVXE", "VFMAXDB", "WFMAXDB"))
+	f.Add(extremumKernel("minVXE", "VFMINDB", "WFMINDB"))
 	if err := os.WriteFile("sum_s390x.s", []byte(f.String()), 0o644); err != nil {
 		fmt.Fprintln(os.Stderr, err)
 		os.Exit(1)
@@ -178,3 +180,39 @@ func gemmKernel() *emit.Function {
 	}
 	return b.Ret().Func()
 }
+
+// extremumKernel builds NAME(a *float64, n int) float64, the maximum (or
+// minimum) of a, for n >= 8. VFMAXDB/VFMINDB and their scalar WF forms come
+// with the vector-enhancements facility 1 (z14), so the Go side gates these
+// on AT_HWCAP's VXRS_EXT. Mode 1 is Java's Math.max/min: a NaN operand gives
+// NaN, and -0 orders below +0, which is NumPy's and Go's builtin max/min rule.
+// The operation is associative under it, so four two-lane accumulators seeded
+// with the first eight elements, folded at the end, give the serial result.
+func extremumKernel(name, vec, sca string) *emit.Function {
+	sig := s390x.Layout([]string{"a", "n"}, []s390x.Type{s390x.Ptr, s390x.Int64},
+		[]string{"ret"}, []s390x.Type{s390x.Float64})
+	b := s390x.NewFunc(name, sig, 0)
+	b.LoadArg("a", "R1").LoadArg("n", "R2").
+		Raw("VLM (R1), V16, V19").
+		Raw("ADD $64, R1").Raw("SUB $8, R2").
+		Raw("loop:").Raw("CMP R2, $8").Raw("BLT fold").
+		Raw("VLM (R1), V0, V3")
+	for r := 0; r < 4; r++ {
+		b.Raw("%s $1, V%d, V%d, V%d", vec, r, 16+r, 16+r)
+	}
+	b.Raw("ADD $64, R1").Raw("SUB $8, R2").Raw("BR loop").
+		Raw("fold:").
+		Raw("%s $1, V17, V16, V16", vec).
+		Raw("%s $1, V19, V18, V18", vec).
+		Raw("%s $1, V18, V16, V0", vec).
+		Raw("VREPG $1, V0, V1").
+		Raw("%s $1, F1, F0, F0", sca). // F0 = extremum of the two lanes
+		Raw("tail:").Raw("CMP R2, $0").Raw("BEQ done").
+		Raw("FMOVD (R1), F1").
+		Raw("%s $1, F1, F0, F0", sca).
+		Raw("ADD $8, R1").Raw("SUB $1, R2").Raw("BR tail").
+		Raw("done:").
+		StoreRet("F0", "ret").
+		Ret()
+	return b.Func()
+}

internal/kernels/fence_s390x_test.go

@@ -2,8 +2,8 @@ package kernels
 
 // forEachFallback runs f again with the scalar kernels used without the vector facility.
 func forEachFallback(f func()) {
-	saved := useVX
-	useVX = false
-	defer func() { useVX = saved }()
+	saved, savedE := useVX, useVXE
+	useVX, useVXE = false, false
+	defer func() { useVX, useVXE = saved, savedE }()
 	f()
 }

internal/kernels/reduce_s390x.go

@@ -10,9 +10,9 @@ import "os"
 //
 //go:generate sh -c "cd asmgen/s390x && go run . && mv sum_s390x.s ../../sum_s390x.s"
 //
-// max/min stay on the scalar reducers: the vector max with NumPy's NaN rule
-// (VFMAXDB, Java mode) needs the z14's vector-enhancements facility, beyond the
-// z13 Go supports.
+// max/min use the vector max with NumPy's NaN rule (VFMAXDB, Java mode) on
+// a z14 or later (vector-enhancements facility 1, AT_HWCAP VXRS_EXT), and the
+// scalar reducers on the z13 that Go also supports.
 
 // HaveReduceSIMD reports that this build ships SIMD kernels; they run when
 // useVX is true.
@@ -46,8 +46,34 @@ func sqrtSIMD(dst, src []float64) {
 	sqrtVX(&dst[0], &src[0], len(dst))
 }
 
-func maxSIMD(a []float64) float64 { return maxUnrolled(a) }
-func minSIMD(a []float64) float64 { return minUnrolled(a) }
+// hwcapVXE is HWCAP_S390_VXRS_EXT (bit 13), the vector-enhancements facility 1
+// of the z14, which brings VFMAXDB/VFMINDB.
+const hwcapVXE = 1 << 13
+
+// useVXE selects the vector max/min. A var so tests can force the scalar path.
+var useVXE = vxeFromAuxv()
+
+func vxeFromAuxv() bool {
+	auxv, err := os.ReadFile("/proc/self/auxv")
+	return err == nil && hwcapFrom(auxv)&hwcapVXE != 0
+}
+
+// maxSIMD / minSIMD are the NaN-propagating extrema of a (non-empty): the
+// VFMAXDB/VFMINDB kernel in Java mode from eight elements up on a z14 or
+// later, the four-chain scalar reducer otherwise.
+func maxSIMD(a []float64) float64 {
+	if !useVXE || len(a) < 8 {
+		return maxUnrolled(a)
+	}
+	return maxVXE(&a[0], len(a))
+}
+
+func minSIMD(a []float64) float64 {
+	if !useVXE || len(a) < 8 {
+		return minUnrolled(a)
+	}
+	return minVXE(&a[0], len(a))
+}
 
 func addBin(dst, a, b []float64) {
 	if !useVX || len(dst) == 0 {
@@ -114,3 +140,9 @@ func mulVX(dst, a, b *float64, n int)
 
 //go:noescape
 func divVX(dst, a, b *float64, n int)
+
+//go:noescape
+func maxVXE(a *float64, n int) float64
+
+//go:noescape
+func minVXE(a *float64, n int) float64

internal/kernels/sum_s390x.s

@@ -331,3 +331,73 @@ TEXT ·gemmMicro4x8VX(SB), NOSPLIT, $0-40
 	VSTM V24, V27, (R1)
 	RET
 
+TEXT ·maxVXE(SB), NOSPLIT, $0-24
+	MOVD a+0(FP), R1
+	MOVD n+8(FP), R2
+	VLM (R1), V16, V19
+	ADD $64, R1
+	SUB $8, R2
+	loop:
+	CMP R2, $8
+	BLT fold
+	VLM (R1), V0, V3
+	VFMAXDB $1, V0, V16, V16
+	VFMAXDB $1, V1, V17, V17
+	VFMAXDB $1, V2, V18, V18
+	VFMAXDB $1, V3, V19, V19
+	ADD $64, R1
+	SUB $8, R2
+	BR loop
+	fold:
+	VFMAXDB $1, V17, V16, V16
+	VFMAXDB $1, V19, V18, V18
+	VFMAXDB $1, V18, V16, V0
+	VREPG $1, V0, V1
+	WFMAXDB $1, F1, F0, F0
+	tail:
+	CMP R2, $0
+	BEQ done
+	FMOVD (R1), F1
+	WFMAXDB $1, F1, F0, F0
+	ADD $8, R1
+	SUB $1, R2
+	BR tail
+	done:
+	FMOVD F0, ret+16(FP)
+	RET
+
+TEXT ·minVXE(SB), NOSPLIT, $0-24
+	MOVD a+0(FP), R1
+	MOVD n+8(FP), R2
+	VLM (R1), V16, V19
+	ADD $64, R1
+	SUB $8, R2
+	loop:
+	CMP R2, $8
+	BLT fold
+	VLM (R1), V0, V3
+	VFMINDB $1, V0, V16, V16
+	VFMINDB $1, V1, V17, V17
+	VFMINDB $1, V2, V18, V18
+	VFMINDB $1, V3, V19, V19
+	ADD $64, R1
+	SUB $8, R2
+	BR loop
+	fold:
+	VFMINDB $1, V17, V16, V16
+	VFMINDB $1, V19, V18, V18
+	VFMINDB $1, V18, V16, V0
+	VREPG $1, V0, V1
+	WFMINDB $1, F1, F0, F0
+	tail:
+	CMP R2, $0
+	BEQ done
+	FMOVD (R1), F1
+	WFMINDB $1, F1, F0, F0
+	ADD $8, R1
+	SUB $1, R2
+	BR tail
+	done:
+	FMOVD F0, ret+16(FP)
+	RET
+

internal/kernels/vx_s390x_test.go

@@ -16,3 +16,16 @@ func TestVXDetected(t *testing.T) {
 
 	})
 }
+
+// TestVXEDetected runs the max/min checks with the vector kernel forced off
+// too, and logs which path this machine takes.
+func TestVXEDetected(t *testing.T) {
+	t.Logf("useVXE = %v (AT_HWCAP)", useVXE)
+	saved := useVXE
+	useVXE = false
+	defer func() { useVXE = saved }()
+	t.Run("scalar", func(t *testing.T) {
+		TestMaxMinSIMD(t)
+		TestMaxMinSIMDSpecialsInVectorLanes(t)
+	})
+}