Repository: cupy/cupy
Pull request: Ability to pass structured data types by value as kernel parameters
This PR is directly related to issue #4828.

Currently, we can only pass by value the following numpy items as kernel parameters: np.int8, np.int16, np.int32, np.int64, np.uint8, np.uint16, np.uint32, np.uint64, np.float16, np.float32, np.float64, np.complex64 and np.comples128 and other numpy builtin datatypes if their itemsize is either 1, 2, 4, 8 or 16 bytes.

This patch aims to enable **pass by value** for any structured data type built with numpy.dtype as kernel parameter.
It enhances current kernel launch capabilities of cupy without changing its interface. It is the responsibility of the user to define a numpy dtype that matches exactly device structure layout.

It provides the ability to:
- pass cuda builtin vector types (char2, double3, int4) as kernel parameters.
- pass "static arrays" of parameters such has `struct double32 { double[32] val };` and other packed structures.
- define and use any cuda composite type as kernel parameter (including structure of structures) such as
 ```cuda
   struct custom_user_struct {
    int4 a;
    char b;
    double c[2];
    short1 d;
    unsigned long long int e[3];
  };
```
Implementation status:
- [x] Patch `cupy/cuda/function.pyx` to enable pass by value to any dtype
- [x] Add one test for each use case in ` tests/cupy_tests/core_tests/test_function.py`: builtin vector type, static array and custom struct
- [x] Document this under `User-Defined Kernels` section
- [x] Throw an error when numpy array argument is not of size 1 (this would be error-prone)
- [x] Add some examples of complex structs 

close #3359, close #4828

Problem statement:
Is there a way to pass a structured buffer to a custom kernel, like arrow cuda buffer or flatbuffers buffer
This isn't actually an issue, just a question and maybe a feature request or at least documentation request, as I didn't find other ways to ask a question, like chat or group. 

I'm looking for a way / example to pass a structured byte buffer to a custom kernel. 

I didn't find any documentation or examples of how to use low level cuda api, for example 

cupy.cuda.MemoryPointer with a custom kernel, ie pass it as an argument, create it from the host and/or device buffer?

I'm looking for an integration with an Arrow and/or Flatbuffers or alike offset based flat structures.


Ability to pass structured types by value as kernel parameters
Currently, it is not possible to pass structures **by value** as kernel arguments. Only some fixed size numpy dtype arguments have been hardcoded yet in `cupy/cuda/function.pyx`. This implies that a cupy kernel launch only succeeds with numpy.generic arguments which are scalar types with itemsize 1, 2, 4, 8 or 16 bytes. This allows us to pass directly by value np.int8 up to np.complex128 numpy data types as kernel arguments. I would like to extend pass-by-value to any numpy dtype.

I am using RawKernels extensively and I would like to be able to pass cuda builtin vector types such as `float3` as kernel parameters. This would make multidimensional code in physics much clearer for example (ie. replace many scalar kernel arguments by vectorized arguments or static arrays of scalars). Depending on device compute capability we can now pass up to 4kB of parameters to cuda kernels, and doing it in 16 bytes chunks is not very efficient.

The definition of the structure layout using numpy.dtype should be the responsibility of the user. This off course implies defining the struct with the same offsets, alignment and itemsize between host and device. `pyopencl` and `pycuda` offer utility functions to convert a device side struct to a numpy dtype with matching layout, but this feature would be overkill for now.

Please note that **this feature is not related to general structured array support in cupy**. It is however already possible to pass a pointer to any array of structures to a RawKernel: simply cast any numpy array to np.int8 data type and let cupy handle np.int8 cuda arrays. Within the kernel, the pointer will be cast to the kernel argument type anyway. This off course only works if host and device structures layout match exactly.

Here is a simple example to illustrate what I would be able to do:
```python
import numpy as np
import cupy as cp

code = r"""
    template <typename T>
    __global__ void affect(T* const __restrict__ output,
                         const T arg0, const T arg1,
                         const size_t nelements) {
        ptrdiff_t i = blockIdx.x * blockDim.x + threadIdx.x;
        if (i >= nelements)
            return;
        output[i] = (i%2 == 0 ? arg0 : arg1);
    }
"""
module = cp.RawModule(code=code,
                      name_expressions=('affect<float3>',),
                      options=('-std=c++11',))
kernel = module.get_function('affect<float3>')

float3 = np.dtype( { 'names': ['x', 'y', 'z'],
                     'formats': [np.float32]*3 } )

arg0 = np.random.rand(3).astype(np.float32).view(float3)
arg1 = np.random.rand(3).astype(np.float32).view(float3)

N = 512
h_output = np.empty(dtype=float3, shape=N)
d_output = cp.asarray(h_output.view(np.float32))

grid = (N,)
block = (256,)
args = (d_output, arg0, arg1, np.uint64(N))
kernel(grid, block, args)

h_output[0::2] = arg0
h_output[1::2] = arg1

np.testing.assert_array_equal(h_output,
                              cp.asnumpy(d_output).view(float3))
```

Hints from issue discussion:
Thanks for the feedback, `cupy.cuda.MemoryPointer` can be passed as an argument to raw kernels but not documented very well. We now have a gitter channel so feel free to join and ask questions when needed. https://gitter.im/cupy/community
This example is now running with PR #4829.
This is very nice capability indeed! I've skimmed over your PR #4829 btw (thanks for quickly putting together a nice work there), but this issue seems to be the more appropriate place for high-level discussions.

> Depending on device compute capability we can now pass up to 4kB of parameters to cuda kernels, and doing it in 16 bytes chunks is not very efficient.

Just to supplement the 4kB thing: This refers to a dedicated constant memory of size 4kB for storing kernel parameters:
https://docs.nvidia.com/cuda/cuda-c-programming-guide/index.html#function-parameters

On the 16-byte limitation: So it indeed looks like this is a hard-coded cap, but actually we've opened a shortcut to bypass it. If you pass an object of type `cupy.cuda.MemoryPointer` (which is the memory object backing any CuPy array, and can be accessed via `arr.data`) to your kernel, then essentially it is a `void*` pointer to a contiguous chunk of GPU memory, meaning it can be a struct of unlimited size. This is actually used in CuPy's internal, for example see this function:
https://github.com/cupy/cupy/blob/2fd13b80ba108bcc12d102fd7fa9d5da0ced7eaf/cupy/core/_cub_reduction.pyx#L412
There I assumed the maximal struct size is of 32 bytes (for a reason not relevant here), which exceeds the 16 byte cap. 

But, the drawback of passing `MemoryPointer` is that your kernel signature has to be modified as
```c++
    template <typename T>
    __global__ void affect(T* const __restrict__ output,
                         const void* arg0, const void* arg1,
                         const size_t nelements)
``` 
and you cast `arg0` and `arg1` to the template type `T` in the kernel body, which may not be pleasant to use (ex: inconvenient for hooking up existing codes), although in terms of performance I think the compiler output should be nearly identical. In addition, this does not support passing values on host, only on device.

> `pyopencl` and `pycuda` offer utility functions to convert a device side struct to a numpy dtype with matching layout, but this feature would be overkill for now.

I think this would be a very useful utility to us (and one that I must have missed!), could you give us a few pointers to that? 🙂

> Please note that **this feature is not related to general structured array support in cupy**.

I'd say it's related, and actually your PR is a good step toward it (#2031). Although there would still be plenty of work needed in CuPy's internal, it's a good basis for thinking about the design.
I am not sure that your example with `void*` and `MemoryPointer` is 100% equivalent. If I understand this correctly you first need to allocate GPU memory with `cupy.cuda.memory.alloc`, wrap it in a `cupy.cuda.MemoryPointer`, find some way to initialize this memory on device side with a RawKernel or by raw host->device copy from some initialized host buffer with matching memory layout, change your kernel arguments to pointer to required type and then dereference those pointers in your actual code:
```cuda
    template <typename T>
    __global__ void affect(T* const __restrict__ output,
                         const T* arg0, const T* arg1,  // <-- void* can be replaced by T* here
                         const size_t nelements)
```
Assuming  compute capability 2.x and higher, kernel arguments (arg0 and arg1 pointers) would be in a dedicated memory bank (actually distinct from device `__constant__` memory space), but underlying data would lie in `__global__` memory. I was not aware of `MemoryPointer`, this indeed makes it possible in the actual version of the library and opens the way for general structured array support in cupy. In order for this to work cupy memory allocators will have to support custom alignement contraints for custom types array allocations.

Just so that it is clear to anyone, the equivalent way to pass `arg0` and `arg1` by value (in `__constant__` memory instead of dedicated kernel argument memory) without PR #4829 would be to use such a kernel:
 ```cuda
    __constant__ T arg0;
    __constant__ T arg1;
    __global__ void affect(T* const __restrict__ output,
                         const size_t nelements)
```
where T cannot be a template anymore. I am however not aware of hostside `__constant__` memory management in cupy (`cudaGetSymbolAddress`,  `cudaGetSymbolSize`, `cudaMemcpyToSymbol`, `cudaMemcpyFromSymbol`).

>  I think this would be a very useful utility to us (and one that I must have missed!), could you give us a few pointers to that

Indeed, the most difficult part of both methods is to create the host side buffer with matching memory layout.

In PyOpenCL, there is a dictionary that tracks device side typenames and corresponding numpy dtype. At module initialization,
this dictionary is populated with scalar and vector builtins. A user can register a new type with `get_or_register_dtype(typename, dtype)`. The most interesting part happens in [match_dtype_to_c_struct ](https://documen.tician.de/pyopencl/array.html?highlight=struct#pyopencl.tools.match_dtype_to_c_struct) where you can provide a typename and a numpy dtype as input. It then proceeds to generate the corresponding C structure declaration (OpenCL code) and calls a kernel to determine device side struct memory layout (see  `cupy/tests/cupy_tests/core_tests/test_function.py::test_custom_user_struct`  in my PR). As an output you get the generated struct code as well as a device compatible numpy dtype.

I think that such a type dataset should be device dependent, even for CUDA. I am not sure about this, but I imagine that devices with different compute capabilities may impose different structure memory layouts. If you point me where I could implement this in CuPy, I may have the time to implement such a feature later this month.

Commit message: Add a test for custom struct examples.

Patch to review:
diff --git a/cupy/cuda/function.pyx b/cupy/cuda/function.pyx
index 8a8d38304a3..fa616be2445 100644
--- a/cupy/cuda/function.pyx
+++ b/cupy/cuda/function.pyx
@@ -90,6 +90,15 @@ cdef class CIntptr(CPointer):
         self.ptr = <void*>&self.val
 
 
+cdef class CNumpyArray(CPointer):
+    cdef:
+        object val
+
+    def __init__(self, v):
+        self.val = v
+        self.ptr = <void*><size_t>v.__array_interface__['data'][0]
+
+
 cdef set _pointer_numpy_types = {numpy.dtype(i).type
                                  for i in '?bhilqBHILQefdFD'}
 
@@ -109,6 +118,20 @@ cdef inline CPointer _pointer(x):
         return x
     if isinstance(x, (TextureObject, SurfaceObject)):
         return CUIntMax(x.ptr)
+    if isinstance(x, numpy.ndarray):
+        # All numpy.ndarray work with CNumpyArray to pass a kernel argument by
+        # value. Here we allow only arrays of size one so that users do not
+        # mistakenly send numpy.ndarrays instead of cupy.ndarrays to kernels.
+        # This may happen if they forget to convert numpy arrays to cupy arrays
+        # prior to kernel call and would pass silently without this check.
+        if (x.size == 1):
+            return CNumpyArray(x)
+        else:
+            msg = ('You are trying to pass a numpy.ndarray of shape {} as a '
+                   'kernel parameter. Only numpy.ndarrays of size one can be '
+                   'passed by value. If you meant to pass a pointer to __glob'
+                   'al__ memory, you need to pass a cupy.ndarray instead.')
+            raise TypeError(msg.format(x.shape))
 
     if type(x) not in _pointer_numpy_types:
         if isinstance(x, int):
diff --git a/docs/source/tutorial/kernel.rst b/docs/source/tutorial/kernel.rst
index e6f9f8388c6..abf3c4cb685 100644
--- a/docs/source/tutorial/kernel.rst
+++ b/docs/source/tutorial/kernel.rst
@@ -263,16 +263,122 @@ Accessing texture (surface) memory in :class:`~cupy.RawKernel` is supported via
     The kernel does not have return values.
     You need to pass both input arrays and output arrays as arguments.
 
-.. note::
-    No validation will be performed by CuPy for arguments passed to the kernel, including types and number of arguments.
-    Especially note that when passing :class:`~cupy.ndarray`, its ``dtype`` should match with the type of the argument declared in the method signature of the CUDA source code (unless you are casting arrays intentionally).
-    For example, ``cupy.float32`` and ``cupy.uint64`` arrays must be passed to the argument typed as ``float*`` and ``unsigned long long*``.
-    For Python primitive types, ``int``, ``float``, ``complex`` and ``bool`` map to ``long long``, ``double``, ``cuDoubleComplex`` and ``bool``, respectively.
-
 .. note::
     When using ``printf()`` in your CUDA kernel, you may need to synchronize the stream to see the output.
     You can use ``cupy.cuda.Stream.null.synchronize()`` if you are using the default stream.
 
+Kernel arguments
+----------------
+Python primitive types and NumPy scalars are passed to the kernel by value.
+Array arguments (pointer arguments) have to be passed as CuPy ndarrays.
+No validation is performed by CuPy for arguments passed to the kernel, including types and number of arguments.
+
+Especially note that when passing a CuPy :class:`~cupy.ndarray`, its ``dtype`` should match with the type of the argument declared in the function signature of the CUDA source code (unless you are casting arrays intentionally). 
+
+As an example, ``cupy.float32`` and ``cupy.uint64`` arrays must be passed to the argument typed as ``float*`` and ``unsigned long long*``, respectively. CuPy does not directly support arrays of non-primitive types such as ``float3``, but nothing prevents you from casting a ``float*`` or ``void*`` to a ``float3*`` in a kernel.
+
+Python primitive types, ``int``, ``float``, ``complex`` and ``bool`` map to ``long long``, ``double``, ``cuDoubleComplex`` and ``bool``, respectively.
+
+NumPy scalars (``numpy.generic``) and NumPy arrays (``numpy.ndarray``) **of size one** 
+are passed to the kernel by value.
+This means that you can pass by value any base NumPy types such as ``numpy.int8`` or ``numpy.float64``, provided the kernel arguments match in size. You can refer to this table to match CuPy/NumPy dtype and CUDA types:
+
++-----------------+-----------------------------------------------+------------------+
+| CuPy/NumPy type | Corresponding kernel types                    | itemsize (bytes) |
++=================+===============================================+==================+
+| bool            | bool                                          | 1                |
++-----------------+-----------------------------------------------+------------------+
+| int8            | char, signed char                             | 1                |
++-----------------+-----------------------------------------------+------------------+
+| int16           | short, signed short                           | 2                |
++-----------------+-----------------------------------------------+------------------+
+| int32           | int, signed int                               | 4                |
++-----------------+-----------------------------------------------+------------------+
+| int64           | long long, signed long long                   | 8                |
++-----------------+-----------------------------------------------+------------------+
+| uint8           | unsigned char                                 | 1                |
++-----------------+-----------------------------------------------+------------------+
+| uint16          | unsigned short                                | 2                |
++-----------------+-----------------------------------------------+------------------+
+| uint32          | unsigned int                                  | 4                |
++-----------------+-----------------------------------------------+------------------+
+| uint64          | unsigned long long                            | 8                |
++-----------------+-----------------------------------------------+------------------+
+| float16         | half                                          | 2                |
++-----------------+-----------------------------------------------+------------------+
+| float32         | float                                         | 4                |
++-----------------+-----------------------------------------------+------------------+
+| float64         | double                                        | 8                |
++-----------------+-----------------------------------------------+------------------+
+| complex64       | float2, cuFloatComplex, complex<float>        | 8                |
++-----------------+-----------------------------------------------+------------------+
+| complex128      | double2, cuDoubleComplex, complex<double>     | 16               |
++-----------------+-----------------------------------------------+------------------+
+
+The CUDA standard guarantees that the size of fundamental types on the host and device always match.
+The itemsize of ``size_t``, ``ptrdiff_t``, ``intptr_t``, ``uintptr_t``, 
+``long``, ``signed long`` and ``unsigned long`` are however platform dependent. 
+To pass any CUDA vector builtins such as ``float3`` or any other user defined structure 
+as kernel arguments (provided it matches the device-side kernel parameter type), see section :ref:`custom_user_structs` below.
+
+.. _custom_user_structs:
+
+Custom user types
+-----------------
+
+Is is possible to use custom types (composite types such as structures and structures of structures) 
+as kernel arguments by defining a custom NumPy dtype.
+When doing this, it is your responsability to match host and device structure memory layout.
+The CUDA standard guarantees that the size of fundamental types on the host and device always match.
+It may however impose device alignment requirements on composite types.
+This means that for composite types, the struct member offsets may be different from what you might expect.
+
+When a kernel argument is passed by value, the CUDA driver will copy exactly ``sizeof(param_type)`` bytes starting from the beginning of the NumPy object data pointer, where ``param_type`` is the parameter type in your kernel. 
+You have to match ``param_type``'s memory layout (ex: size, alignment and struct padding/packing) 
+by defining a corresponding `NumPy dtype <https://numpy.org/doc/stable/reference/arrays.dtypes.html>`_.
+
+For builtin CUDA vector types such as ``int2`` and ``double4`` and other packed structures with 
+named members you can directly define such NumPy dtypes as the following:
+
+.. doctest::
+
+    >>> import numpy as np
+    >>> names = ['x', 'y', 'z']
+    >>> types = [np.float32]*3
+    >>> float3 = np.dtype({'names': names, 'formats': types})
+    >>> arg = np.random.rand(3).astype(np.float32).view(float3)
+    >>> print(arg)  # doctest: +SKIP
+    [(0.9940819, 0.62873816, 0.8953669)]
+    >>> arg['x'] = 42.0
+    >>> print(arg)  # doctest: +SKIP
+    [(42., 0.62873816, 0.8953669)]
+
+Here ``arg`` can be used directly as a kernel argument.
+When there is no need to name fields you may prefer this syntax to define packed structures such as 
+vectors or matrices:
+
+.. doctest::
+
+    >>> import numpy as np
+    >>> float5x5 = np.dtype({'names': ['dummy'], 'formats': [(np.float32,(5,5))]}) 
+    >>> arg = np.random.rand(25).astype(np.float32).view(float5x5)
+    >>> print(arg.itemsize)
+    100
+
+Here ``arg`` represents a 100-byte scalar (i.e. a NumPy array of size 1)
+that can be passed by value to any kernel.
+Kernel parameters are passed by value in a dedicated 4kB memory bank which has its own cache with broadcast.
+Upper bound for total kernel parameters size is thus 4kB
+(see this `link <https://docs.nvidia.com/cuda/cuda-c-programming-guide/index.html#function-parameters>`_). 
+It may be important to note that this dedicated memory bank is not shared with the device ``__constant__`` memory space.
+
+For now, CuPy offers no helper routines to create user defined composite types. 
+Such composite types can however be built recursively using NumPy dtype `offsets` and `itemsize` capabilities,
+see ``cupy/examples/user_structs`` for examples of advanced usage. 
+
+.. warning::
+    You cannot directly pass static arrays as kernel arguments with the ``type arg[N]`` syntax where N is a compile time constant. The signature of ``__global__ void kernel(float arg[5])`` is seen as ``__global__ void kernel(float* arg)`` by the compiler. If you want to pass five floats to the kernel by value you need to define a custom structure ``struct float5 { float val[5]; };`` and modify the kernel signature to ``__global__ void kernel(float5 arg)``.
+
 
 Raw modules
 -----------
diff --git a/examples/custom_struct/README.md b/examples/custom_struct/README.md
new file mode 100644
index 00000000000..de814502b85
--- /dev/null
+++ b/examples/custom_struct/README.md
@@ -0,0 +1,11 @@
+# Custom user structure examples
+
+This folder contains examples of custom user structures in `cupy.RawKernel` (see [https://docs.cupy.dev/en/stable/tutorial/kernel.html](https://docs.cupy.dev/en/stable/tutorial/kernel.html) for corresponding documentation).
+
+This folder provides three scripts ranked by increasing complexity:
+
+1. `builtins_vectors.py` shows how to use CUDA builtin vectors such as `float4` both as scalar parameter (pass by value from host) and array parameter in RawKernels.
+2. `packed_matrix.py` demonstrates how to create and use templated packed structures in RawModules.
+3. `complex_struct.py` illustrates the possibility to recursively build complex NumPy dtypes matching device structure memory layout.
+
+All examples can be run as simple python scripts: `python3.x example_name.py`.
diff --git a/examples/custom_struct/builtin_vectors.py b/examples/custom_struct/builtin_vectors.py
new file mode 100644
index 00000000000..667d6e71307
--- /dev/null
+++ b/examples/custom_struct/builtin_vectors.py
@@ -0,0 +1,48 @@
+import sys
+import numpy
+import cupy
+
+code = '''
+__device__ double3 operator+(const double3& lhs, const double3& rhs) {
+    return make_double3(lhs.x + rhs.x,
+                        lhs.y + rhs.y,
+                        lhs.z + rhs.z);
+}
+
+extern "C" __global__ void sum_kernel(const double3* lhs,
+                                            double3  rhs,
+                                            double3* out) {
+  int i = threadIdx.x;
+  out[i] = lhs[i] + rhs;
+}
+'''
+
+double3 = numpy.dtype(
+    {
+        'names': ['x', 'y', 'z'],
+        'formats': [numpy.float64]*3
+    }
+)
+
+
+def main():
+    N = 8
+
+    # The kernel computes out = lhs+rhs where lhs and rhs are double3 vectors.
+    # lhs is an array of N such vectors and rhs is double3 kernel parameter.
+
+    lhs = cupy.random.rand(3*N, dtype=numpy.float64).reshape(N, 3)
+    rhs = numpy.random.rand(3).astype(numpy.float64)
+    out = cupy.empty_like(lhs)
+
+    kernel = cupy.RawKernel(code, 'sum_kernel')
+    args = (lhs, rhs.view(double3), out)
+    kernel((1,), (N,), args)
+
+    expected = lhs + cupy.asarray(rhs[None, :])
+    cupy.testing.assert_array_equal(expected, out)
+    print("Kernel output matches expected value.")
+
+
+if __name__ == '__main__':
+    sys.exit(main())
diff --git a/examples/custom_struct/complex_struct.py b/examples/custom_struct/complex_struct.py
new file mode 100644
index 00000000000..227cf8b3154
--- /dev/null
+++ b/examples/custom_struct/complex_struct.py
@@ -0,0 +1,132 @@
+import sys
+import numpy
+import cupy
+
+struct_definition = '''
+struct complex_struct {
+    int4 a;
+    char b;
+    double c[2];
+    short1 d;
+    unsigned long long int e[3];
+};
+'''
+
+struct_layout_code = '''
+{struct_definition}
+
+extern "C" __global__ void get_struct_layout(
+                                unsigned long long *itemsize,
+                                unsigned long long *sizes,
+                                unsigned long long *offsets) {{
+    const complex_struct* ptr = NULL;
+
+    itemsize[0] = sizeof(complex_struct);
+
+    sizes[0] = sizeof(ptr->a);
+    sizes[1] = sizeof(ptr->b);
+    sizes[2] = sizeof(ptr->c);
+    sizes[3] = sizeof(ptr->d);
+    sizes[4] = sizeof(ptr->e);
+
+    offsets[0] = (unsigned long long)&ptr->a;
+    offsets[1] = (unsigned long long)&ptr->b;
+    offsets[2] = (unsigned long long)&ptr->c;
+    offsets[3] = (unsigned long long)&ptr->d;
+    offsets[4] = (unsigned long long)&ptr->e;
+}}
+'''.format(struct_definition=struct_definition)
+
+
+kernel_code = '''
+{struct_definition}
+
+extern "C" __global__ void test_kernel(const complex_struct s,
+                                       double* out) {{
+    int i = threadIdx.x;
+    double sum = 0.0;
+    sum += s.a.x + s.a.y + s.a.z + s.a.w;
+    sum += s.b;
+    sum += s.c[0] + s.c[1];
+    sum += s.d.x;
+    sum += s.e[0] + s.e[1] + s.e[2];
+    out[i] = i * sum;
+}}
+'''.format(struct_definition=struct_definition)
+
+
+def make_packed(basetype, N, itemsize):
+    # A small utility function to make packed structs
+    # Can represent simple packed vectors such as float4 or double[3].
+    assert 0 < N <= 4, N
+    names = list('xyzw')[:N]
+    formats = [basetype]*N
+    return numpy.dtype(dict(names=names,
+                            formats=formats,
+                            itemsize=itemsize))
+
+
+def main():
+    # This program demonstrate how to build a hostside
+    # representation of device structure 'complex_struct'
+    # defined in variable 'struct_definition' that can be
+    # used as a RawKernel argument.
+
+    # First step is to determine structure memory layout
+    #  itemsize -> overall struct size
+    #  sizes    -> individual struct member sizes, determined with sizeof
+    #  offsets  -> individual struct member offsets, determined with offsetof
+    # Results (in terms of bytes) are copied to host after kernel launch.
+    # Note that 'complex_struct' has 5 members named a, b, c, d and e.
+    itemsize = cupy.ndarray(shape=(1,), dtype=numpy.uint64)
+    sizes = cupy.ndarray(shape=(5,), dtype=numpy.uint64)
+    offsets = cupy.ndarray(shape=(5,), dtype=numpy.uint64)
+
+    kernel = cupy.RawKernel(struct_layout_code, 'get_struct_layout')
+    kernel((1,), (1,), (itemsize, sizes, offsets))
+
+    (itemsize, sizes, offsets) = map(cupy.asnumpy, (itemsize, sizes, offsets))
+    print("Overall structure itemsize: {} bytes".format(itemsize.item()))
+    print("Structure members itemsize: {}".format(sizes))
+    print("Structure members offsets: {}".format(offsets))
+
+    # Second step: build a numpy dtype for each struct member
+    atype = make_packed(numpy.int32,   4, sizes[0])
+    btype = make_packed(numpy.int8,    1, sizes[1])
+    ctype = make_packed(numpy.float64, 2, sizes[2])
+    dtype = make_packed(numpy.int16,   1, sizes[3])
+    etype = make_packed(numpy.uint64,  3, sizes[4])
+
+    # Third step: create the complex struct representation with
+    #  the right offsets
+    names = list('abcde')
+    formats = [atype, btype, ctype, dtype, etype]
+    complex_struct = numpy.dtype(dict(names=names,
+                                      formats=formats,
+                                      offsets=offsets,
+                                      itemsize=itemsize.item()))
+
+    # Build a complex_struct kernel argument
+    s = numpy.empty(shape=(1,), dtype=complex_struct)
+    s['a'] = numpy.arange(0, 4).astype(numpy.int32).view(atype)
+    s['b'] = numpy.arange(4, 5).astype(numpy.int8).view(btype)
+    s['c'] = numpy.arange(5, 7).astype(numpy.float64).view(ctype)
+    s['d'] = numpy.arange(7, 8).astype(numpy.int16).view(dtype)
+    s['e'] = numpy.arange(8, 11).astype(numpy.uint64).view(etype)
+    print("Complex structure value:\n  {}".format(s))
+
+    # Setup test kernel
+    N = 8
+    out = cupy.empty(shape=(N,), dtype=numpy.float64)
+    kernel = cupy.RawKernel(kernel_code, 'test_kernel')
+    kernel((1,), (N,), (s, out))
+
+    # the sum of all members of our complex struct instance is 55.0
+    expected = cupy.arange(N) * 55.0
+
+    cupy.testing.assert_array_almost_equal(expected, out)
+    print("Kernel output matches expected value.")
+
+
+if __name__ == '__main__':
+    sys.exit(main())
diff --git a/examples/custom_struct/packed_matrix.py b/examples/custom_struct/packed_matrix.py
new file mode 100644
index 00000000000..beb7de45278
--- /dev/null
+++ b/examples/custom_struct/packed_matrix.py
@@ -0,0 +1,91 @@
+import sys
+import numpy
+import cupy
+
+code = '''
+template<typename T>
+struct Matrix {
+    T value[4][4];
+
+    __device__ T& operator() (int i, int j) {
+        return this->value[i][j];
+    }
+
+    __device__ const T& operator() (int i, int j) const {
+        return this->value[i][j];
+    }
+};
+
+template<typename T>
+__device__ Matrix<T> operator+ (const Matrix<T>& lhs, const Matrix<T>& rhs) {
+    Matrix<T> res;
+    for (int i = 0; i<4; i++) {
+        for (int j = 0; j<4; j++) {
+            res(i,j) = lhs(i,j) + rhs(i,j);
+        }
+    }
+    return res;
+}
+
+template<typename T>
+__device__ Matrix<T> operator* (const Matrix<T>& lhs, const Matrix<T>& rhs) {
+    Matrix<T> res;
+    for (int i = 0; i<4; i++) {
+        for (int j = 0; j<4; j++) {
+            res(i,j) = T(0);
+            for (int k = 0; k<4; k++) {
+                res(i,j) += lhs(i,k) * rhs(k,j);
+            }
+        }
+    }
+    return res;
+}
+
+template<typename T>
+__global__ void kernel(const Matrix<T>* A,
+                       const Matrix<T>* B,
+                       const Matrix<T>  C,
+                             Matrix<T>* out) {
+  int i = threadIdx.x;
+  out[i] = A[i] * B[i] + C;
+}
+'''
+
+
+def main():
+    N = 8
+    module = cupy.RawModule(code=code, options=('-std=c++03',),
+                            name_expressions=('kernel<float>',
+                                              'kernel<double>'))
+
+    # The kernel computes out = A*B+C where A, B and C are 4x4 matrices.
+    # A and B are arrays of N such matrices and C is a matrix kernel parameter.
+
+    for (ctype, dtype) in zip(('float', 'double'),
+                              (numpy.float32, numpy.float64)):
+
+        A = cupy.random.rand(16*N, dtype=dtype).reshape(N, 4, 4)
+        B = cupy.random.rand(16*N, dtype=dtype).reshape(N, 4, 4)
+        C = numpy.random.rand(16).astype(dtype).reshape(4, 4)
+        out = cupy.empty_like(A)
+
+        Matrix = numpy.dtype(
+            {
+                'names': ['value'],
+                'formats': [(dtype, (4, 4))]
+            }
+        )
+
+        kernel = module.get_function('kernel<{}>'.format(ctype))
+        args = (A, B, C.ravel().view(Matrix), out)
+        kernel((1,), (N,), args)
+
+        expected = cupy.matmul(A, B) + cupy.asarray(C[None, :, :])
+
+        cupy.testing.assert_array_almost_equal(expected, out)
+        print("Kernel output matches expected value for "
+              "type '{}'.".format(ctype))
+
+
+if __name__ == '__main__':
+    sys.exit(main())
diff --git a/tests/cupy_tests/core_tests/test_function.py b/tests/cupy_tests/core_tests/test_function.py
index 7765e37ce7e..f2970e5f2f1 100644
--- a/tests/cupy_tests/core_tests/test_function.py
+++ b/tests/cupy_tests/core_tests/test_function.py
@@ -1,14 +1,19 @@
 import unittest
+import pytest
 
 import numpy
 
 import cupy
+from cupy.core import core
 from cupy.cuda import compiler
+from cupy.cuda import runtime
 from cupy import testing
 
 
 def _compile_func(kernel_name, code):
-    mod = compiler.compile_with_cache(code)
+    # workaround for hipRTC
+    extra_source = core._get_header_source() if runtime.is_hip else None
+    mod = compiler.compile_with_cache(code, extra_source=extra_source)
     return mod.get_function(kernel_name)
 
 
@@ -54,3 +59,160 @@ def test_numpy_scalar(self):
 
         expected = a_cpu + b
         testing.assert_array_equal(x, expected)
+
+    def test_numpy_dtype(self):
+        code = '''
+extern "C" __global__ void test_kernel(const double* a, double3 b, double* x) {
+  int i = threadIdx.x;
+  x[i] = a[i] + b.x + b.y + b.z;
+}
+'''
+
+        double3 = numpy.dtype({'names': ['x', 'y', 'z'],
+                               'formats': [numpy.float64]*3})
+        a_cpu = numpy.arange(24, dtype=numpy.float64)
+        a = cupy.array(a_cpu)
+        b = numpy.random.rand(3).view(double3)
+        x = cupy.empty_like(a)
+
+        func = _compile_func('test_kernel', code)
+
+        func.linear_launch(a.size, (a, b, x))
+
+        expected = a_cpu + b['x'] + b['y'] + b['z']
+        testing.assert_array_equal(x, expected)
+
+    def test_static_array(self):
+        code = '''
+struct double5 {
+    double value[5];
+    __device__ const double& operator[](size_t i) const { return value[i]; }
+};
+
+extern "C" __global__ void test_kernel(const double* a, double5 b, double* x) {
+  int i = threadIdx.x;
+  x[i] = a[i] + b[0] + b[1] + b[2] + b[3] + b[4];
+}
+'''
+
+        a_cpu = numpy.arange(24, dtype=numpy.float64)
+        a = cupy.array(a_cpu)
+        x = cupy.empty_like(a)
+
+        func = _compile_func('test_kernel', code)
+
+        # We cannot pass np.ndarray kernel arguments of size > 1
+        b = numpy.arange(5).astype(numpy.float64)
+        with pytest.raises(TypeError):
+            func.linear_launch(a.size, (a, b, x))
+
+        double5 = numpy.dtype({
+            'names': ['dummy'],
+            'formats': [(numpy.float64, (5,))]
+        })
+        func.linear_launch(a.size, (a, b.view(double5), x))
+
+        expected = a_cpu + b.sum()
+        testing.assert_array_equal(x, expected)
+
+    def test_custom_user_struct(self):
+        struct_definition = '''
+struct custom_user_struct {
+    int4 a;
+    char b;
+    double c[2];
+    short1 d;
+    unsigned long long int e[3];
+};
+'''
+
+        # first step is to determine struct memory layout
+        struct_layout_code = '''
+{struct_definition}
+extern "C" __global__ void get_struct_layout(
+                                unsigned long long *itemsize,
+                                unsigned long long *sizes,
+                                unsigned long long *offsets) {{
+    const custom_user_struct* ptr = NULL;
+
+    itemsize[0] = sizeof(custom_user_struct);
+
+    sizes[0] = sizeof(ptr->a);
+    sizes[1] = sizeof(ptr->b);
+    sizes[2] = sizeof(ptr->c);
+    sizes[3] = sizeof(ptr->d);
+    sizes[4] = sizeof(ptr->e);
+
+    offsets[0] = (unsigned long long)&ptr->a;
+    offsets[1] = (unsigned long long)&ptr->b;
+    offsets[2] = (unsigned long long)&ptr->c;
+    offsets[3] = (unsigned long long)&ptr->d;
+    offsets[4] = (unsigned long long)&ptr->e;
+}}
+'''.format(struct_definition=struct_definition)
+
+        itemsize = cupy.ndarray(shape=(1,), dtype=numpy.uint64)
+        sizes = cupy.ndarray(shape=(5,), dtype=numpy.uint64)
+        offsets = cupy.ndarray(shape=(5,), dtype=numpy.uint64)
+        func = _compile_func('get_struct_layout', struct_layout_code)
+        func.linear_launch(1, (itemsize, sizes, offsets))
+
+        # Build structure data type recursively
+        names = list('abcde')
+        itemsize = cupy.asnumpy(itemsize).item()
+        sizes = cupy.asnumpy(sizes).tolist()
+        offsets = cupy.asnumpy(offsets).tolist()
+
+        def make_packed(basetype, N, itemsize):
+            assert 0 < N <= 4, N
+            names = list('xyzw')[:N]
+            formats = [basetype]*N
+            return numpy.dtype(dict(names=names,
+                                    formats=formats,
+                                    itemsize=itemsize))
+
+        # structure member data types
+        int4 = make_packed(numpy.int32, 4, sizes[0])
+        char = make_packed(numpy.int8, 1, sizes[1])
+        double2 = make_packed(numpy.float64, 2, sizes[2])
+        short1 = make_packed(numpy.int16, 1, sizes[3])
+        ulong3 = make_packed(numpy.uint64, 3, sizes[4])
+
+        formats = [int4, char, double2, short1, ulong3]
+        struct_dtype = numpy.dtype(dict(names=names,
+                                        formats=formats,
+                                        offsets=offsets,
+                                        itemsize=itemsize))
+
+        s = numpy.empty(shape=(1,), dtype=struct_dtype)
+        s['a'] = numpy.arange(0, 4).astype(numpy.int32).view(int4)
+        s['b'] = numpy.arange(4, 5).astype(numpy.int8).view(char)
+        s['c'] = numpy.arange(5, 7).astype(numpy.float64).view(double2)
+        s['d'] = numpy.arange(7, 8).astype(numpy.int16).view(short1)
+        s['e'] = numpy.arange(8, 11).astype(numpy.uint64).view(ulong3)
+
+        # test kernel code
+        code = '''
+{struct_definition}
+extern "C" __global__ void test_kernel(const double* a,
+                                       custom_user_struct s,
+                                       double* x) {{
+    int i = threadIdx.x;
+    double sum = s.a.x + s.a.y + s.a.z + s.a.w;
+    sum += s.b;
+    sum += s.c[0] + s.c[1];
+    sum += s.d.x;
+    sum += s.e[0] + s.e[1] + s.e[2];
+    x[i] = a[i] + sum;
+}}
+'''.format(struct_definition=struct_definition)
+
+        a_cpu = numpy.arange(24, dtype=numpy.float64)
+        a = cupy.array(a_cpu)
+        x = cupy.empty_like(a)
+
+        func = _compile_func('test_kernel', code)
+        func.linear_launch(a.size, (a, s, x))
+
+        expected = a_cpu + 55.0
+        testing.assert_array_equal(x, expected)
diff --git a/tests/example_tests/test_custom_struct.py b/tests/example_tests/test_custom_struct.py
new file mode 100644
index 00000000000..1d4e5a662bc
--- /dev/null
+++ b/tests/example_tests/test_custom_struct.py
@@ -0,0 +1,34 @@
+import re
+import unittest
+
+from example_tests import example_test
+
+
+class TestCustomStruct(unittest.TestCase):
+
+    def test_builtin_vectors(self):
+        output = example_test.run_example('custom_struct/builtin_vectors.py')
+        assert re.match(
+            r"Kernel output matches expected value.",
+            output.decode('utf-8')
+        )
+
+    def test_packed_matrix(self):
+        output = example_test.run_example('custom_struct/packed_matrix.py')
+        assert re.match(
+            r"Kernel output matches expected value for type 'float'.\n"
+            r"Kernel output matches expected value for type 'double'.",
+            output.decode('utf-8'),
+        )
+
+    def test_complex_struct(self):
+        output = example_test.run_example('custom_struct/complex_struct.py')
+        assert re.match(
+            r"Overall structure itemsize: \d+ bytes\n"
+            r"Structure members itemsize: \[(\s*\d+){5}]\n"
+            r"Structure members offsets: \[(\s*\d+){5}]\n"
+            r"Complex structure value:\n"
+            r"\s+\[.*\]\n"
+            r"Kernel output matches expected value.",
+            output.decode('utf-8'),
+        )


Code review comments:
diff_hunk: @@ -0,0 +1,91 @@
+import sys
+import numpy
+import cupy
+
+code = '''
+template<typename T>
+struct Matrix {
+    T value[4][4];
+
+    __device__ T& operator() (int i, int j) {
+        return this->value[i][j];
+    }
+
+    __device__ const T& operator() (int i, int j) const {
+        return this->value[i][j];
+    }
+};
+
+template<typename T>
+__device__ Matrix<T> operator+ (const Matrix<T>& lhs, const Matrix<T>& rhs) {
+    Matrix<T> res;
+    for (int i = 0; i<4; i++) {
+        for (int j = 0; j<4; j++) {
+            res(i,j) = lhs(i,j) + rhs(i,j);
+        }
+    }
+    return res;
+}
+
+template<typename T>
+__device__ Matrix<T> operator* (const Matrix<T>& lhs, const Matrix<T>& rhs) {
+    Matrix<T> res;
+    for (int i = 0; i<4; i++) {
+        for (int j = 0; j<4; j++) {
+            res(i,j) = T(0);
+            for (int k = 0; k<4; k++) {
+                res(i,j) += lhs(i,k) * rhs(k,j);
+            }
+        }
+    }
+    return res;
+}
+
+template<typename T>
+__global__ void kernel(const Matrix<T>* A,
+                       const Matrix<T>* B,
+                       const Matrix<T>  C,
+                             Matrix<T>* out) {
+  int i = threadIdx.x;
+  out[i] = A[i] * B[i] + C;
+}
+'''
+
+
+def main():
+    N = 8
+    module = cupy.RawModule(code=code, options=('-std=c++03',),
path: examples/custom_struct/packed_matrix.py
text: @user1:
Apparently we hit the same issue in two different PRs...See https://github.com/cupy/cupy/pull/4922#issuecomment-803718555.
```suggestion
    module = cupy.RawModule(code=code, options=('-std=c++11',),
```

@author:
I was not aware that C++03 was CUDA 11 specific, hope CUDA 12 will bring support for C++98 !