The AI CUDA Engineer 👷

import torch
import torch.nn as nn
import torch.nn.functional as F


def module_fn(
    x: torch.Tensor,
    stride: int,
    padding: int,
    conv_transpose: torch.Tensor,
    conv_transpose_bias: torch.Tensor,
) -> torch.Tensor:
    """
    Applies a 3D transposed convolution operation followed by two max pooling layers and a sum operation.

    Args:
        x (torch.Tensor): Input tensor of shape (batch_size, in_channels, depth, height, width)
        stride (int): Stride of the transposed convolution
        padding (int): Padding of the transposed convolution
        conv_transpose (torch.Tensor): Transposed convolution weight tensor
        conv_transpose_bias (torch.Tensor): Bias tensor for transposed convolution

    Returns:
        torch.Tensor: Output tensor after applying transposed convolution, max pooling and sum reduction
    """
    x = F.conv_transpose3d(
        x, conv_transpose, bias=conv_transpose_bias, stride=stride, padding=padding
    )
    x = F.max_pool3d(x, kernel_size=2)
    x = F.max_pool3d(x, kernel_size=3)
    x = torch.sum(x, dim=1, keepdim=True)
    return x


class Model(nn.Module):
    """
    Model that performs a 3D transposed convolution, followed by two max pooling layers and a sum operation.
    """

    def __init__(self, in_channels, out_channels, kernel_size, stride, padding):
        super(Model, self).__init__()
        conv = nn.ConvTranspose3d(in_channels, out_channels, kernel_size)
        self.conv_transpose_parameter = nn.Parameter(conv.weight)
        self.conv_transpose_bias = nn.Parameter(conv.bias)

    def forward(self, x, stride, padding, fn=module_fn):
        return fn(
            x, stride, padding, self.conv_transpose_parameter, self.conv_transpose_bias
        )


batch_size = 16
in_channels = 8
out_channels = 16
depth, height, width = 16, 32, 32
kernel_size = 3
stride = 2
padding = 1


def get_inputs():
    return [torch.randn(batch_size, in_channels, depth, height, width), stride, padding]


def get_init_inputs():
    return [in_channels, out_channels, kernel_size, stride, padding]

#include <torch/extension.h>
#include <ATen/ATen.h>
#include <cuda.h>
#include <cuda_runtime.h>
#include <cfloat>

// Device function: 2x2x2 max pooling on the input tensor
__device__ inline float pool2x2x2(const float* input, int n, int c,
                                    int start_d1, int start_h1, int start_w1,
                                    int D1, int H1, int W1, int C) {
    float max_val = -FLT_MAX;
    #pragma unroll
    for (int dd = 0; dd < 2; dd++) {
        int d1 = start_d1 + dd;
        if (d1 < D1) {
            #pragma unroll
            for (int hh = 0; hh < 2; hh++) {
                int h1 = start_h1 + hh;
                if (h1 < H1) {
                    #pragma unroll
                    for (int ww = 0; ww < 2; ww++) {
                        int w1 = start_w1 + ww;
                        if (w1 < W1) {
                            int idx = ((n * C + c) * D1 + d1) * H1 * W1 + h1 * W1 + w1;
                            max_val = max(max_val, input[idx]);
                        }
                    }
                }
            }
        }
    }
    return max_val;
}

// Device function: 3x3x3 pooling over the results of the 2x2x2 pooling
__device__ inline float pool3x3x3(const float* input, int n, int c,
                                    int start_d2, int start_h2, int start_w2,
                                    int D1, int H1, int W1, int C) {
    float max_val = -FLT_MAX;
    // The intermediate dimensions from 2x2x2 pooling are D2 = D1/2, H2 = H1/2, W2 = W1/2
    #pragma unroll
    for (int d_off = 0; d_off < 3; d_off++) {
        int d2 = start_d2 + d_off;
        if (d2 < (D1 / 2)) {
            #pragma unroll
            for (int h_off = 0; h_off < 3; h_off++) {
                int h2 = start_h2 + h_off;
                if (h2 < (H1 / 2)) {
                    #pragma unroll
                    for (int w_off = 0; w_off < 3; w_off++) {
                        int w2 = start_w2 + w_off;
                        if (w2 < (W1 / 2)) {
                            // Compute starting indices for the 2x2x2 pooling window
                            int start_d1 = d2 * 2;
                            int start_h1 = h2 * 2;
                            int start_w1 = w2 * 2;
                            float temp = pool2x2x2(input, n, c, start_d1, start_h1, start_w1, D1, H1, W1, C);
                            max_val = max(max_val, temp);
                        }
                    }
                }
            }
        }
    }
    return max_val;
}

// Main kernel: Each thread computes one element of the final pooled output
__global__ void modular_maxpool_kernel(
    const float* __restrict__ input,
    float* __restrict__ output,
    const int N, const int C,
    const int D1, const int H1, const int W1,
    const int D3, const int H3, const int W3) {

    int idx = blockIdx.x * blockDim.x + threadIdx.x;
    if (idx >= N * C * D3 * H3 * W3) return;

    // Decode output index (n, c, d3, h3, w3)
    int w3 = idx % W3;
    int h3 = (idx / W3) % H3;
    int d3 = (idx / (W3 * H3)) % D3;
    int c  = (idx / (W3 * H3 * D3)) % C;
    int n  = idx / (W3 * H3 * D3 * C);

    // Calculate starting indices for the 3x3x3 pooling window on the intermediate tensor
    // The intermediate tensor dimensions are: D2 = D1/2, H2 = H1/2, W2 = W1/2.
    // The second maxpool has kernel size 3 and stride 3, so we compute:
    int start_d2 = d3 * 3;
    int start_h2 = h3 * 3;
    int start_w2 = w3 * 3;

    // Compute the final pooled value using the modular device function
    float val = pool3x3x3(input, n, c, start_d2, start_h2, start_w2, D1, H1, W1, C);
    output[idx] = val;
}

// The forward function applies a conv_transpose3d operation, then the modular max pooling.
// It sums over the channel dimension at the end.

torch::Tensor forward(
    torch::Tensor x,
    int64_t stride,
    int64_t padding,
    torch::Tensor conv_transpose,
    torch::Tensor conv_transpose_bias) {

    x = x.contiguous();
    conv_transpose = conv_transpose.contiguous();
    conv_transpose_bias = conv_transpose_bias.contiguous();

    TORCH_CHECK(x.is_cuda(), "Input x must be a CUDA tensor");
    TORCH_CHECK(conv_transpose.is_cuda(), "conv_transpose must be a CUDA tensor");
    TORCH_CHECK(conv_transpose_bias.is_cuda(), "conv_transpose_bias must be a CUDA tensor");

    // Apply the transposed convolution operation
    x = at::conv_transpose3d(
        x,
        conv_transpose,
        conv_transpose_bias,
        {stride, stride, stride},
        {padding, padding, padding}
    );

    // Get dimensions after conv_transpose: [N, C, D1, H1, W1]
    auto sizes = x.sizes();
    const int N = sizes[0];
    const int C = sizes[1];
    const int D1 = sizes[2];
    const int H1 = sizes[3];
    const int W1 = sizes[4];

    // The two-stage pooling reduces dimensions as follows:
    // First max pool (2x2x2): D2 = D1 / 2, H2 = H1 / 2, W2 = W1 / 2
    // Second max pool (3x3x3): D3 = D2 / 3 = D1 / 6, H3 = H1 / 6, W3 = W1 / 6
    const int D3 = D1 / 6;
    const int H3 = H1 / 6;
    const int W3 = W1 / 6;

    auto output = torch::empty({N, C, D3, H3, W3}, x.options());

    // Launch the kernel
    const int total_elements = N * C * D3 * H3 * W3;
    const int threads = 256;
    const int blocks = (total_elements + threads - 1) / threads;

    modular_maxpool_kernel<<<blocks, threads>>>(
        x.data_ptr<float>(),
        output.data_ptr<float>(),
        N, C, D1, H1, W1,
        D3, H3, W3
    );

    // Sum over channels (dim=1), keeping the dimension
    return output.sum(1, /*keepdim=*/true);
}

PYBIND11_MODULE(TORCH_EXTENSION_NAME, m) {
    m.def("forward", &forward, "Modular forward pass with device functions for max pooling");
}