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,
    output_padding: int,
    conv_transpose: torch.Tensor,
    conv_transpose_bias: torch.Tensor,
    multiplier: float,
) -> torch.Tensor:
    """
    Applies transposed convolution, scalar multiplication, and multiple global average pooling operations.

    Args:
        x (torch.Tensor): Input tensor of shape (batch_size, in_channels, height, width)
        stride (int): Stride of the transposed convolution
        padding (int): Padding of the transposed convolution
        output_padding (int): Additional size added to output shape
        conv_transpose (torch.Tensor): Transposed convolution weight tensor
        conv_transpose_bias (torch.Tensor): Bias tensor for transposed convolution
        multiplier (float): Scalar multiplier value

    Returns:
        torch.Tensor: Scalar output after applying operations
    """
    x = F.conv_transpose2d(
        x,
        conv_transpose,
        bias=conv_transpose_bias,
        stride=stride,
        padding=padding,
        output_padding=output_padding,
    )
    x = x * multiplier
    x = torch.mean(x, dim=[2, 3], keepdim=True)
    x = torch.mean(x, dim=[2, 3], keepdim=True)
    x = torch.mean(x)
    return x


class Model(nn.Module):
    """
    Model that performs a transposed convolution, multiplies by a scalar, applies global average pooling,
    another global average pooling, and then calculates the mean.
    """

    def __init__(
        self,
        in_channels,
        out_channels,
        kernel_size,
        stride,
        padding,
        output_padding,
        multiplier,
    ):
        super(Model, self).__init__()
        conv = nn.ConvTranspose2d(
            in_channels,
            out_channels,
            kernel_size,
            stride=stride,
            padding=padding,
            output_padding=output_padding,
        )
        self.conv_transpose_parameter = nn.Parameter(conv.weight)
        self.conv_transpose_bias = nn.Parameter(
            conv.bias
            + torch.randn(
                conv.bias.shape, device=conv.bias.device, dtype=conv.bias.dtype
            )
            * 0.02
        )
        self.multiplier = multiplier

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


batch_size = 128
in_channels = 3
out_channels = 16
height, width = 32, 32
kernel_size = 3
stride = 2
padding = 1
output_padding = 1
multiplier = 0.5


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


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

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

// Kernel to compute the mean of each (batch, channel) slice using vectorized loads with __ldg()
// and then atomically accumulate the slice means into a global accumulator.

// The kernel uses 128-bit (float4) vectorized loads when the number of elements is divisible by 4,
// otherwise it falls back to scalar loads with __ldg().

// Each block is responsible for one (batch, channel) slice.
// The final result is computed as the average over all (batch, channel) means.

template <unsigned int blockSize>
__global__ void vectorized_ldg_mean_kernel(
    const float* __restrict__ input,
    float* __restrict__ global_accum,
    int H,
    int W,
    int C
) {
    extern __shared__ float shared[]; // shared memory for block reduction
    int num_elements = H * W;
    int batch = blockIdx.x / C;
    int channel = blockIdx.x % C;
    int input_offset = (batch * C + channel) * num_elements;
    float sum = 0.0f;

    // Check if the slice size is divisible by 4 to use vectorized loads
    if ((num_elements & 3) == 0) {
        int num_vec = num_elements >> 2; // equivalent to num_elements / 4
        // Cast the input pointer to float4 pointer for 128-bit aligned loads
        const float4* in_vec = reinterpret_cast<const float4*>(input + input_offset);
        for (int i = threadIdx.x; i < num_vec; i += blockDim.x) {
            // Use __ldg for read-only global memory access
            float4 v = __ldg(&in_vec[i]);
            sum += v.x + v.y + v.z + v.w;
        }
    } else {
        for (int i = threadIdx.x; i < num_elements; i += blockDim.x) {
            sum += __ldg(&input[input_offset + i]);
        }
    }

    // Write the partial sum to shared memory
    shared[threadIdx.x] = sum;
    __syncthreads();

    // Intra-block reduction in shared memory
    if (blockSize >= 512) {
        if (threadIdx.x < 256) {
            shared[threadIdx.x] += shared[threadIdx.x + 256];
        }
        __syncthreads();
    }
    if (blockSize >= 256) {
        if (threadIdx.x < 128) {
            shared[threadIdx.x] += shared[threadIdx.x + 128];
        }
        __syncthreads();
    }
    if (blockSize >= 128) {
        if (threadIdx.x < 64) {
            shared[threadIdx.x] += shared[threadIdx.x + 64];
        }
        __syncthreads();
    }
    if (threadIdx.x < 32) {
        volatile float* vsmem = shared;
        vsmem[threadIdx.x] += vsmem[threadIdx.x + 32];
        vsmem[threadIdx.x] += vsmem[threadIdx.x + 16];
        vsmem[threadIdx.x] += vsmem[threadIdx.x + 8];
        vsmem[threadIdx.x] += vsmem[threadIdx.x + 4];
        vsmem[threadIdx.x] += vsmem[threadIdx.x + 2];
        vsmem[threadIdx.x] += vsmem[threadIdx.x + 1];
    }

    // Thread 0 computes the mean for this (batch, channel) slice and atomically adds it to the global accumulator
    if (threadIdx.x == 0) {
        float slice_mean = shared[0] / static_cast<float>(num_elements);
        atomicAdd(global_accum, slice_mean);
    }
}

at::Tensor module_fn(
    at::Tensor x,
    int64_t stride,
    int64_t padding,
    int64_t output_padding,
    at::Tensor conv_transpose,
    at::Tensor conv_transpose_bias,
    double multiplier
) {
    // Perform transposed convolution
    at::Tensor y = at::conv_transpose2d(
        x,
        conv_transpose,
        conv_transpose_bias,
        {stride, stride},
        {padding, padding},
        {output_padding, output_padding},
        1,
        {1, 1}
    );
    
    // Scale the output
    y = y * multiplier;
    
    // Get dimensions (N, C, H, W)
    auto dims = y.sizes();
    int N = dims[0];
    int C = dims[1];
    int H = dims[2];
    int W = dims[3];
    
    // Allocate a scalar accumulator on the device initialized to zero
    auto options = torch::TensorOptions().device(y.device()).dtype(y.dtype());
    at::Tensor accum = torch::zeros({1}, options);
    
    // Launch one block per (batch, channel) slice
    constexpr int blockSize = 256;
    int numBlocks = N * C;
    size_t sharedMemSize = blockSize * sizeof(float);

    vectorized_ldg_mean_kernel<blockSize><<<numBlocks, blockSize, sharedMemSize>>>(
        y.data_ptr<float>(),
        accum.data_ptr<float>(),
        H, W, C
    );
    
    // Compute the final overall mean: each block contributed its slice mean, so average over N*C slices
    accum = accum / static_cast<float>(N * C);
    return accum;
}

PYBIND11_MODULE(TORCH_EXTENSION_NAME, m) {
    m.def("forward", &module_fn, "Module function");
}