Fix f ordering

CHANGELOG.md

@@ -7,6 +7,10 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 ## [Unreleased]
 
+### Fixed
+- Preserve Fortran order when copying CuPy arrays to the host through
+  `to_numpy`, `host_call`, `evaluate_on_host`, and `PyccelKernel` conversions.
+
 ### Added
 - `kernel_testing.emulated_launches()` runs every `CudaKernel` launch in a block
   on the CPU, on the host buffers of the fake CuPy arrays, so code that launches

docs/source/api.md

@@ -206,6 +206,10 @@ Converts to a host-side NumPy array. CuPy arrays are copied from device to
 host. Other array-like inputs are passed through `numpy.asarray`; NumPy arrays
 may therefore be returned as-is rather than copied.
 
+F-contiguous CuPy inputs keep F order on the host; other CuPy inputs become
+C-contiguous. Arbitrary device strides are not preserved. See
+[Array ordering and strides](guides/array-ordering.md).
+
 ### `to_cupy(array)`
 
 Converts an array-like input to a CuPy array. Raises `ImportError` if CuPy or

docs/source/guides/array-ordering.md

@@ -0,0 +1,139 @@
+# Array ordering and strides
+
+CuNumpy supports C-ordered (row-major) and F-ordered (column-major) arrays on
+both NumPy and CuPy. Choose the layout when allocating an array, and keep it
+consistent with the compiled kernels that consume it.
+
+For a two-dimensional array `a[row, column]`, C order stores each row together;
+F order stores each column together. Indexing and shape stay the same. F order
+can help code that processes columns of a particle buffer, but measure the
+whole workload, including kernels and copies, before choosing a layout.
+
+## Allocate and inspect
+
+```python
+import cunumpy as xp
+
+markers = xp.zeros((100, 6), dtype=xp.float64, order="F")
+assert markers.flags.f_contiguous
+
+positions = markers[:, 1:4]
+assert positions.flags.f_contiguous
+assert positions.shape == (100, 3)
+```
+
+A block of complete columns of an F-contiguous array remains F-contiguous.
+A block of complete rows of a C-contiguous array remains C-contiguous.
+Stepped slices such as `markers[::2, :]` generally have neither layout.
+Check `a.flags.c_contiguous`, `a.flags.f_contiguous`, and `a.strides` when
+debugging a kernel boundary. Strides are measured in bytes. Arrays with
+singleton dimensions or no elements can satisfy both contiguity flags.
+
+Use `xp.asfortranarray(a)` to require F order or `xp.ascontiguousarray(a)` to
+require C order. These may allocate a copy; a transpose changes the shape and
+indexing meaning and is not a substitute for an order conversion.
+
+## Host and device conversions
+
+| Operation | Layout behavior |
+| --- | --- |
+| `xp.to_numpy(cupy_array)` | F-contiguous input becomes F-contiguous host storage; other input becomes C-contiguous host storage. |
+| `xp.to_numpy(numpy_array)` | Uses `numpy.asarray`; an ordinary NumPy array keeps its layout and strides without copying. |
+| `xp.to_cupy(numpy_array)` | Uses `cupy.asarray`, preserving C or F order for contiguous input. |
+| `xp.to_cunumpy(a)` | Uses the conversion for the active backend. |
+| `xp.host_call`, `xp.evaluate_on_host` | Device arguments follow `to_numpy`; returned NumPy arrays follow `to_cupy`. |
+| `xp.kernels.PyccelKernel` | Device arguments use the same host conversion as `to_numpy`; declared outputs are copied back into the original device arrays. |
+
+The device-to-host conversion uses `array.get(order="A")`. This preserves
+F-contiguous layout, including complete column blocks, but does not preserve
+arbitrary strides or shared storage between distinct views. Calling CuPy's
+`array.get()` directly defaults to C order; see the
+[CuPy array documentation](https://docs.cupy.dev/en/stable/reference/generated/cupy.ndarray.html).
+
+The host function itself controls the layout of any new arrays it creates.
+Returning a C-ordered result from `host_call` does not make that result F-ordered
+just because an input was F-ordered.
+
+```{note}
+Preserving F order in CuPy-to-host conversions is an unreleased fix. Earlier
+releases used `array.get()` and produced C-ordered host copies. Applications
+depending on F order at a host-kernel boundary need a release containing this fix.
+```
+
+## Copies, masks, and assignment
+
+`a.copy()` defaults to C order on both NumPy and CuPy. Use `a.copy(order="K")`
+to retain the layout of a C- or F-contiguous array, or `a.copy(order="F")` to
+require F order. `order="K"` is not a promise to preserve arbitrary strides.
+
+Boolean row masks and fancy indexing create new arrays. Do not rely on those
+results retaining F order. For example, a boolean row selection from an
+F-ordered marker buffer produces a C-ordered result on CuPy. Normalize an
+independent result before passing it to a kernel requiring F order:
+
+```python
+alive = markers[:, 0] >= 0
+selected = xp.asfortranarray(markers[alive])
+assert selected.flags.f_contiguous
+```
+
+Assignment into an existing buffer keeps the destination's layout:
+
+```python
+markers[...] = markers.copy()  # C-ordered temporary, F-ordered destination
+assert markers.flags.f_contiguous
+```
+
+Likewise, `xp.algorithms.compact_by_mask(alive, markers)` gathers selected rows
+and writes them into the front of the original buffer, preserving that
+buffer's layout. Only the returned count of leading rows is valid. However,
+`markers[:n_kept]` generally loses F contiguity when it trims rows: the columns
+still have the original buffer's stride. Pass the full buffer plus a valid-row
+count to a kernel designed for that interface, or make an F-contiguous copy
+of the prefix for a kernel that requires one. A stride-aware CUDA array view
+can consume the prefix directly.
+
+Check or normalize arrays created by concatenation, file readers, and other
+libraries before passing them across a layout-sensitive boundary.
+
+## Compiled kernel boundaries
+
+Pyccel array signatures specify an expected order. A two-dimensional
+`"float[:, :]"` argument expects C order; use `"float[:, :](order=F)"` for an
+F-ordered argument. Every caller must supply the declared layout, including
+temporary arrays and array attributes in argument objects.
+`PyccelKernel` preserves F-contiguous device arguments when moving them to
+the host; it does not inspect signatures or fix mismatches for you. The same
+requirement applies to a `Kernel` using its host fallback without a CUDA version.
+See [Host kernels with GPU data](../kernels/pyccel-kernel.md).
+
+CUDA `Array2D<T>` and other `ArrayND<T>` view arguments carry shape and strides,
+so indexing through the view supports C order, F order, and strided slices.
+The same applies to array view fields in `CudaStruct` and to emulated launches.
+Raw pointer arguments and `CArrayND<T>` views require C-contiguous storage.
+Hand-written indexing such as `data[row * n_columns + column]` assumes C order;
+use the array view's indexing to handle both layouts. See
+[Index macros and array views](../kernels/cuda-kernel.md).
+
+Some convenience helpers deliberately require C order:
+
+| Helper | Contract |
+| --- | --- |
+| `xp.as_device_array` | Returns a C-contiguous device array, copying if needed. |
+| `xp.kernels.as_kernel_array` | Returns a C-contiguous array on the side of `like`. |
+| `xp.kernels.kernel_output` | Yields a C-contiguous working buffer and copies updates back into the original output when needed. |
+
+For an F-order kernel, allocate or normalize its arguments explicitly instead
+of using these helpers to prepare F-ordered inputs.
+
+## Check values and layout separately
+
+Numerical equality does not check memory order. Tests for layout-sensitive
+code should assert the host argument's contiguity inside the called function,
+as well as checking values. Use a two-dimensional shape with both dimensions
+greater than one so C and F order are distinguishable.
+
+The fake CuPy supports these conversion checks without a GPU. Its
+`kernel_testing.host_buffer(a)` exposes the existing NumPy storage without a
+copy, preserving layout and strides. It does not exercise the `to_numpy`
+transfer path; test that path separately. See [Testing kernels](../kernels/testing.md).

docs/source/guides/data-movement.md

@@ -21,6 +21,10 @@ implicitly: every transfer is a visible function call.
   `to_numpy()`, on CuPy like `to_cupy()`.
 * None of them modify the source array or change the active backend.
 
+CuPy-to-host conversions preserve F-contiguous layout; other device inputs
+become C-contiguous host arrays. See [Array ordering and strides](array-ordering.md)
+for conversion rules, copies, and compiled kernel requirements.
+
 ## Transfer at boundaries, not in loops
 
 Load or generate data, move it to the device once, run the whole computation

docs/source/index.md

@@ -47,6 +47,7 @@ quickstart
 guides/backends
 guides/portable-code
 guides/data-movement
+guides/array-ordering
 guides/gpu-devices
 guides/execution-helpers
 guides/mpi

docs/source/kernels/pyccel-kernel.md

@@ -37,6 +37,11 @@ with xp.use_backend("cupy"):
     smooth_kernel(field, out)  # out is updated on the device
 ```
 
+Device-to-host conversion preserves F-contiguous arrays. The host kernel must
+accept that order: the wrapper does not adapt arrays to a compiled signature.
+See [Array ordering and strides](../guides/array-ordering.md) for Pyccel signatures
+and the operations that can change layout.
+
 ## Declare the outputs
 
 The wrapper cannot know which arguments the function writes. By default it

pyproject.toml

@@ -5,7 +5,7 @@ requires = [ "setuptools", "wheel" ]
 
 [project]
 name = "cunumpy"
-version = "0.6.1"
+version = "0.6.2"
 description = "Simple wrapper for numpy and cupy. Replace `import numpy as np` with `import cunumpy as xp`."
 readme = "README.md"
 keywords = [ "python" ]

src/cunumpy/xp.py

@@ -269,7 +269,7 @@ def synchronize() -> None:
 def _to_numpy(array: Any) -> np.ndarray:
     """`to_numpy` without transfer counting, for internal use."""
     if get_array_backend(array) == "cupy":
-        return array.get()
+        return array.get(order="A")
 
     return np.asarray(array)
 
@@ -289,6 +289,7 @@ def to_numpy(array: Any) -> np.ndarray:
 
     A CuPy array is copied to the host, which `count_transfers()` counts as a
     ``to_host`` transfer; anything else is passed through `numpy.asarray`.
+    Fortran-contiguous CuPy arrays keep F order; other CuPy arrays use C order.
     """
     result = _to_numpy(array)
     if _COUNTERS and get_array_backend(array) == "cupy":

tests/unit/test_cupy.py

@@ -2,6 +2,50 @@
 import pytest
 
 import cunumpy as xp
+from cunumpy.kernels import PyccelKernel
+
+
+@pytest.mark.parametrize("order", ["C", "F"])
+@pytest.mark.parametrize("column_block", [False, True])
+@pytest.mark.parametrize(
+    "conversion", ["to_numpy", "host_call", "evaluate_on_host", "pyccel_kernel"]
+)
+def test_host_conversion_preserves_order(order, column_block, conversion):
+    if not xp.cupy_available():
+        pytest.skip("CuPy not installed or not functional")
+
+    expected = np.array(np.arange(30.0).reshape(5, 6), order=order)
+    with xp.use_backend("cupy"):
+        device = xp.to_cupy(expected)
+        if column_block:
+            device = device[:, 1:4]
+            expected = expected[:, 1:4]
+
+        def check_host(array):
+            assert isinstance(array, np.ndarray)
+            np.testing.assert_array_equal(array, expected)
+            assert array.flags.f_contiguous == (order == "F")
+            assert array.flags.c_contiguous == (order == "C")
+            return array
+
+        if conversion == "to_numpy":
+            check_host(xp.to_numpy(device))
+        else:
+            if conversion == "host_call":
+                result = xp.host_call(check_host, device)
+            elif conversion == "evaluate_on_host":
+
+                class Model:
+                    @xp.evaluate_on_host
+                    def evaluate(self, array):
+                        return check_host(array)
+
+                result = Model().evaluate(device)
+            else:
+                result = PyccelKernel(check_host)(device)
+
+            assert xp.is_gpu(result)
+            check_host(xp.to_numpy(result))
 
 
 def test_to_cupy_available():

tests/unit/test_execution_transfers.py

@@ -31,10 +31,10 @@ def set(self, host):
             cp.copies.append("set")
             np.copyto(self.array, host)
 
-        def get(self, out=None, stream=None):
+        def get(self, out=None, stream=None, order="C"):
             cp.copies.append(("get", stream))
             if out is None:
-                return self.array.copy()
+                return self.array.copy(order=order)
             np.copyto(out, self.array)
             return out
 

tests/unit/test_transfers.py

@@ -42,8 +42,8 @@ def shape(self):
     def dtype(self):
         return self.data.dtype
 
-    def get(self):
-        return self.data.copy()
+    def get(self, order="C"):
+        return self.data.copy(order=order)
 
     def __setitem__(self, key, value):
         self.data[key] = value.data if isinstance(value, _FakeDeviceArray) else value