fft: apply plans on cursors (.home()), the complement of extents

A plan holds the shape (extents + precomputed tables); it is now applied on
a cursor -- base + strides, no shape of its own:

  plan(A.home());   // primitive: cursor gives base+strides, plan gives shape
  plan(A);          // sugar: checks A's extents == plan, then plan(A.home())

plan(cursor) rebuilds a strided view from the plan's extents and the cursor
(via layout_t's size = nelems/stride invariant) and drives the unchanged
orchestration; plan(array) delegates to it. Both entry points are SFINAE-
disjoint (cursor = base()/strides() but not multi-like). Returns the cursor
by value (safe when built from a temporary .home()).

This is the GPU-facing execution primitive: a device cursor (device pointer
base + strides) flows through plan(cursor) unchanged, so the memory-space
dispatch can hang off it (see fft.NOTES.md section 8). fft.hpp now includes
array_ref.hpp for layout_t/subarray and compiles standalone.

Perf unchanged (reconstruction is one-time per execute). Verified: full suite
+ new cursor cases under GCC -Werror strict, Clang, _GLIBCXX_DEBUG, and libc++
hardening + ASan/UBSan; clang-tidy and clang-format clean.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>

Author: alfC

include/boost/multi/algorithms/fft.NOTES.md

@@ -20,13 +20,27 @@ are computed once and reused across repeated transforms):
 // reusable plan: tables + scratch built once, executed many times
 multi::fft_plan<complex, 2> plan{multi::extents_t<2>{1024, 1024}, multi::fft_forward};
 multi::fft_plan             plan2{A, multi::fft_forward};   // CTAD from a prototype array
-plan(A);            // in place; A: any array/subarray with the planned sizes,
-plan.execute(B);    // of ANY strided layout, repeatedly, no re-allocation
+
+// the plan holds the *extents*; it is applied on a *cursor* -- base + strides,
+// the complement of extents -- obtained from any array/subarray with `.home()`:
+plan(A.home());     // primitive form: cursor supplies base+strides, plan supplies shape
+plan(A);            // sugar: checks A.extensions() == plan shape, then applies on A.home()
+plan.execute(B);    // same, in place, any strided layout, repeatedly, no re-allocation
 
 // one-shot convenience (builds a throw-away plan)
 multi::fft_inplace(A, sign);   // also fft_inplace_forward / fft_inplace_backward
 ```
 
+The plan/cursor split is deliberate: a plan is *shape + precomputed tables*, and
+a cursor (`.home()`) is *where the data lives* (base pointer + strides), with no
+shape of its own. `plan(cursor)` trusts the cursor to have the planned shape (a
+cursor carries no sizes to check); `plan(array)` is the checked convenience that
+extracts `array.home()`. This is also the GPU-facing shape: a cursor is exactly
+what a device kernel receives (base + strides), so the same execution primitive
+carries over (see section 8). Internally `plan(cursor)` rebuilds a strided view
+from the plan's extents + the cursor via `layout_t`'s size = nelems/stride
+invariant; that view then drives the unchanged orchestration.
+
 Properties:
 
 - **No external dependency** (no FFTW/cuFFT). Pure C++17 + standard library.
@@ -414,7 +428,10 @@ layout-agnostic; and plans already separate "build tables once" (host) from
 2. **`fft_memory_space<Ptr>` trait** dispatching host (raw `T*`, current
    paths) vs device (thrust/Multi-CUDA fancy pointers), specialized in an
    opt-in adaptor header so the core stays CUDA-free. The current
-   host-iterator gather fallback must *never* run on device pointers.
+   host-iterator gather fallback must *never* run on device pointers. The
+   entry point is already the right one: `plan(cursor)` accepts any pointer
+   type in `cursor.base()`, so a device cursor flows in unchanged and the
+   trait picks the device orchestration.
 3. **Device table mirror**: keep host `std::vector` tables; add a lazily
    uploaded POD `engine_view` (raw device pointers + sizes) per engine.
 4. **Device stage launchers**: one `__global__` kernel per stage kind around

include/boost/multi/algorithms/fft.hpp

@@ -56,6 +56,8 @@
 #ifndef BOOST_MULTI_ALGORITHMS_FFT_HPP
 #define BOOST_MULTI_ALGORITHMS_FFT_HPP
 
+#include <boost/multi/array_ref.hpp>  // for layout_t and subarray (cursor -> strided view reconstruction)
+
 #include <algorithm>    // for copy, fill, min, max, find_if
 #include <array>        // for plan sizes
 #include <cassert>      // for assert
@@ -1129,6 +1131,56 @@ void fft_apply_last(ViewND&& view, fft_engine<T> const& eng) {  // NOLINT(cppcor
 	}
 }
 
+// A Multi cursor (`.home()`) is base + strides with no extents. These helpers
+// recover its rank and rebuild a full strided view once the extents (which the
+// plan owns) are supplied -- the "extents (plan) + cursor (target)" split.
+template<class Cursor>
+using fft_cursor_strides_t = std::decay_t<decltype(std::declval<Cursor>().strides())>;
+
+// A cursor exposes base()/strides() but, unlike an array/subarray, no sizes().
+template<class C, class = void> struct fft_is_cursor_like : std::false_type {};
+template<class C>
+struct fft_is_cursor_like<C, std::void_t<decltype(std::declval<C>().base()), decltype(std::declval<C>().strides()), std::enable_if_t<!fft_is_multi_like<C>::value>>> : std::true_type {};
+
+template<class Cursor>
+inline constexpr std::ptrdiff_t fft_cursor_rank = static_cast<std::ptrdiff_t>(std::tuple_size_v<fft_cursor_strides_t<Cursor>>);
+
+// Build layout_t<D> from extents and strides (all offsets 0, as for a home
+// cursor): each level carries nelems = size*stride, which is the invariant
+// layout_t uses to recover size() = nelems/stride and extent() bounds.
+template<std::ptrdiff_t D>
+auto fft_layout_from(std::array<std::size_t, static_cast<std::size_t>(D)> const& ext, std::array<std::ptrdiff_t, static_cast<std::size_t>(D)> const& str) -> multi::layout_t<D> {
+	if constexpr(D == 0) {
+		return multi::layout_t<0>{};
+	} else {
+		std::array<std::size_t, static_cast<std::size_t>(D) - 1>    sub_ext{};
+		std::array<std::ptrdiff_t, static_cast<std::size_t>(D) - 1> sub_str{};
+		for(std::size_t i = 0; i != static_cast<std::size_t>(D) - 1; ++i) {
+			sub_ext.at(i) = ext.at(i + 1);
+			sub_str.at(i) = str.at(i + 1);
+		}
+		return multi::layout_t<D>{
+			fft_layout_from<D - 1>(sub_ext, sub_str),
+			str[0], 0,
+			static_cast<std::ptrdiff_t>(ext[0]) * str[0]
+		};  // NOLINT(cppcoreguidelines-pro-bounds-constant-array-index) [0] is valid for D>=1
+	}
+}
+
+template<class Cursor, std::size_t... Is>
+auto fft_strides_array(Cursor const& cur, std::index_sequence<Is...> /*unused*/)
+	-> std::array<std::ptrdiff_t, sizeof...(Is)> {
+	return {{static_cast<std::ptrdiff_t>(cur.template stride<static_cast<multi::dimensionality_type>(Is)>())...}};
+}
+
+// (cursor, extents) -> strided subarray sharing the cursor's memory.
+template<class T, std::ptrdiff_t D, class Cursor>
+auto fft_view_from_cursor(Cursor const& cur, std::array<std::size_t, static_cast<std::size_t>(D)> const& ext) {
+	auto const str = fft_strides_array(cur, std::make_index_sequence<static_cast<std::size_t>(D)>{});
+	using ptr_type = typename Cursor::element_ptr;
+	return multi::subarray<T, D, ptr_type>{fft_layout_from<D>(ext, str), cur.base()};
+}
+
 }  // end namespace detail
 
 // Reusable multidimensional FFT plan: precomputes twiddle tables, stage
@@ -1202,32 +1254,57 @@ class fft_plan {
 	explicit fft_plan(MultiSubArray const& arr, int sign = fft_forward)
 	: sizes_{to_sizes_(arr.sizes(), std::make_index_sequence<static_cast<std::size_t>(D)>{})}, sign_{sign} { init_(); }
 
-	auto sign() const -> int { return sign_; }
-
-	// Execute the plan on `arr` in place. `arr` must have the planned sizes;
-	// its layout (strides, subarray-ness) is free to differ between calls.
-	template<class MultiSubArray>
-	auto operator()(MultiSubArray&& arr) const -> MultiSubArray&& {  // NOLINT(cppcoreguidelines-missing-std-forward)
-		static_assert(std::decay_t<MultiSubArray>::dimensionality == D, "array rank must match the plan");
-		assert(matches_(arr.sizes(), std::make_index_sequence<static_cast<std::size_t>(D)>{}));
+ private:
+	// The axis walk, shared by the cursor and array entry points. `view` is a
+	// strided view of the planned shape; transformed in place.
+	template<class View>
+	void apply_(View&& view) const {  // NOLINT(cppcoreguidelines-missing-std-forward)
 		if constexpr(D == 1) {
-			detail::fft_apply_last(arr(), engine_<0>());
+			detail::fft_apply_last(view, engine_<0>());
 		} else {
 			// Transform the last two axes together, slab by slab (cache
 			// locality), then the remaining axes 0 .. D-3 by static recursion
 			// over rotated views (each axis bound to its engine at compile
 			// time; no runtime axis parameter anywhere).
-			detail::fft_apply_last_pair(arr(), engine_<D - 1>(), engine_<D - 2>());
+			detail::fft_apply_last_pair(view, engine_<D - 1>(), engine_<D - 2>());
 			if constexpr(D >= 3) {
-				transform_middle_<1>(arr().rotated());
+				transform_middle_<1>(view.rotated());
 			}
 		}
+	}
+
+ public:
+	auto sign() const -> int { return sign_; }
+
+	auto sizes() const -> std::array<std::size_t, static_cast<std::size_t>(D)> const& { return sizes_; }
+
+	// Apply the plan on a cursor (`A.home()`) -- base + strides, no extents.
+	// The plan supplies the extents, so the cursor's rank must match and its
+	// shape is *trusted* to be the planned one (a cursor carries no sizes to
+	// check). This is the primitive execution form; the array overload below
+	// delegates to it. Transforms in place; returns the cursor.
+	template<class Cursor, std::enable_if_t<detail::fft_is_cursor_like<std::decay_t<Cursor>>::value, int> = 0>  // NOLINT(modernize-use-constraints) C++17
+	auto operator()(Cursor const& home) const -> Cursor {
+		static_assert(detail::fft_cursor_rank<std::decay_t<Cursor>> == D, "cursor rank must match the plan");
+		auto view = detail::fft_view_from_cursor<T, D>(home, sizes_);
+		apply_(view);
+		return home;  // cursors are value types (base + strides); returned by value
+	}
+
+	// Execute the plan on `arr` in place. `arr` must have the planned sizes;
+	// its layout (strides, subarray-ness) is free to differ between calls.
+	// Delegates to the cursor form via `arr.home()`.
+	template<class MultiSubArray, std::enable_if_t<detail::fft_is_multi_like<std::decay_t<MultiSubArray>>::value, int> = 0>  // NOLINT(modernize-use-constraints) C++17
+	auto operator()(MultiSubArray&& arr) const -> MultiSubArray&& {                                                          // NOLINT(cppcoreguidelines-missing-std-forward)
+		static_assert(std::decay_t<MultiSubArray>::dimensionality == D, "array rank must match the plan");
+		assert(matches_(arr.sizes(), std::make_index_sequence<static_cast<std::size_t>(D)>{}));
+		operator()(arr.home());
 		return std::forward<MultiSubArray>(arr);
 	}
 
-	template<class MultiSubArray>
-	auto execute(MultiSubArray&& arr) const -> MultiSubArray&& {
-		return operator()(std::forward<MultiSubArray>(arr));
+	template<class Target>
+	auto execute(Target&& tgt) const -> decltype(auto) {
+		return operator()(std::forward<Target>(tgt));
 	}
 };
 
@@ -1241,20 +1318,20 @@ fft_plan(MultiSubArray const&) -> fft_plan<typename MultiSubArray::element, Mult
 // throw-away plan; for repeated transforms of the same shape, build an
 // fft_plan once and execute it).
 template<class MultiSubArray>
-auto fft_inplace(MultiSubArray&& A, int sign = fft_forward) -> MultiSubArray&& {
+auto fft_inplace(MultiSubArray&& arr, int sign = fft_forward) -> MultiSubArray&& {
 	using array_type = std::decay_t<MultiSubArray>;
-	fft_plan<typename array_type::element, array_type::dimensionality> const plan{A, sign};
-	return plan(std::forward<MultiSubArray>(A));
+	fft_plan<typename array_type::element, array_type::dimensionality> const plan{arr, sign};
+	return plan(std::forward<MultiSubArray>(arr));
 }
 
 template<class MultiSubArray>
-auto fft_inplace_forward(MultiSubArray&& A) -> MultiSubArray&& {
-	return fft_inplace(std::forward<MultiSubArray>(A), fft_forward);
+auto fft_inplace_forward(MultiSubArray&& arr) -> MultiSubArray&& {
+	return fft_inplace(std::forward<MultiSubArray>(arr), fft_forward);
 }
 
 template<class MultiSubArray>
-auto fft_inplace_backward(MultiSubArray&& A) -> MultiSubArray&& {
-	return fft_inplace(std::forward<MultiSubArray>(A), fft_backward);
+auto fft_inplace_backward(MultiSubArray&& arr) -> MultiSubArray&& {
+	return fft_inplace(std::forward<MultiSubArray>(arr), fft_backward);
 }
 
 }  // end namespace boost::multi

test/algorithms_fft.cpp

@@ -327,6 +327,48 @@ auto main() -> int {  // NOLINT(readability-function-cognitive-complexity,bugpro
 		BOOST_TEST( m / static_cast<double>(n) < tol );
 	}
 
+	// plan applied on a cursor (.home()): extents live in the plan, the cursor
+	// supplies base + strides
+	{
+		multi::array<complex, 2> arr({6, 10}, complex{});
+		for(int i = 0; i != 6; ++i) {
+			for(int j = 0; j != 10; ++j) {
+				arr[i][j] = complex{static_cast<double>(i - j), static_cast<double>((i * j) % 7)};
+			}
+		}
+		auto reference = arr;
+
+		multi::fft_plan<complex, 2> const plan{
+			multi::extents_t<2>{6, 10},
+			multi::fft_forward
+		};
+		plan(arr.home());  // apply on the cursor directly
+		plan(reference);   // apply on the array (delegates to .home())
+		BOOST_TEST( max_abs_diff(arr.elements(), reference.elements()) < tol );
+	}
+
+	// cursor application on a strided sub-block matches the same plan on a contiguous copy
+	{
+		multi::array<complex, 3> big({6, 7, 8}, complex{});
+		int                      c = 0;
+		for(auto& e : big.elements()) {
+			e = complex{static_cast<double>(c % 11) - 5.0, static_cast<double>(c % 7) - 3.0};
+			++c;
+		}
+		auto&& blk = big({1, 5}, {2, 6}, {1, 7});  // 4 x 4 x 6 strided view
+
+		multi::array<complex, 3> flat{blk};  // contiguous copy of the same values
+
+		multi::fft_plan<complex, 3> const plan{
+			multi::extents_t<3>{4, 4, 6},
+			multi::fft_forward
+		};
+		plan(flat.home());  // contiguous cursor
+		plan(blk.home());   // strided cursor, same plan
+
+		BOOST_TEST( max_abs_diff(blk.elements(), flat.elements()) < tol );
+	}
+
 	// 3D round-trip
 	{
 		multi::array<complex, 3> arr({3, 4, 5}, complex{});