cmd/compile: archsimd MulAdd cannot accumulate in place or fold a broadcast memory operand (no VFMADD231PS, no .BCST), though the assembler encodes both

[scttfrdmn]

Go version

go1.26.5 (toolchain host darwin/arm64, cross-compiling and executing on linux/amd64), GOEXPERIMENT=simd.

What did you do?

Wrote the inner loop of an SGEMM microkernel — one row of a 1×32 tile, i.e.
c[0:32] += a[k] * b[k][0:32] accumulated over k — the canonical shape for
MulAdd on packed panels:

package main

import "simd/archsimd"

//go:noinline
func row(a, b, c []float32, kc int) {
	c0 := archsimd.LoadFloat32x16Slice(c[0:16])
	c1 := archsimd.LoadFloat32x16Slice(c[16:32])
	ap, bp := a[:kc], b[:kc*32]
	for len(ap) >= 1 && len(bp) >= 32 {
		av := archsimd.BroadcastFloat32x16(ap[0])
		c0 = av.MulAdd(archsimd.LoadFloat32x16Slice(bp[0:16]), c0)
		c1 = av.MulAdd(archsimd.LoadFloat32x16Slice(bp[16:32]), c1)
		ap, bp = ap[1:], bp[32:]
	}
	c0.StoreSlice(c[0:16])
	c1.StoreSlice(c[16:32])
}

What did you expect to see?

Two instructions of arithmetic per k-step, plus the two B-panel loads:

VMOVDQU64        (DI), Z3
VMOVDQU64      64(DI), Z4
VFMADD231PS.BCST (AX), Z3, Z0     // Z0 += Z3 * broadcast(float32 at (AX))
VFMADD231PS.BCST (AX), Z4, Z1

Go's assembler already emits exactly that, including the EVEX embedded
broadcast, so this is not an ISA or an assembler limitation:

$ cat x_amd64.s
TEXT ·probe(SB), NOSPLIT, $0-0
	VFMADD231PS.BCST 12(SI), Z1, Z0
	RET

$ GOARCH=amd64 GOOS=linux go build -o probe . && llvm-objdump -d probe | grep -A1 probe.abi0
000000000047a7c0 <main.probe.abi0>:
  47a7c0: 62 f2 75 58 b8 46 03  vfmadd231ps 0xc(%rsi){1to16}, %zmm1, %zmm0 # zmm0 = (zmm1 * mem) + zmm0

$ echo 'vfmadd231ps 12(%rsi){1to16}, %zmm1, %zmm0' | llvm-mc -triple=x86_64 --show-encoding
	vfmadd231ps 12(%rsi){1to16}, %zmm1, %zmm0 # encoding: [0x62,0xf2,0x75,0x58,0xb8,0x46,0x03]

(Byte-identical. go tool objdump cannot decode EVEX, which is why llvm-objdump
is used above.)

What did you see instead?

Nine instructions of arithmetic path per k-step for two FMAs (full steady-state
loop is 27 instructions; the rest is loop control and bounds-check-free panel
re-slicing):

00088 (main.go:14)  VMOVSS        (AX), X2      // a[k]
00129 (broadcast)   VBROADCASTSS  X2, Z2        // ... in two instructions
00135 (load)        VMOVDQU64     (DI), Z3
00141 (load)        VMOVDQU64   64(DI), Z4
00148 (main.go:15)  VMOVDQU64     Z2, Z5        // copy: 213 clobbers arg0
00154 (main.go:15)  VFMADD213PS   Z0, Z3, Z2    // Z2 = Z2*Z3 + Z0
00160 (main.go:16)  VFMADD213PS   Z1, Z4, Z5
00179 (main.go:13)  VMOVDQU64     Z5, Z1        // move results back to the
00185 (main.go:13)  VMOVDQU64     Z2, Z0        //   accumulators' registers

Three independent things block the expected form, which is why I think this is one
report rather than three:

1. There is no 231- or 132-shaped FMA SSA op for vector types. All 20
VFMADD* entries in _gen/simdAMD64ops.go are 213, plus their load and
Masked variants. 213 has resultInArg0: true and writes its first
multiplicand, so acc = a.MulAdd(b, acc) can never land in acc's register: it
needs a live scratch register beyond the working set, always, and a copy to get the
value home. The 231 shape is not foreign to the backend — VFMADD231SS/SD
exist at _gen/AMD64Ops.go:809–810 and are what math.FMA uses — it just does
not exist for vectors.

2. The one load-merging rule that exists can only fold the addend, which in a
GEMM is the accumulator. From _gen/simdAMD64.rules:2774:

(VFMADD213PS512 x y l:(VMOVDQUload512 {sym} [off] ptr mem)) && canMergeLoad(v, l) && clobber(l)
    => (VFMADD213PS512load {sym} [off] x y ptr mem)

213 computes arg0*arg1 + arg2, so the folded operand is arg2 — the addend.
x.MulAdd(y, z) is x*y + z, so in an accumulating loop z is the accumulator,
and folding it means reading the accumulator from memory every iteration. That is
the one operand a GEMM kernel must keep in a register, so the rule cannot fire
usefully here even though it is the right rule for other shapes. Had a 231 form
existed, its memory operand would have been a multiplicand, which is exactly what
a packed panel is.

3. Nothing in the SSA generators emits the .BCST suffix. grep -rn BCST src/cmd/compile/internal/ssa/_gen/ is empty, while src/cmd/internal/obj/x86/
(evex.go, asm6.go) supports it. So the embedded broadcast — the piece that
removes the VMOVSS/VBROADCASTSS pair entirely — is reachable from hand-written
assembly and from nothing else.

Why it matters

Measured on the shipped 2×32 microkernel of a pure-Go float32 BLAS subset
(https://github.com/scttfrdmn/keel), steady-state loop, 16 FMAs per pass:

body	insns	insns/FMA	fraction of measured peak on an i9-9960X
as emitted today	74	4.625	46.0% (measured)
with a 231 form (folds the copies)	66	4.125	51.7% (projected)
with 231 + .BCST	50	3.125	68.2% (projected)
with 231 + .BCST + no statement-anchor NOPs	46	2.875	74.2% (projected)

The denominator is that host's own measured single-core register-only FMA
saturation rate, taken in the same benchmark invocation, not a formula. On that
part the kernel is at its front end's limit, not its FMA units': the two shipped
tile shapes' throughputs stand in the inverse ratio of their instruction counts
(1.308 measured against 1.351 predicted) and both derive the same ~4.2
instructions per cycle despite differing 35% in instruction count. So on a
2-FMA/cycle, ~4-wide-issue machine, instruction count is the performance, and 32
of these 74 instructions are the four groups above.

Two Zen hosts with wider front ends (7950X3D, Ryzen AI MAX+ 395) reach 96.6% and
64.2% of their measured peaks with the same source, so this is not a
whole-fleet regression — it is a specific, quantified gap between what the
assembler can encode and what the intrinsics can reach, and it binds hardest on
exactly the parts with two full-width 512-bit FMA units.

Note that the 231 form alone is not sufficient: folding the accumulate but not
the broadcast leaves 66 instructions and is projected at 51.7%, still short. The
embedded broadcast is the larger of the two levers.

To be clear about scope: this is not "archsimd is missing an operation".
MulAdd is exactly the right source-level primitive and the source above is the
natural way to write the loop. It is a lowering gap, and each of the three parts
looks independently fixable — a 231 op set with a load variant that folds a
multiplicand, and a rule that recognizes a splat feeding an FMA and emits .BCST.

Related, and measured in the same loop: BroadcastFloat32x16 is an emulated op
(SetElem + Broadcast1To16) rather than a VBROADCASTSS, which is the
VMOVSS+VBROADCASTSS pair above. Filed separately.

cc #73787 #78979

[gabyhelp]

Related Issues

• simd/archsimd: AVX instructions do not make use of memory operands #78159 (closed)
• cmd/compile: GOEXPERIMENT=simd allocates 512-bit vector values from 15 of the 32 zmm registers #80828
• cmd/compile: no automatic use of fused multiply-add on amd64 even with GOAMD64=v3 #71204 (closed)
• simd/archsimd: BroadcastFloat32x16 is emulated, and every inlined non-intrinsic wrapper costs a 1-byte anchor NOP in the loop body #80830
• cmd/compile: compiler doesn't optmize fused multiply-subtract when using GOAMD64=v3 #80278 (closed)
• cmd/compile: code generates unnecessary mov-s #65192
• cmd/internal/obj/x86: add fma #8037 (closed)
• simd/archsimd: missing combined load/op optimization #79408 (closed)
• simd/archsimd: archsimd.LoadUint8x32Slice() emits nops #76968 (closed)
• x/arch: incorrect disassembly for x86 FMA instruction #41043

_{(Emoji vote if this was helpful or unhelpful; more detailed feedback welcome in this discussion.)}

[scttfrdmn]

Measured impact of this issue on a real float32 GEMM microkernel, on Sapphire Rapids
(2 FMA ports), go1.26.5, GOEXPERIMENT=simd, GOMAXPROCS=1.

The kernel. 4 rows x 32 columns, eight archsimd.Float32x16 accumulators, packed
panels, no calls in the K-loop. The steady-state loop body is 50 instructions for 8
FMAs = 6.25 instructions per FMA. Of those 50, 8 are VMOVDQU64 copy-backs of the
accumulators, one per accumulator per iteration.

They are there because of the operand order this issue describes. The emitted form is

VFMADD213PS Z7, Z12, Z8      // Z8 = Z12*Z8 + Z7

so the sum lands in the register holding the broadcast operand, not in the accumulator.
With no 231 form reachable from the intrinsic layer, every accumulator has to be copied
back each iteration.

The denominator. The same binary's FMA-saturation loop — register-only, audited at
2.25 instructions per FMA — measures 232.3 GFLOP/s +/- 0.0% on this host. The
microkernel measures 79.35 GFLOP/s, i.e. 34.16% of that. Taking the retirement
model attainable = 2.25 / (insns per FMA), 6.25 predicts 36.0% and the kernel achieves
94.9% of that prediction, so the instruction count explains the rate.

The delta. Removing only the 8 copy-backs — this issue's in-place accumulate, nothing
else — takes the body from 50 to 42 instructions, 6.25 to 5.25 insns/FMA, and the
attainable ceiling from 36.0% to 42.9% of FMA peak. That is a +19% relative
improvement in the ceiling for this kernel from this issue alone.

Folding the broadcast as a .BCST memory operand — the second half of this issue — would
remove 4 more instructions: 38, 4.75 insns/FMA, 47.4%. Worth noting for anyone
comparing with #78159, which is closed as completed: memory operands do fold here, but
the load-merge rule folds memory into the addend, which is the operand the accumulator
needs to occupy. So the fold landing and this issue being open are consistent, and the
fold as it exists cannot substitute for the 231 form.

A control, since percent-of-peak on a shared cloud host invites the obvious objection.
Three readings, two guest sizes including the whole 2-socket host:

guest	peak GFLOP/s	kernel GFLOP/s	kernel / peak
c7i.4xlarge (8 cores)	228.9	78.32	34.216%
c7i.4xlarge (8 cores)	240.2	82.06	34.163%
c7i.48xlarge (96 cores, 2 sockets)	232.3	79.35	34.158%

The peaks span 4.9% and the fraction holds to 0.057 points (0.17% relative), with the
whole-host reading in the middle of the peak range rather than at an end. Numerator and
denominator are measured in the same process on the same silicon, so the ratio divides the
clock out; a neighbour effect would have to reproduce itself to two parts in a thousand
three times.

One more datum that may help prioritise, because it is microarchitecture-specific. The
same accumulator pressure inverts which shape wins depending on FMA port count. A 2-row x
32-column shape has fewer instructions per FMA (4.625) but its FMAs form a 4-long serial
dependency chain; the 4-row shape has 6.25 insns/FMA and eight independent chains. On a
1-FMA/cycle host the 2-row shape is the faster of the two; on this 2-FMA/cycle host the
4-row shape wins, 34.16% versus 26.5% of the same peak. So the copy-back tax is paid
hardest exactly where the extra accumulators are needed to cover latency x throughput, and
the register budget those accumulators come from is #80828's 15 of 32.

Context, not a request: this is from a pure-Go float32 BLAS subset built on
GOEXPERIMENT=simd, and stripping this issue plus #80830's anchor NOPs from the model
still lands at 50.0% of FMA peak, against a 55% target — so this is not the whole gap on
our side, just the largest single piece of it that is not ours.

[scttfrdmn]

Correction to the last paragraph of my comment above, plus the figure it should have
carried. The measured delta on this issue is unaffected — 8 copy-backs, 50 to 42
instructions, 36.0% to 42.9% of FMA peak — but the mechanism I offered for the shape
comparison was wrong, and the corrected one argues more strongly for this issue.

First, the number I quoted as a bare percentage: 2x32 measures 61.45 GFLOP/s against
the same 232.3 GFLOP/s peak = 26.45%, same run and same binary as the 4x32 row
(c7i.48xlarge, keel rev 75bfaf5).

A clock-free comparison. Multiplying each shape's audited insns/FMA by its achieved
rate gives its instruction throughput relative to the FMA-saturation loop measured in the
same process, which divides the clock out and needs no pipeline constants at all:

shape	insns/FMA	GFLOP/s	% of the peak loop's own instruction rate
4x32, 8 accumulators	6.250	79.35	94.9%
2x32, 4 accumulators	4.625	61.45	54.4%

4x32 — the shape that wins on this host — runs at 94.9% of the instruction throughput the
same silicon demonstrably sustains. That is the claim worth making here: the 8 copy-backs
sit on the one axis that is measurably binding for the winning shape. So the ceiling
arithmetic above (36.0% to 42.9%) is the gain to expect, not an upper bound that some other
limit would absorb.

What I got wrong. I attributed the shape inversion to accumulator pressure — "the
copy-back tax is paid hardest exactly where the extra accumulators are needed to cover
latency x throughput." Accumulator count binds neither shape:

• 4x32 has 8 accumulators at 4-cycle latency on 2 ports, so its dependency term is 4.00
  cycles — exactly tied with its port term, against a measured body of 11.71 cycles. A >= P*L
  makes that term non-binding at any unroll, because the FMA count cancels out of (F/A)*L
  versus F/P.
• 2x32's 4-long serial chain floors its body at 16.00 cycles. It measures 30.24. The
  chain is 1.89x too loose to be the limit, and only 54% of its instruction supply is used,
  so 2x32 is bounded by neither term. Something outside the static model stalls it; the
  candidate the model does supply is the operand feed, which is exact arithmetic —
  1/MR + 1/V is 1.00 vectors/FMA for 2x32 against 0.75 for 4x32.

So the inversion establishes something narrower than I claimed: 4x32 wins because it
reaches its instruction ceiling, not because 8 accumulators cover latency x throughput. A
shape carrying more accumulators would be strictly worse here on both terms that bind, which
points the lever at instruction count — this issue — rather than at register pressure.

The only pipeline constant used above is 2 FMA ports, for the clock in 232.3 / 64. Sapphire
Rapids' issue width is not used anywhere in it, and I have not sourced it; the
instruction-throughput column needs neither.