ToT: conj() on nested-tile contractions; real-plain x complex-ToT products via one real gemm

[kshitij-05]

Follow-up to #573 (merged); rebased onto master.

Problem

1. conj() on a ToT contraction does not compile. conj(A * B) carries a ComplexConjugate<void> contraction factor; the ToT per-cell multiply-add ops in ContEngine static_cast that factor to the element type, Tensor::gemm forwards it as the BLAS alpha, and the ContractReduce<…, ComplexConjugate<…>> specializations name their result type through a value-returning gemm that nested tiles do not have.
2. A complex ToT contracted with a real plain tensor (A("i,j;a") * t("j,k"), Tensor<Tensor<complex<double>>> × Tensor<double>) skips the ToT × plain-tensor GEMM fast path, which requires identical element types, and runs the per-cell AXPY loop.

Changes

• detail::elem_factor<Numeric>(factor): the per-element multiplier for a contraction factor — the factor for numeric factors, 1 for ComplexConjugate<…> (conjugation and scale are applied to the finished result by ContractReduce's finalization, as for non-nested tiles). Used by the ToT per-cell ops in ContEngine and by Tensor::gemm's alpha.
• ContractReduce<…, ComplexConjugate<void>> / <…, ComplexConjugate<Scalar>>: result_type is Result, as in the primary template, and the nested-tile accumulate (which was an abort(); // not yet implemented stub) and the partial-result reduce are shared with the primary template through ContractReduceBase::accumulate_nested / reduce_results.
• ToT × real-plain products: when the plain element type is the real part of the inner element type, the complex slabs are viewed as real matrices with the inner extent doubled (re/im interleaved) and the existing strided GEMM runs in real arithmetic with alpha = beta = 1; results are identical to the complex product.
• Tensor::gemm's interleaved real-GEMM path takes its layout requirement (std::complex<Vr> is two contiguous Vr, alignof at least Vr's) and the ToT operand's inner-scalar match as terms of the interleaved predicate (detail::is_interleaved_real_view_v, shared by the tot_x_t and t_x_tot blocks), so an ABI that breaks the layout falls back to the per-cell loop instead of failing to compile or reinterpreting into UB.
• ContractReduce<…, ComplexConjugate<…>>::result_type is Result, as in the primary template, instead of the value type of a gemm(left, right, 1, helper) call (identical for homogeneous tiles; public-header change).
• Review round: the factor is normalized once in ContEngine (elem_scale() for the per-cell ops, outer_factor() for the outer tile op) instead of per consumer, so S*conj(A*B), -conj(A*B), the inner-Hadamard and complex-plain-operand forms, and batched (fused-mode) conj(A*B) all apply the conjugation and scale exactly once; arena-cell (Tensor<ArenaTensor<complex>>) conj(A*B) compiles (ArenaTensor gains operator*= for ComplexConjugate factors); nested-cell conjugation in the finalization is in place; Tensor::gemm rejects a ComplexConjugate alpha at compile time again; an integer-literal scale (2 * conj(a*b)) compiles.
• Tests (tot_conj suite in tests/tot_expressions.cpp, real and complex rows): conj(a), permuted conj(a), conj(a) * b, a * conj(b), conj(a * b) (outer contraction with inner outer product), conj(a) * b with an inner contraction, and the mixed-type ToT × real-plain product, all against explicit references.

Testing

tot_expressions (incl. the extended tot_conj: scaled/negated/integer-scaled conj, inner contraction and inner Hadamard, complex plain operand, multi-tile, fused mode, arena cells, interleaved real-gemm path), einsum, arena_strided_gemm: tot_expressions 362/362 cases (3486 assertions; tot_conj 22 cases / 808 assertions), einsum 59/59 (516), arena_strided_gemm 60/60 (1869), all against current master (post-#573, #581).

Not addressed

• Mixed real-ToT × complex-ToT products (the strided kernels require one element type for all three operands).
• Two review notes for a follow-up: the permuted conj finalization on owning nested tiles materializes the tile twice (conjugated copy, then permute), and the interleaved-gemm preamble is duplicated between the tot_x_t and t_x_tot blocks of Tensor::gemm.

[Copilot]

🟡 Changes recommended

The new interleaved real-GEMM fast path relies on reinterpret_cast assumptions about std::complex<T> representation and needs compile-time layout/alignment guards to avoid potential UB.

Once you've addressed the issues Copilot identified, you can request another Copilot review.

Pull request overview

This PR extends the ToT (tensor-of-tensor) contraction machinery to correctly support conj() around nested-tile contractions and adds a fast path for complex-ToT × real-plain tensor products by running a single real GEMM on an interleaved (re,im) view.

Changes:

• Introduces detail::elem_factor<Numeric>(factor) to extract the per-element contraction multiplier (numeric factor, or 1 for ComplexConjugate<...> factors) and wires it into ToT per-cell ops and Tensor::gemm alpha handling.
• Refactors ContractReduce to share nested-tile accumulation and partial-result reduction logic across the primary template and ComplexConjugate<...> specializations (removing the previous nested-tile abort() stub).
• Adds a mixed-type ToT × real-plain product optimization by viewing complex slabs as interleaved real matrices (doubling inner extent) and adds a focused tot_conj test suite covering these cases.

File summaries

File	Description
tests/tot_expressions.cpp	Adds tot_conj tests covering unary conj, conj around products, inner+outer contractions, and complex-ToT × real-plain products.
src/TiledArray/tile_op/contract_reduce.h	Shares nested accumulation + reduction logic and enables nested ToT accumulation for ComplexConjugate<...> factors.
src/TiledArray/tensor/tensor.h	Uses elem_factor for GEMM alpha and adds the interleaved real-GEMM fast path for complex-ToT × real-plain tensor scale products.
src/TiledArray/tensor/complex.h	Adds is_complex_conjugate_v and elem_factor helper for contraction-factor handling.
src/TiledArray/expressions/cont_engine.h	Switches ToT per-cell contraction factor application to elem_factor so ComplexConjugate<...> factors no longer require element-type casts.

Review details

Suppressed comments (1)

src/TiledArray/tensor/tensor.h:3531

• Same interleaved real-GEMM reinterpret_cast is used in the mirrored T * ToT scale path. Add the same compile-time layout/alignment guard for std::complex here as well, so the reinterpret_cast to Ur* is only used when the representation assumptions hold.

      constexpr bool same_type = std::is_same_v<Ur, Real>;
      constexpr bool interleaved =
          !same_type && std::is_same_v<std::complex<Ur>, Real>;
      constexpr integer cw = interleaved ? 2 : 1;  // reals per inner element
      if constexpr (same_type || interleaved) {

• Files reviewed: 5/5 changed files
• Comments generated: 1
• Review effort level: Lite

---

💡 Add a code-review agent skill or configure MCP servers for context-aware, tailored reviews. Learn more in the docs.

[evaleev]

Reviewed f399b490f..4fd5dcaa8 on top of #573. I built each changed header into probe programs (linked against a debug libtiledarray) and ran them, so the findings below are confirmed by execution rather than by reading.

Copilot's comment: already addressed

Copilot's single finding (guard the reinterpret_cast<Vr*> on std::complex<Vr> layout) is answered by 4fd5dcaa8, which adds the static_assert at both sites — including the mirrored t x ToT path at tensor.h:3729 that Copilot listed as a suppressed comment. Nothing left to do there.

One nit on the guard: alignof(Real) == alignof(Vr) is stricter than the cast needs (over-alignment of complex<double> is harmless), and because it is a static_assert rather than a term in the interleaved predicate, an ABI that over-aligns complex would make TA fail to compile rather than fall back to the AXPY loop. Folding the condition into constexpr bool interleaved degrades gracefully instead.

Blocking: two silent wrong answers, both new

Both cases were compile errors before this PR, so the PR turns hard failures into quiet incorrect results.

1. A ComplexConjugate<Scalar> prefactor is dropped on ToT contractions

complex.h:310 + cont_engine.h:435,498.

elem_factor returns Numeric(1) for any ComplexConjugate<...>, deferring conjugation and scale to ContractReduce's finalization. But the nested branch builds the outer op as op_type(..., scalar_type(1), ...), so finalization runs conj_to(temp, factor().factor()) with factor() == 1. The magnitude is applied nowhere.

On Tensor<Tensor<complex<double>>>:

2.0 * conj(a("i,j;a") * b("j,k;b"))   got (-0.5,-0)  ref (-1,-0)     # factor 2 gone
    -conj(a("i,j;a") * b("j,k;b"))    got (-0.5,-0)  ref ( 0.5, 0)   # sign gone

-conj(A*B) is an ordinary expression, not a corner case. Either propagate factor_ into the op's scalar_type for the nested branch, or have elem_factor return the ComplexConjugate<Scalar> magnitude.

2. conj(ToT * complex-plain) double-conjugates the plain operand

cont_engine.h:2006.

fallback_op bakes the raw factor_ in per element (right * factor => conj(right) for ComplexConjugate<void>), and finalization then conjugates the whole result again. The fused arena lambda is factor-free and correct, so the two paths now disagree. Reproduced with an inner permutation, which makes make_contraction_arena_plan return nullopt and select fallback_op:

c("i,k;b,a") = conj(a("i,j;a,b") * t("j,k"));   // t complex plain
got (-0.5,-7.1)  ref (-3.5,-5.9)

(-0.5,-7.1) is exactly sum_j conj(a_j) * t_j — t never gets conjugated. make_fused_hadamard_scaled_lambda(this->factor_) at cont_engine.h:1257/1899 has the same shape and was also left on the raw factor.

The stated goal is not reached for the arena tile type

conj(a("i,j;x") * b("j,k;y")) still does not compile for Tensor<ArenaTensor<complex<double>>> — the production tile type this series exists to optimize (per #573's commit message: MPQC Kramers PNS-CCD). The static_cast<T>(factor) pattern this PR replaced with elem_factor in Tensor::gemm survives in the arena mirror at arena_tensor.h:599, reached via arena_einsum.h:1110 <- fused_contraction_inplace <- init_inner_tile_op. Everything else on arena tiles works:

expression	owning Tensor<Tensor<cplx>>	Tensor<ArenaTensor<cplx>>
conj(a)	works	works
conj(a) * b, a * conj(b)	works	works
conj(a * b)	works	compile error

Not a regression (it fails on master too), but it means the first half of the PR title holds only for owning tiles.

Test coverage gaps

• The interleaved real-GEMM kernel ships untested. Both blocks gate on is_tensor_view_v<U> && is_tensor_view_v<value_type>, but conj_test_params uses owning Tensor<Tensor<complex<double>>> (the file comment says so). tot_times_real_plain therefore validates only the per-cell AXPY fallback. I wrote the missing case — arena complex ToT x real plain, checked against an owning reference — and it is correct: 16 strided GEMM runs, 0 fallbacks, max error 0 (verified via TA_GEMM_TIMING=1). The kernel is fine; it just has no test, and that case is worth adding.
• No test covers a scaled or negated conj on a ToT (finding 1), or conj combined with view/arena inner cells — reduce_results now routes those through arena_tot_add_to for the ComplexConjugate specializations, which is new and unexercised.

Minor

• tensor.h:3552,3729: the reinterpret_cast<const Vr*>(lc0[0].data()) is applied in the same_type branch too, and interleaved is decided from V and the result's inner scalar — the left operand's inner scalar (typename U::value_type) is never checked. I could not construct a reachable mismatch (the mixed real-ToT -> complex-result assignment does not compile), so this is hardening rather than a live bug: adding std::is_same_v<typename U::value_type, Real> to the predicate restores what the natural-typed pointer used to guarantee.
• contract_reduce.h:491,620: result_type changed from decltype(gemm(declval<Left>(), declval<Right>(), 1, helper)) to Result. For homogeneous cases these coincide and the change looks right (it aligns with ContractReduceBase), but it is a public-header behavior change worth calling out per the repo's PR-review guidance.
• contract_reduce.h:578,707: using TiledArray::empty; now appears twice (function scope + else block). Harmless; drop the inner one.
• Genuinely good: sharing accumulate_nested / reduce_results across the primary template and both specializations removes the abort() stub and keeps the three in sync by construction.

Verdict

Findings 1 and 2 need fixing before merge — both produce wrong numbers with no diagnostic on expressions a user would naturally write. The interleaved GEMM is correct but needs the arena test that actually exercises it.

[evaleev]

Requesting changes on the two blockers from the detailed review above. Both were compile errors before this PR and are silent wrong answers after it, on expressions a user would naturally write:

1. ComplexConjugate<Scalar> prefactor dropped on ToT contractions (complex.h:310, cont_engine.h:435,498) — elem_factor returns 1 for any ComplexConjugate<...>, but the nested branch builds the op with scalar_type(1), so finalization conjugates with factor 1 and the magnitude is applied nowhere. -conj(a("i,j;a") * b("j,k;b")) returns the unnegated result.

2. conj(ToT * complex-plain) double-conjugates the plain operand (cont_engine.h:2006) — fallback_op applies the raw factor_ per element and finalization conjugates the whole result. With an inner permutation (arena plan bails to nullopt) the result is sum_j conj(a_j) * t_j instead of conj(sum_j a_j * t_j). Same shape at cont_engine.h:1257/1899.

The interleaved real-GEMM path is correct as written — it just needs the arena test that actually exercises it, since the added tot_times_real_plain case uses owning tiles and never enters that code.

[Copilot]

🔵 Needs a closer look

It changes low-level contraction/GEMM paths and factor handling across multiple ToT execution modes, so it warrants final human review despite strong test additions.

Review details

• Files reviewed: 7/7 changed files
• Comments generated: 1
• Review effort level: Lite