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_runtime.h>

#define WARP_SIZE 32

// Optimized kernel with balanced workload distribution
// Grid mapping: grid.x = N * D, grid.y = H, grid.z = W
// Each block computes the sum over channels for one spatial location of one (n,d) pair.

template <int BLOCK_SIZE>
__global__ void sum_channels_kernel(const float* __restrict__ input, float* __restrict__ output,
                                      int N, int C, int D, int H, int W) {
    // Compute indices from grid: combine N and D into grid.x
    int nd = blockIdx.x;            // nd ranges from 0 to (N*D - 1)
    int n = nd / D;                 // n index
    int d = nd % D;                 // d index
    int h = blockIdx.y;             // h index
    int w = blockIdx.z;             // w index

    // Compute the base offset for (n, d, h, w) in a contiguous tensor [N, C, D, H, W]
    // Offset for element [n, c, d, h, w] = n * (C*D*H*W) + c * (D*H*W) + d * (H*W) + h * W + w
    int base = n * (C * D * H * W) + d * (H * W) + h * W + w;

    float thread_sum = 0.0f;

    // Loop over channels with stride equal to BLOCK_SIZE
    for (int c = threadIdx.x; c < C; c += BLOCK_SIZE) {
        thread_sum += input[base + c * (D * H * W)];
    }

    // Warp-level reduction using shuffle intrinsics
    unsigned int mask = 0xffffffff;
    for (int offset = WARP_SIZE / 2; offset > 0; offset /= 2) {
        thread_sum += __shfl_down_sync(mask, thread_sum, offset);
    }

    // Use shared memory to accumulate results from different warps if needed
    __shared__ float warp_sums[BLOCK_SIZE / WARP_SIZE];
    int lane = threadIdx.x % WARP_SIZE;
    int warp_id = threadIdx.x / WARP_SIZE;
    if (lane == 0) {
        warp_sums[warp_id] = thread_sum;
    }
    __syncthreads();

    // Final reduction from warp sums: let the first warp perform the reduction
    float block_sum = 0.0f;
    if (threadIdx.x < (BLOCK_SIZE / WARP_SIZE)) {
        block_sum = warp_sums[threadIdx.x];
    }
    if (threadIdx.x < WARP_SIZE) {
        for (int offset = (BLOCK_SIZE / WARP_SIZE) / 2; offset > 0; offset /= 2) {
            block_sum += __shfl_down_sync(mask, block_sum, offset);
        }
    }
    
    if (threadIdx.x == 0) {
        // Output tensor has shape [N, 1, D, H, W] so index accordingly
        int out_index = n * (D * H * W) + d * (H * W) + h * W + w;
        output[out_index] = block_sum;
    }
}

// Forward function that applies conv_transpose3d, pooling, then channel sum reduction
torch::Tensor forward(
    torch::Tensor x,
    int64_t stride,
    int64_t padding,
    torch::Tensor conv_transpose,
    torch::Tensor conv_transpose_bias) {

    // Ensure inputs are contiguous
    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 transposed convolution
    x = at::conv_transpose3d(
        x,
        conv_transpose,
        conv_transpose_bias,
        {stride, stride, stride},
        {padding, padding, padding}
    );

    // Apply max pooling operations
    x = at::max_pool3d(x, {2, 2, 2});
    x = at::max_pool3d(x, {3, 3, 3});

    // Prepare output tensor with summed channels (result shape: [N, 1, D, H, W])
    const auto sizes = x.sizes(); // sizes: [N, C, D, H, W]
    auto output = torch::zeros({sizes[0], 1, sizes[2], sizes[3], sizes[4]}, x.options());

    // Launch kernel: use grid dimensions mapping (N*D, H, W) to cover the output spatial domain
    const int BLOCK_SIZE = 256;
    dim3 grid(sizes[0] * sizes[2], sizes[3], sizes[4]);
    sum_channels_kernel<BLOCK_SIZE><<<grid, BLOCK_SIZE>>>(
        x.data_ptr<float>(),
        output.data_ptr<float>(),
        sizes[0],
        sizes[1],
        sizes[2],
        sizes[3],
        sizes[4]
    );

    return output;
}

PYBIND11_MODULE(TORCH_EXTENSION_NAME, m) {
    m.def("forward", &forward, "Optimized 78_ConvTranspose3d_Max_Max_Sum kernel with balanced workload distribution");
}