Skip to main content

strided_einsum2/
contiguous.rs

1//! GEMM-ready operand types and preparation functions for contiguous data.
2//!
3//! These types encapsulate the logic for preparing strided operands for GEMM:
4//! checking fusability, copying to col-major buffers when needed, and managing
5//! the writeback for borrowed output operands.
6
7use crate::ScalarBase;
8use std::any::{Any, TypeId};
9use std::cell::RefCell;
10use std::collections::HashMap;
11use std::mem::MaybeUninit;
12use strided_kernel::MaybeSendSync;
13use strided_view::{RawStridedMut, RawStridedRef, StridedArray, StridedView, StridedViewMut};
14
15/// GEMM-ready input operand with contiguous data.
16pub struct ContiguousOperand<T: Copy + 'static> {
17    ptr: *const T,
18    row_stride: isize,
19    col_stride: isize,
20    batch_strides: Vec<isize>,
21    conj: bool,
22    /// Owns the buffer if a copy was made or input was consumed.
23    pub(crate) _buf: Option<StridedArray<T>>,
24    buf_is_pooled: bool,
25}
26
27/// GEMM-ready output operand with contiguous data.
28pub struct ContiguousOperandMut<T: Copy + 'static> {
29    ptr: *mut T,
30    row_stride: isize,
31    col_stride: isize,
32    batch_strides: Vec<isize>,
33    /// Whether the caller must copy the buffer back to the original destination
34    /// after GEMM completes (true only for borrowed non-contiguous C).
35    needs_writeback: bool,
36    /// Owns the buffer if a copy was made.
37    pub(crate) _buf: Option<StridedArray<T>>,
38    buf_is_pooled: bool,
39}
40
41/// Overwrite-only C operand. The destination borrow is retained until the
42/// provider has returned successfully; no initialized C view is constructed
43/// while the provider is running.
44#[allow(dead_code)]
45pub(crate) struct UninitContiguousOperand<'a, 'b, T: Copy + 'static> {
46    destination: &'a mut RawStridedMut<'b, MaybeUninit<T>>,
47    ptr: *mut MaybeUninit<T>,
48    row_stride: isize,
49    col_stride: isize,
50    batch_strides: Vec<isize>,
51    temp: Option<StridedArray<MaybeUninit<T>>>,
52    writeback: Option<strided_kernel::CopyPlan>,
53}
54
55thread_local! {
56    static BUFFER_POOL: RefCell<HashMap<TypeId, Box<dyn Any>>> = RefCell::new(HashMap::new());
57}
58
59const MAX_POOL_PER_TYPE: usize = 16;
60const MAX_POOLED_BYTES: usize = 64 * 1024 * 1024;
61
62fn take_pooled_vec_uninit<T: Copy + 'static>(len: usize) -> Vec<T> {
63    BUFFER_POOL.with(|pool| {
64        let mut pool = pool.borrow_mut();
65        let entry = pool
66            .entry(TypeId::of::<T>())
67            .or_insert_with(|| Box::new(Vec::<Vec<T>>::new()));
68        let vecs = entry
69            .downcast_mut::<Vec<Vec<T>>>()
70            .expect("buffer pool type mismatch");
71
72        let mut best_idx = None;
73        let mut best_cap = usize::MAX;
74        for (idx, v) in vecs.iter().enumerate() {
75            let cap = v.capacity();
76            if cap >= len && cap < best_cap {
77                best_idx = Some(idx);
78                best_cap = cap;
79            }
80        }
81
82        let mut data = best_idx
83            .map(|idx| vecs.swap_remove(idx))
84            .unwrap_or_else(|| Vec::with_capacity(len));
85        if data.capacity() < len {
86            data.reserve(len - data.capacity());
87        }
88        unsafe { data.set_len(len) };
89        data
90    })
91}
92
93fn return_pooled_vec<T: Copy + 'static>(mut data: Vec<T>) {
94    let bytes = data.capacity().saturating_mul(std::mem::size_of::<T>());
95    if bytes == 0 || bytes > MAX_POOLED_BYTES {
96        return;
97    }
98    data.clear();
99    BUFFER_POOL.with(|pool| {
100        let mut pool = pool.borrow_mut();
101        let entry = pool
102            .entry(TypeId::of::<T>())
103            .or_insert_with(|| Box::new(Vec::<Vec<T>>::new()));
104        let vecs = entry
105            .downcast_mut::<Vec<Vec<T>>>()
106            .expect("buffer pool type mismatch");
107        if vecs.len() >= MAX_POOL_PER_TYPE {
108            if let Some((min_idx, min_cap)) = vecs
109                .iter()
110                .enumerate()
111                .map(|(i, v)| (i, v.capacity()))
112                .min_by_key(|(_, cap)| *cap)
113            {
114                if min_cap < data.capacity() {
115                    vecs.swap_remove(min_idx);
116                    vecs.push(data);
117                }
118            }
119        } else {
120            vecs.push(data);
121        }
122    });
123}
124
125fn alloc_col_major_uninit_with_pool<T: Copy + 'static>(
126    dims: &[usize],
127) -> strided_view::Result<(StridedArray<T>, bool)> {
128    let total = dims
129        .iter()
130        .try_fold(1usize, |total, &dim| total.checked_mul(dim))
131        .ok_or(strided_view::StridedError::OffsetOverflow)?
132        .max(1);
133    let bytes = total.saturating_mul(std::mem::size_of::<T>());
134    if bytes == 0 || bytes > MAX_POOLED_BYTES {
135        return Ok((alloc_col_major_uninit(dims)?, false));
136    }
137    let data = take_pooled_vec_uninit::<T>(total);
138    let arr = unsafe { StridedArray::col_major_from_buffer_uninit(data, dims) };
139    Ok((arr, true))
140}
141
142/// Allocate a col-major buffer, optionally reusing from the thread-local pool.
143fn alloc_maybe_pooled<T: Copy + 'static>(
144    dims: &[usize],
145    use_pool: bool,
146) -> strided_view::Result<(StridedArray<T>, bool)> {
147    if use_pool {
148        alloc_col_major_uninit_with_pool(dims)
149    } else {
150        Ok((alloc_col_major_uninit(dims)?, false))
151    }
152}
153
154#[cfg(test)]
155fn pooled_count_for_type<T: 'static>() -> usize {
156    BUFFER_POOL.with(|pool| {
157        let mut pool = pool.borrow_mut();
158        let Some(entry) = pool.get_mut(&TypeId::of::<T>()) else {
159            return 0;
160        };
161        entry
162            .downcast_mut::<Vec<Vec<T>>>()
163            .map_or(0, |vecs| vecs.len())
164    })
165}
166
167impl<T: Copy + 'static> ContiguousOperand<T> {
168    /// Raw const pointer to the operand data at the base offset.
169    #[inline]
170    pub fn ptr(&self) -> *const T {
171        self.ptr
172    }
173
174    /// Row (lo-group) stride for the fused 2D matrix.
175    #[inline]
176    pub fn row_stride(&self) -> isize {
177        self.row_stride
178    }
179
180    /// Column (sum/ro-group) stride for the fused 2D matrix.
181    #[inline]
182    pub fn col_stride(&self) -> isize {
183        self.col_stride
184    }
185
186    /// Batch dimension strides.
187    #[inline]
188    pub fn batch_strides(&self) -> &[isize] {
189        &self.batch_strides
190    }
191
192    /// Whether this operand requires conjugation.
193    #[inline]
194    pub fn conj(&self) -> bool {
195        self.conj
196    }
197
198    /// Returns `true` if this operand owns a buffer (copy was made or ownership transferred).
199    #[cfg(test)]
200    #[inline]
201    pub(crate) fn has_buf(&self) -> bool {
202        self._buf.is_some()
203    }
204}
205
206impl<T: Copy + 'static> ContiguousOperandMut<T> {
207    /// Raw mutable pointer to the operand data at the base offset.
208    #[inline]
209    pub fn ptr(&self) -> *mut T {
210        self.ptr
211    }
212
213    /// Row (lo-group) stride for the fused 2D matrix.
214    #[inline]
215    pub fn row_stride(&self) -> isize {
216        self.row_stride
217    }
218
219    /// Column (ro-group) stride for the fused 2D matrix.
220    #[inline]
221    pub fn col_stride(&self) -> isize {
222        self.col_stride
223    }
224
225    /// Batch dimension strides.
226    #[inline]
227    pub fn batch_strides(&self) -> &[isize] {
228        &self.batch_strides
229    }
230
231    /// Returns `true` if this operand owns a buffer (copy was made).
232    #[cfg(test)]
233    #[inline]
234    pub(crate) fn has_buf(&self) -> bool {
235        self._buf.is_some()
236    }
237
238    /// Returns `true` if the caller must copy the buffer back to the original
239    /// destination after GEMM completes.
240    #[cfg(test)]
241    #[inline]
242    pub(crate) fn needs_writeback(&self) -> bool {
243        self.needs_writeback
244    }
245}
246
247impl<T: Copy + Send + Sync> ContiguousOperandMut<T> {
248    /// After GEMM: copy the internal buffer back to `dest` if needed.
249    ///
250    /// This is a no-op when the GEMM wrote directly to the destination
251    /// (contiguous case or owned output).
252    pub fn finalize_into(self, dest: &mut StridedViewMut<T>) -> crate::Result<()> {
253        if self.needs_writeback {
254            if let Some(ref buf) = self._buf {
255                strided_perm::copy_into(dest, &buf.view())?;
256            }
257        }
258        Ok(())
259    }
260
261    /// After GEMM: copy the internal buffer back to a raw destination if needed.
262    ///
263    /// This mirrors [`Self::finalize_into`] for prepared replay paths that use
264    /// borrowed raw layout metadata instead of owned-metadata views.
265    pub fn finalize_raw_into(self, dest: &mut RawStridedMut<'_, T>) -> crate::Result<()> {
266        if self.needs_writeback {
267            if let Some(ref buf) = self._buf {
268                let mut dest_view = dest.as_view_mut();
269                strided_perm::copy_into(&mut dest_view, &buf.view())?;
270            }
271        }
272        Ok(())
273    }
274}
275
276impl<T: Copy + 'static> Drop for ContiguousOperand<T> {
277    fn drop(&mut self) {
278        if self.buf_is_pooled {
279            if let Some(arr) = self._buf.take() {
280                return_pooled_vec(arr.into_data());
281            }
282        }
283    }
284}
285
286impl<T: Copy + 'static> Drop for ContiguousOperandMut<T> {
287    fn drop(&mut self) {
288        if self.buf_is_pooled {
289            if let Some(arr) = self._buf.take() {
290                return_pooled_vec(arr.into_data());
291            }
292        }
293    }
294}
295
296/// Result of checking whether dimension groups are contiguous enough for GEMM.
297struct ContiguityCheck {
298    fused_g1: Option<(usize, isize)>,
299    fused_g2: Option<(usize, isize)>,
300    needs_copy: bool,
301}
302
303/// Try to fuse a GEMM dimension group while preserving canonical col-major
304/// logical order. GEMM sees each group as a single 1D index, so independently
305/// fusing row-major and col-major groups would pair different logical indices.
306fn try_fuse_col_major_group(dims: &[usize], strides: &[isize]) -> Option<(usize, isize)> {
307    if dims.len() != strides.len() {
308        return None;
309    }
310    let total = dims
311        .iter()
312        .try_fold(1usize, |acc, &dim| acc.checked_mul(dim))?;
313    if dims.is_empty() {
314        return Some((1, 0));
315    }
316
317    let mut base_stride = None;
318    let mut expected_stride = None;
319    for (&dim, &stride) in dims.iter().zip(strides.iter()) {
320        if dim <= 1 {
321            continue;
322        }
323        if stride == 0 {
324            return None;
325        }
326        if let Some(expected) = expected_stride {
327            if stride != expected {
328                return None;
329            }
330        } else {
331            base_stride = Some(stride);
332        }
333        let dim = isize::try_from(dim).ok()?;
334        expected_stride = Some(stride.checked_mul(dim)?);
335    }
336
337    let stride = base_stride.unwrap_or_else(|| {
338        strides
339            .iter()
340            .copied()
341            .min_by_key(|stride| stride.unsigned_abs())
342            .unwrap_or(0)
343    });
344    Some((total, stride))
345}
346
347/// Check if two dimension groups are fusable (contiguous) for GEMM.
348///
349/// The fused logical dimension order must match the canonical axis order used
350/// by the plan. When `requires_unit_stride` is true (e.g., CBLAS backend), also
351/// checks that at least one of the fused strides is 0 or 1.
352fn check_contiguity(
353    group1_dims: &[usize],
354    group1_strides: &[isize],
355    group2_dims: &[usize],
356    group2_strides: &[isize],
357    requires_unit_stride: bool,
358) -> ContiguityCheck {
359    let fused_g1 = try_fuse_col_major_group(group1_dims, group1_strides);
360    let fused_g2 = try_fuse_col_major_group(group2_dims, group2_strides);
361
362    let mut needs_copy = fused_g1.is_none() || fused_g2.is_none();
363
364    if requires_unit_stride && !needs_copy {
365        let (_, rs) = fused_g1.unwrap();
366        let (_, cs) = fused_g2.unwrap();
367        if rs != 0 && rs != 1 && cs != 0 && cs != 1 {
368            needs_copy = true;
369        }
370    }
371
372    ContiguityCheck {
373        fused_g1,
374        fused_g2,
375        needs_copy,
376    }
377}
378
379/// Compute col-major layout parameters from a freshly-allocated col-major buffer.
380///
381/// Returns `(row_stride, col_stride, batch_strides)`.
382fn col_major_layout(
383    buf: &StridedArray<impl Copy>,
384    n_group1: usize,
385    n_inner: usize,
386) -> (isize, isize, Vec<isize>) {
387    let m: usize = buf.dims()[..n_group1].iter().product::<usize>().max(1);
388    let row_stride = if m == 0 { 0 } else { 1isize };
389    let col_stride = m as isize;
390    let batch_strides = buf.strides()[n_inner..].to_vec();
391    (row_stride, col_stride, batch_strides)
392}
393
394/// Allocate a column-major StridedArray with uninitialized data.
395///
396/// With batch-last canonical order `[inner..., batch...]`, pure column-major
397/// naturally gives batch dims the largest strides — each batch slice is a
398/// contiguous column-major matrix.
399pub(crate) fn alloc_col_major_uninit<T: Copy>(
400    dims: &[usize],
401) -> strided_view::Result<StridedArray<T>> {
402    let total = dims
403        .iter()
404        .try_fold(1usize, |total, &dim| total.checked_mul(dim))
405        .ok_or(strided_view::StridedError::OffsetOverflow)?
406        .max(1);
407    // SAFETY: `T: Copy` guarantees no drop glue, so leaving elements
408    // uninitialised is safe. Each caller must establish initialization before
409    // exposing the array as initialized; uninitialized-output callers use
410    // `T = MaybeUninit<U>` and finalize only after complete overwrite.
411    let mut data = Vec::with_capacity(total);
412    unsafe { data.set_len(total) };
413
414    // Pure column-major: stride 1 for first dim, each subsequent dim
415    // has stride = previous stride * previous dim size.
416    let mut strides = vec![0isize; dims.len()];
417    if !dims.is_empty() {
418        strides[0] = 1;
419        for i in 1..dims.len() {
420            strides[i] = strides[i - 1] * dims[i - 1] as isize;
421        }
422    }
423
424    Ok(StridedArray::from_parts(data, dims, &strides, 0)?)
425}
426
427/// Prepare a borrowed input view for GEMM.
428///
429/// Expects batch-last canonical order: `[group1..., group2..., batch...]`.
430/// Checks if the two inner dimension groups are fusable.
431/// If not, copies to a contiguous col-major buffer.
432///
433/// - `requires_unit_stride`: backend needs at least one unit stride (e.g. CBLAS).
434/// - `use_pool`: reuse thread-local buffers to avoid repeated allocation.
435/// - `materialize_conj_fn`: when `Some(f)` and `conj == true`, applies `f` to each
436///   element during copy (for backends that cannot pass conj flags to GEMM).
437pub fn prepare_input_view<T: ScalarBase + 'static>(
438    view: &StridedView<T>,
439    n_group1: usize,
440    n_group2: usize,
441    conj: bool,
442    requires_unit_stride: bool,
443    use_pool: bool,
444    materialize_conj_fn: Option<fn(T) -> T>,
445) -> crate::Result<ContiguousOperand<T>> {
446    let dims = view.dims();
447    let strides = view.strides();
448    let n_inner = n_group1 + n_group2;
449
450    // For backends that cannot pass conjugation flags to GEMM (e.g., CBLAS),
451    // materialize conj into the data before the GEMM call.
452    if let Some(conj_fn) = materialize_conj_fn {
453        if conj {
454            let (mut buf, buf_is_pooled) = alloc_maybe_pooled(dims, use_pool)?;
455            strided_kernel::map_into(&mut buf.view_mut(), view, conj_fn)?;
456            let ptr = buf.view().ptr();
457            let (row_stride, col_stride, batch_strides) = col_major_layout(&buf, n_group1, n_inner);
458            return Ok(ContiguousOperand {
459                ptr,
460                row_stride,
461                col_stride,
462                batch_strides,
463                conj: false,
464                _buf: Some(buf),
465                buf_is_pooled,
466            });
467        }
468    }
469
470    let check = check_contiguity(
471        &dims[..n_group1],
472        &strides[..n_group1],
473        &dims[n_group1..n_inner],
474        &strides[n_group1..n_inner],
475        requires_unit_stride,
476    );
477
478    if check.needs_copy {
479        let (mut buf, buf_is_pooled) = alloc_maybe_pooled(dims, use_pool)?;
480        strided_kernel::copy_into_col_major(&mut buf.view_mut(), view)?;
481        let ptr = buf.view().ptr();
482        let (row_stride, col_stride, batch_strides) = col_major_layout(&buf, n_group1, n_inner);
483        Ok(ContiguousOperand {
484            ptr,
485            row_stride,
486            col_stride,
487            batch_strides,
488            conj,
489            _buf: Some(buf),
490            buf_is_pooled,
491        })
492    } else {
493        let (_, rs) = check.fused_g1.unwrap();
494        let (_, cs) = check.fused_g2.unwrap();
495        Ok(ContiguousOperand {
496            ptr: view.ptr(),
497            row_stride: rs,
498            col_stride: cs,
499            batch_strides: strides[n_inner..].to_vec(),
500            conj,
501            _buf: None,
502            buf_is_pooled: false,
503        })
504    }
505}
506
507/// Prepare a borrowed raw input layout for GEMM.
508///
509/// This is the raw-layout counterpart to [`prepare_input_view`]. Direct
510/// fusable layouts do not construct `StridedView` wrappers; a view is created
511/// only when a copy/materialization path needs to call a generic strided kernel.
512pub fn prepare_input_raw<T: ScalarBase + 'static>(
513    view: &RawStridedRef<'_, T>,
514    n_group1: usize,
515    n_group2: usize,
516    conj: bool,
517    requires_unit_stride: bool,
518    use_pool: bool,
519    materialize_conj_fn: Option<fn(T) -> T>,
520) -> crate::Result<ContiguousOperand<T>> {
521    let dims = view.dims();
522    let strides = view.strides();
523    let n_inner = n_group1 + n_group2;
524
525    if let Some(conj_fn) = materialize_conj_fn {
526        if conj {
527            let (mut buf, buf_is_pooled) = alloc_maybe_pooled(dims, use_pool)?;
528            strided_kernel::map_into(&mut buf.view_mut(), &view.as_view(), conj_fn)?;
529            let ptr = buf.view().ptr();
530            let (row_stride, col_stride, batch_strides) = col_major_layout(&buf, n_group1, n_inner);
531            return Ok(ContiguousOperand {
532                ptr,
533                row_stride,
534                col_stride,
535                batch_strides,
536                conj: false,
537                _buf: Some(buf),
538                buf_is_pooled,
539            });
540        }
541    }
542
543    let check = check_contiguity(
544        &dims[..n_group1],
545        &strides[..n_group1],
546        &dims[n_group1..n_inner],
547        &strides[n_group1..n_inner],
548        requires_unit_stride,
549    );
550
551    if check.needs_copy {
552        let (mut buf, buf_is_pooled) = alloc_maybe_pooled(dims, use_pool)?;
553        strided_kernel::copy_into_col_major(&mut buf.view_mut(), &view.as_view())?;
554        let ptr = buf.view().ptr();
555        let (row_stride, col_stride, batch_strides) = col_major_layout(&buf, n_group1, n_inner);
556        Ok(ContiguousOperand {
557            ptr,
558            row_stride,
559            col_stride,
560            batch_strides,
561            conj,
562            _buf: Some(buf),
563            buf_is_pooled,
564        })
565    } else {
566        let (_, rs) = check.fused_g1.unwrap();
567        let (_, cs) = check.fused_g2.unwrap();
568        Ok(ContiguousOperand {
569            ptr: view.ptr(),
570            row_stride: rs,
571            col_stride: cs,
572            batch_strides: strides[n_inner..].to_vec(),
573            conj,
574            _buf: None,
575            buf_is_pooled: false,
576        })
577    }
578}
579
580/// Prepare an owned input array for GEMM.
581///
582/// Expects batch-last canonical order: `[group1..., group2..., batch...]`.
583/// If already contiguous after dimension grouping, transfers ownership without copying.
584/// Otherwise, copies to a new col-major buffer.
585///
586/// Parameters are the same as [`prepare_input_view`].
587pub fn prepare_input_owned<T: ScalarBase + 'static>(
588    arr: StridedArray<T>,
589    n_group1: usize,
590    n_group2: usize,
591    conj: bool,
592    requires_unit_stride: bool,
593    use_pool: bool,
594    materialize_conj_fn: Option<fn(T) -> T>,
595) -> crate::Result<ContiguousOperand<T>> {
596    let dims = arr.dims().to_vec();
597    let strides = arr.strides().to_vec();
598    let n_inner = n_group1 + n_group2;
599
600    // For backends that cannot pass conjugation flags to GEMM,
601    // materialize conj into the data before the GEMM call.
602    if let Some(conj_fn) = materialize_conj_fn {
603        if conj {
604            let (mut buf, buf_is_pooled) = alloc_maybe_pooled(&dims, use_pool)?;
605            strided_kernel::map_into(&mut buf.view_mut(), &arr.view(), conj_fn)?;
606            let ptr = buf.view().ptr();
607            let (row_stride, col_stride, batch_strides) = col_major_layout(&buf, n_group1, n_inner);
608            return Ok(ContiguousOperand {
609                ptr,
610                row_stride,
611                col_stride,
612                batch_strides,
613                conj: false,
614                _buf: Some(buf),
615                buf_is_pooled,
616            });
617        }
618    }
619
620    let check = check_contiguity(
621        &dims[..n_group1],
622        &strides[..n_group1],
623        &dims[n_group1..n_inner],
624        &strides[n_group1..n_inner],
625        requires_unit_stride,
626    );
627
628    if check.needs_copy {
629        let (mut buf, buf_is_pooled) = alloc_maybe_pooled(&dims, use_pool)?;
630        strided_kernel::copy_into_col_major(&mut buf.view_mut(), &arr.view())?;
631        let ptr = buf.view().ptr();
632        let (row_stride, col_stride, batch_strides) = col_major_layout(&buf, n_group1, n_inner);
633        Ok(ContiguousOperand {
634            ptr,
635            row_stride,
636            col_stride,
637            batch_strides,
638            conj,
639            _buf: Some(buf),
640            buf_is_pooled,
641        })
642    } else {
643        let (_, rs) = check.fused_g1.unwrap();
644        let (_, cs) = check.fused_g2.unwrap();
645        let ptr = arr.view().ptr();
646        Ok(ContiguousOperand {
647            ptr,
648            row_stride: rs,
649            col_stride: cs,
650            batch_strides: strides[n_inner..].to_vec(),
651            conj,
652            _buf: Some(arr),
653            buf_is_pooled: false,
654        })
655    }
656}
657
658/// Prepare a borrowed mutable output view for GEMM.
659///
660/// Expects batch-last canonical order: `[group1..., group2..., batch...]`.
661/// Checks if the two inner dimension groups (lo, ro) are fusable.
662/// If not, allocates a col-major buffer and copies the existing data into it
663/// when `beta` is non-zero (so the GEMM accumulation is correct).
664///
665/// After GEMM, call [`ContiguousOperandMut::finalize_into`] with the original
666/// view to copy results back if needed.
667///
668/// # Safety contract
669///
670/// When inner dims are fusable (no copy needed), the returned `ContiguousOperandMut`
671/// holds a raw pointer into `view`'s data. The caller must ensure `view` outlives
672/// the returned operand and that no aliasing mutable references exist during GEMM.
673pub fn prepare_output_view<T: ScalarBase + 'static>(
674    view: &mut StridedViewMut<T>,
675    n_group1: usize,
676    n_group2: usize,
677    beta: T,
678    requires_unit_stride: bool,
679    use_pool: bool,
680) -> crate::Result<ContiguousOperandMut<T>> {
681    let dims = view.dims().to_vec();
682    let strides = view.strides().to_vec();
683    let n_inner = n_group1 + n_group2;
684
685    let check = check_contiguity(
686        &dims[..n_group1],
687        &strides[..n_group1],
688        &dims[n_group1..n_inner],
689        &strides[n_group1..n_inner],
690        requires_unit_stride,
691    );
692
693    if check.needs_copy {
694        let (mut buf, buf_is_pooled) = alloc_maybe_pooled(&dims, use_pool)?;
695        if beta != T::zero() {
696            strided_kernel::copy_into_col_major(&mut buf.view_mut(), &view.as_view())?;
697        }
698        let ptr = buf.view_mut().as_mut_ptr();
699        let (row_stride, col_stride, batch_strides) = col_major_layout(&buf, n_group1, n_inner);
700        Ok(ContiguousOperandMut {
701            ptr,
702            row_stride,
703            col_stride,
704            batch_strides,
705            needs_writeback: true,
706            _buf: Some(buf),
707            buf_is_pooled,
708        })
709    } else {
710        let (_, rs) = check.fused_g1.unwrap();
711        let (_, cs) = check.fused_g2.unwrap();
712        Ok(ContiguousOperandMut {
713            ptr: view.as_mut_ptr(),
714            row_stride: rs,
715            col_stride: cs,
716            batch_strides: strides[n_inner..].to_vec(),
717            needs_writeback: false,
718            _buf: None,
719            buf_is_pooled: false,
720        })
721    }
722}
723
724/// Prepare a borrowed raw output layout for GEMM.
725///
726/// This is the raw-layout counterpart to [`prepare_output_view`]. Direct
727/// fusable layouts do not construct `StridedViewMut` wrappers; a view is
728/// created only when a copy/writeback path needs a generic strided kernel.
729pub fn prepare_output_raw<T: ScalarBase + 'static>(
730    view: &mut RawStridedMut<'_, T>,
731    n_group1: usize,
732    n_group2: usize,
733    beta: T,
734    requires_unit_stride: bool,
735    use_pool: bool,
736) -> crate::Result<ContiguousOperandMut<T>> {
737    let dims = view.dims().to_vec();
738    let strides = view.strides().to_vec();
739    let n_inner = n_group1 + n_group2;
740
741    let check = check_contiguity(
742        &dims[..n_group1],
743        &strides[..n_group1],
744        &dims[n_group1..n_inner],
745        &strides[n_group1..n_inner],
746        requires_unit_stride,
747    );
748
749    if check.needs_copy {
750        let (mut buf, buf_is_pooled) = alloc_maybe_pooled(&dims, use_pool)?;
751        if beta != T::zero() {
752            strided_kernel::copy_into_col_major(&mut buf.view_mut(), &view.as_view())?;
753        }
754        let ptr = buf.view_mut().as_mut_ptr();
755        let (row_stride, col_stride, batch_strides) = col_major_layout(&buf, n_group1, n_inner);
756        Ok(ContiguousOperandMut {
757            ptr,
758            row_stride,
759            col_stride,
760            batch_strides,
761            needs_writeback: true,
762            _buf: Some(buf),
763            buf_is_pooled,
764        })
765    } else {
766        let (_, rs) = check.fused_g1.unwrap();
767        let (_, cs) = check.fused_g2.unwrap();
768        Ok(ContiguousOperandMut {
769            ptr: view.as_mut_ptr(),
770            row_stride: rs,
771            col_stride: cs,
772            batch_strides: strides[n_inner..].to_vec(),
773            needs_writeback: false,
774            _buf: None,
775            buf_is_pooled: false,
776        })
777    }
778}
779
780#[allow(dead_code)]
781impl<'a, 'b, T: Copy + MaybeSendSync + 'static> UninitContiguousOperand<'a, 'b, T> {
782    #[inline]
783    pub(crate) fn ptr(&self) -> *mut MaybeUninit<T> {
784        self.ptr
785    }
786
787    #[inline]
788    pub(crate) fn row_stride(&self) -> isize {
789        self.row_stride
790    }
791
792    #[inline]
793    pub(crate) fn col_stride(&self) -> isize {
794        self.col_stride
795    }
796
797    #[inline]
798    pub(crate) fn batch_strides(&self) -> &[isize] {
799        &self.batch_strides
800    }
801
802    /// Finalize the destination only after a successful overwrite provider.
803    pub(crate) fn finalize(self) -> crate::Result<()> {
804        let Self {
805            destination,
806            temp,
807            writeback,
808            ..
809        } = self;
810        let (Some(temp), Some(writeback)) = (temp, writeback) else {
811            return Ok(());
812        };
813        // The temporary is dense and the provider has returned success, so
814        // every element is initialized. This conversion is intentionally
815        // after provider success and is never used for direct C storage.
816        let dims = temp.dims().to_vec();
817        let strides = temp.strides().to_vec();
818        let data = temp.into_data();
819        let len = data.len();
820        let cap = data.capacity();
821        let ptr = data.as_ptr().cast_mut().cast::<T>();
822        std::mem::forget(data);
823        let initialized = unsafe {
824            StridedArray::from_parts(Vec::from_raw_parts(ptr, len, cap), &dims, &strides, 0)
825        }?;
826        let source = RawStridedRef::new(
827            initialized.data(),
828            initialized.dims(),
829            initialized.strides(),
830            initialized.view().offset(),
831        )?;
832        writeback.execute_uninit(destination, &source)?;
833        Ok(())
834    }
835}
836
837/// Prepare an overwrite-only C operand. A non-fusable destination receives a
838/// dense `MaybeUninit` temporary and a compiled uninitialized writeback plan.
839#[allow(dead_code)]
840pub(crate) fn prepare_output_raw_uninit<'a, 'b, T: ScalarBase + MaybeSendSync + 'static>(
841    destination: &'a mut RawStridedMut<'b, MaybeUninit<T>>,
842    n_group1: usize,
843    n_group2: usize,
844    requires_unit_stride: bool,
845) -> crate::Result<UninitContiguousOperand<'a, 'b, T>> {
846    let dims = destination.dims().to_vec();
847    let strides = destination.strides().to_vec();
848    let n_inner = n_group1 + n_group2;
849    let check = check_contiguity(
850        &dims[..n_group1],
851        &strides[..n_group1],
852        &dims[n_group1..n_inner],
853        &strides[n_group1..n_inner],
854        requires_unit_stride,
855    );
856    if !check.needs_copy {
857        let Some((_, row_stride)) = check.fused_g1 else {
858            return Err(strided_view::StridedError::PlanLayoutMismatch.into());
859        };
860        let Some((_, col_stride)) = check.fused_g2 else {
861            return Err(strided_view::StridedError::PlanLayoutMismatch.into());
862        };
863        let ptr = destination.as_mut_ptr();
864        return Ok(UninitContiguousOperand {
865            destination,
866            ptr: ptr.cast(),
867            row_stride,
868            col_stride,
869            batch_strides: strides[n_inner..].to_vec(),
870            temp: None,
871            writeback: None,
872        });
873    }
874
875    let temp = alloc_col_major_uninit::<MaybeUninit<T>>(&dims)?;
876    let ptr = temp.view().ptr().cast_mut();
877    let (row_stride, col_stride, batch_strides) = col_major_layout(&temp, n_group1, n_inner);
878    let writeback = strided_kernel::CopyPlan::compile(&dims, &strides, temp.strides())?;
879    Ok(UninitContiguousOperand {
880        destination,
881        ptr,
882        row_stride,
883        col_stride,
884        batch_strides,
885        temp: Some(temp),
886        writeback: Some(writeback),
887    })
888}
889
890#[cfg(test)]
891mod tests_generic_backend {
892    use super::*;
893    use crate::backend::{Backend, NaiveBackend};
894
895    #[test]
896    fn test_input_for_backend_contiguous() {
897        let a = StridedArray::<f64>::col_major(&[2, 3]);
898        let view = a.view();
899        let op = prepare_input_view(
900            &view,
901            1,
902            1,
903            false,
904            <NaiveBackend as Backend<f64>>::REQUIRES_UNIT_STRIDE,
905            false,
906            None,
907        )
908        .unwrap();
909        assert!(op._buf.is_none());
910        assert_eq!(op.row_stride(), 1);
911        assert_eq!(op.col_stride(), 2);
912        assert!(!op.conj());
913    }
914
915    #[test]
916    fn test_input_for_backend_non_contiguous() {
917        let data = vec![0.0f64; 100];
918        let a = StridedArray::<f64>::from_parts(data, &[2, 3, 4], &[20, 4, 1], 0).unwrap();
919        let view = a.view();
920        let op = prepare_input_view(
921            &view,
922            2,
923            1,
924            false,
925            <NaiveBackend as Backend<f64>>::REQUIRES_UNIT_STRIDE,
926            false,
927            None,
928        )
929        .unwrap();
930        assert!(op._buf.is_some());
931        assert_eq!(op.row_stride(), 1);
932        assert_eq!(op.col_stride(), 6);
933    }
934
935    #[test]
936    fn test_output_for_backend_contiguous() {
937        let mut c = StridedArray::<f64>::col_major(&[2, 3]);
938        let mut view = c.view_mut();
939        let op = prepare_output_view(
940            &mut view,
941            1,
942            1,
943            0.0,
944            <NaiveBackend as Backend<f64>>::REQUIRES_UNIT_STRIDE,
945            false,
946        )
947        .unwrap();
948        assert!(!op.needs_writeback);
949        assert!(op._buf.is_none());
950        assert_eq!(op.row_stride(), 1);
951        assert_eq!(op.col_stride(), 2);
952    }
953
954    #[test]
955    fn test_output_for_backend_non_contiguous_beta_zero() {
956        let data = vec![0.0f64; 100];
957        let mut c = StridedArray::<f64>::from_parts(data, &[2, 3, 4], &[20, 4, 1], 0).unwrap();
958        let mut view = c.view_mut();
959        let op = prepare_output_view(
960            &mut view,
961            2,
962            1,
963            0.0,
964            <NaiveBackend as Backend<f64>>::REQUIRES_UNIT_STRIDE,
965            false,
966        )
967        .unwrap();
968        assert!(op.needs_writeback);
969        assert!(op._buf.is_some());
970        assert_eq!(op.row_stride(), 1);
971        assert_eq!(op.col_stride(), 6);
972    }
973
974    #[test]
975    fn test_output_for_backend_non_contiguous_beta_nonzero_and_finalize() {
976        let mut data = vec![0.0f64; 30];
977        data[0] = 10.0;
978        data[1] = 20.0;
979        data[10] = 40.0;
980        let mut c = StridedArray::<f64>::from_parts(data, &[2, 3, 1], &[10, 1, 1], 0).unwrap();
981        let mut view = c.view_mut();
982        let op = prepare_output_view(
983            &mut view,
984            2,
985            1,
986            1.0,
987            <NaiveBackend as Backend<f64>>::REQUIRES_UNIT_STRIDE,
988            false,
989        )
990        .unwrap();
991        assert!(op.needs_writeback);
992        let buf = op._buf.as_ref().unwrap();
993        assert_eq!(buf.get(&[0, 0, 0]), 10.0);
994        assert_eq!(buf.get(&[0, 1, 0]), 20.0);
995        assert_eq!(buf.get(&[1, 0, 0]), 40.0);
996        op.finalize_into(&mut view).unwrap();
997    }
998}
999
1000#[cfg(test)]
1001mod tests {
1002    use super::*;
1003    use crate::backend::{ActiveBackend, Backend};
1004
1005    // Helper to construct prepare_input_view params matching the active backend.
1006    const UNIT_STRIDE: bool = <ActiveBackend as Backend<f64>>::REQUIRES_UNIT_STRIDE;
1007
1008    #[test]
1009    fn test_borrowed_contiguous_no_copy() {
1010        let a = StridedArray::<f64>::col_major(&[2, 3]);
1011        let view = a.view();
1012
1013        let op = prepare_input_view(&view, 1, 1, false, UNIT_STRIDE, true, None).unwrap();
1014
1015        assert!(!op.has_buf());
1016        assert_eq!(op.row_stride(), 1);
1017        assert_eq!(op.col_stride(), 2);
1018        assert!(!op.conj());
1019    }
1020
1021    #[test]
1022    fn test_borrowed_transposed_matrix_no_copy() {
1023        let data = vec![0.0f64; 6];
1024        let a_t = StridedArray::<f64>::from_parts(data, &[2, 3], &[3, 1], 0).unwrap();
1025        let view = a_t.view();
1026
1027        let op = prepare_input_view(&view, 1, 1, false, UNIT_STRIDE, true, None).unwrap();
1028
1029        assert!(!op.has_buf());
1030        assert_eq!(op.row_stride(), 3);
1031        assert_eq!(op.col_stride(), 1);
1032    }
1033
1034    #[test]
1035    fn test_borrowed_batched_transposed_matrix_no_copy() {
1036        let data = vec![0.0f64; 2 * 3 * 5];
1037        let a_t = StridedArray::<f64>::from_parts(data, &[2, 3, 5], &[3, 1, 6], 0).unwrap();
1038        let view = a_t.view();
1039
1040        let op = prepare_input_view(&view, 1, 1, false, UNIT_STRIDE, true, None).unwrap();
1041
1042        assert!(!op.has_buf());
1043        assert_eq!(op.row_stride(), 3);
1044        assert_eq!(op.col_stride(), 1);
1045        assert_eq!(op.batch_strides(), &[6]);
1046    }
1047
1048    #[test]
1049    fn test_borrowed_non_contiguous_copies() {
1050        let data = vec![0.0f64; 100];
1051        let a = StridedArray::<f64>::from_parts(data, &[2, 3, 4], &[20, 4, 1], 0).unwrap();
1052        let view = a.view();
1053
1054        let op = prepare_input_view(&view, 2, 1, false, UNIT_STRIDE, true, None).unwrap();
1055
1056        assert!(op.has_buf());
1057        assert_eq!(op.row_stride(), 1);
1058        assert_eq!(op.col_stride(), 6);
1059    }
1060
1061    #[test]
1062    fn test_owned_contiguous_no_copy() {
1063        let a = StridedArray::<f64>::col_major(&[2, 3]);
1064
1065        let op = prepare_input_owned(a, 1, 1, false, UNIT_STRIDE, true, None).unwrap();
1066
1067        assert!(op.has_buf());
1068        assert_eq!(op.row_stride(), 1);
1069        assert_eq!(op.col_stride(), 2);
1070    }
1071
1072    #[test]
1073    fn test_owned_non_contiguous_copies() {
1074        let data = vec![0.0f64; 100];
1075        let a = StridedArray::<f64>::from_parts(data, &[2, 3, 4], &[20, 4, 1], 0).unwrap();
1076
1077        let op = prepare_input_owned(a, 2, 1, false, UNIT_STRIDE, true, None).unwrap();
1078
1079        assert!(op.has_buf());
1080        assert_eq!(op.row_stride(), 1);
1081        assert_eq!(op.col_stride(), 6);
1082    }
1083
1084    #[test]
1085    fn test_output_view_contiguous() {
1086        let mut c = StridedArray::<f64>::col_major(&[2, 3]);
1087        let mut view = c.view_mut();
1088
1089        let op = prepare_output_view(&mut view, 1, 1, 0.0, UNIT_STRIDE, true).unwrap();
1090
1091        assert!(!op.needs_writeback());
1092        assert!(!op.has_buf());
1093        assert_eq!(op.row_stride(), 1);
1094        assert_eq!(op.col_stride(), 2);
1095    }
1096
1097    #[test]
1098    fn test_output_view_non_contiguous_beta_zero() {
1099        let data = vec![0.0f64; 100];
1100        let mut c = StridedArray::<f64>::from_parts(data, &[2, 3, 4], &[20, 4, 1], 0).unwrap();
1101        let mut view = c.view_mut();
1102
1103        let op = prepare_output_view(&mut view, 2, 1, 0.0, UNIT_STRIDE, true).unwrap();
1104
1105        assert!(op.needs_writeback());
1106        assert!(op.has_buf());
1107        assert_eq!(op.row_stride(), 1);
1108        assert_eq!(op.col_stride(), 6);
1109    }
1110
1111    #[test]
1112    fn test_output_view_non_contiguous_beta_nonzero_and_finalize() {
1113        let mut data = vec![0.0f64; 30];
1114        data[0] = 10.0;
1115        data[1] = 20.0;
1116        data[2] = 30.0;
1117        data[10] = 40.0;
1118        data[11] = 50.0;
1119        data[12] = 60.0;
1120        let mut c = StridedArray::<f64>::from_parts(data, &[2, 3, 1], &[10, 1, 1], 0).unwrap();
1121
1122        assert_eq!(c.get(&[0, 0, 0]), 10.0);
1123        assert_eq!(c.get(&[1, 1, 0]), 50.0);
1124
1125        let mut view = c.view_mut();
1126
1127        let mut op = prepare_output_view(&mut view, 2, 1, 1.0, UNIT_STRIDE, true).unwrap();
1128
1129        assert!(op.needs_writeback());
1130        assert!(op.has_buf());
1131
1132        let buf = op._buf.as_ref().unwrap();
1133        assert_eq!(buf.get(&[0, 0, 0]), 10.0);
1134        assert_eq!(buf.get(&[1, 1, 0]), 50.0);
1135
1136        {
1137            let result_data = vec![100.0f64; 6];
1138            let result =
1139                StridedArray::<f64>::from_parts(result_data, &[2, 3, 1], &[3, 1, 1], 0).unwrap();
1140            strided_kernel::copy_into(&mut op._buf.as_mut().unwrap().view_mut(), &result.view())
1141                .unwrap();
1142            op.ptr = op._buf.as_mut().unwrap().view_mut().as_mut_ptr();
1143        }
1144
1145        op.finalize_into(&mut view).unwrap();
1146
1147        assert_eq!(c.get(&[0, 0, 0]), 100.0);
1148        assert_eq!(c.get(&[0, 1, 0]), 100.0);
1149        assert_eq!(c.get(&[0, 2, 0]), 100.0);
1150        assert_eq!(c.get(&[1, 0, 0]), 100.0);
1151        assert_eq!(c.get(&[1, 1, 0]), 100.0);
1152        assert_eq!(c.get(&[1, 2, 0]), 100.0);
1153    }
1154
1155    #[test]
1156    fn test_prepare_input_view_temp_buffer_is_recycled() {
1157        let before = pooled_count_for_type::<f64>();
1158        let data = vec![0.0f64; 100];
1159        let a = StridedArray::<f64>::from_parts(data, &[2, 3, 4], &[20, 4, 1], 0).unwrap();
1160        let view = a.view();
1161
1162        {
1163            let op = prepare_input_view(&view, 2, 1, false, UNIT_STRIDE, true, None).unwrap();
1164            assert!(op.has_buf());
1165        }
1166
1167        let after = pooled_count_for_type::<f64>();
1168        assert!(after >= before.saturating_add(1));
1169    }
1170
1171    #[test]
1172    fn test_output_raw_uninit_direct_finalize() {
1173        let mut storage = vec![MaybeUninit::<f64>::uninit(); 6];
1174        let mut raw = RawStridedMut::new(&mut storage, &[2, 3], &[1, 2], 0).unwrap();
1175        let op = prepare_output_raw_uninit(&mut raw, 1, 1, false).unwrap();
1176
1177        assert!(op.temp.is_none());
1178        assert!(op.writeback.is_none());
1179        for j in 0..3 {
1180            for i in 0..2 {
1181                let offset = i as isize * op.row_stride() + j as isize * op.col_stride();
1182                // SAFETY: the direct layout is injective and every logical
1183                // destination element is written exactly once before finalize.
1184                unsafe {
1185                    op.ptr()
1186                        .offset(offset)
1187                        .write(MaybeUninit::new((10 + i + 2 * j) as f64));
1188                }
1189            }
1190        }
1191        op.finalize().unwrap();
1192        drop(raw);
1193
1194        let values: Vec<f64> = storage
1195            .into_iter()
1196            .map(|x| unsafe { x.assume_init() })
1197            .collect();
1198        assert_eq!(values, vec![10.0, 11.0, 12.0, 13.0, 14.0, 15.0]);
1199    }
1200
1201    #[test]
1202    fn test_output_raw_uninit_temp_finalize_writes_back() {
1203        let mut storage = vec![MaybeUninit::<f64>::uninit(); 30];
1204        let mut raw = RawStridedMut::new(&mut storage, &[2, 3, 1], &[10, 1, 1], 0).unwrap();
1205        let op = prepare_output_raw_uninit(&mut raw, 2, 1, true).unwrap();
1206
1207        assert!(op.temp.is_some());
1208        assert!(op.writeback.is_some());
1209        for i in 0..6 {
1210            let offset = i as isize * op.row_stride();
1211            // SAFETY: the dense temporary has one initialized value for every
1212            // logical output element before writeback.
1213            unsafe {
1214                op.ptr()
1215                    .offset(offset)
1216                    .write(MaybeUninit::new((20 + i) as f64));
1217            }
1218        }
1219        op.finalize().unwrap();
1220        drop(raw);
1221
1222        // SAFETY: finalize wrote every asserted destination element before
1223        // these reads.
1224        unsafe {
1225            assert_eq!(storage[0].assume_init_ref(), &20.0);
1226            assert_eq!(storage[1].assume_init_ref(), &22.0);
1227            assert_eq!(storage[2].assume_init_ref(), &24.0);
1228            assert_eq!(storage[10].assume_init_ref(), &21.0);
1229            assert_eq!(storage[11].assume_init_ref(), &23.0);
1230            assert_eq!(storage[12].assume_init_ref(), &25.0);
1231        }
1232    }
1233}