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

// Fused kernel: performs elementwise multiplication then reduction over spatial dims
// Each block processes one (batch, channel) slice of the tensor

template <int BLOCK_SIZE>
__global__ void fused_reduction_kernel(
    const float* __restrict__ input,
    float multiplier,
    int N,
    int C,
    int H,
    int W,
    float* __restrict__ global_sum
) {
    extern __shared__ float sdata[];
    const int tid = threadIdx.x;

    // Each block handles one (batch, channel) pair
    const int bid = blockIdx.x; // range: 0 to (N*C - 1)
    const int batch_idx = bid / C;
    const int channel_idx = bid % C;
    const int spatial_size = H * W;

    // Get pointer to the beginning of the (batch, channel) slice
    const float* input_ptr = input + (batch_idx * C * spatial_size) + (channel_idx * spatial_size);

    // Each thread processes multiple elements along the spatial dimension
    float sum = 0.0f;
    for (int i = tid; i < spatial_size; i += BLOCK_SIZE) {
        // Fuse multiplication by 'multiplier' with accumulation
        sum += input_ptr[i] * multiplier;
    }

    // Store the local sum in shared memory
    sdata[tid] = sum;
    __syncthreads();

    // Perform reduction in shared memory
    if (BLOCK_SIZE >= 512) { 
        if (tid < 256) { sdata[tid] += sdata[tid + 256]; } 
        __syncthreads(); 
    }
    if (BLOCK_SIZE >= 256) { 
        if (tid < 128) { sdata[tid] += sdata[tid + 128]; } 
        __syncthreads(); 
    }
    if (BLOCK_SIZE >= 128) { 
        if (tid < 64) { sdata[tid] += sdata[tid + 64]; } 
        __syncthreads(); 
    }

    // Warp-level reduction (no __syncthreads needed within a warp)
    if (tid < 32) {
        volatile float* vsmem = sdata;
        if (BLOCK_SIZE >= 64) { vsmem[tid] += vsmem[tid + 32]; }
        if (BLOCK_SIZE >= 32) { vsmem[tid] += vsmem[tid + 16]; }
        if (BLOCK_SIZE >= 16) { vsmem[tid] += vsmem[tid + 8]; }
        if (BLOCK_SIZE >= 8)  { vsmem[tid] += vsmem[tid + 4]; }
        if (BLOCK_SIZE >= 4)  { vsmem[tid] += vsmem[tid + 2]; }
        if (BLOCK_SIZE >= 2)  { vsmem[tid] += vsmem[tid + 1]; }
    }

    // Thread 0 atomically adds the block's result to the global accumulator
    if (tid == 0) {
        atomicAdd(global_sum, sdata[0]);
    }
}

// Module function: performs transposed convolution, then fuses elementwise multiplication and reduction
// to compute the overall mean of the output tensor
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
) {
    // Apply transposed convolution (optimized via cuDNN)
    at::Tensor y = at::conv_transpose2d(
        x,
        conv_transpose,
        conv_transpose_bias,
        {stride, stride},
        {padding, padding},
        {output_padding, output_padding},
        1,
        {1, 1}
    );

    // Retrieve dimensions of y: (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 device tensor to hold the global sum (initialized to 0)
    auto options = torch::TensorOptions().device(y.device()).dtype(torch::kFloat32);
    at::Tensor global_sum_tensor = torch::zeros({1}, options);

    // Launch the fused reduction kernel
    constexpr int BLOCK_SIZE = 256;
    int grid = N * C;  // one block per (batch, channel) slice
    int shared_mem_size = BLOCK_SIZE * sizeof(float);
    
    fused_reduction_kernel<BLOCK_SIZE><<<grid, BLOCK_SIZE, shared_mem_size>>>(
        y.data_ptr<float>(),
        static_cast<float>(multiplier),
        N, C, H, W,
        global_sum_tensor.data_ptr<float>()
    );

    // Synchronize to ensure kernel completion
    cudaDeviceSynchronize();

    // Compute the final mean: global sum divided by total number of elements
    float total_elements = static_cast<float>(N * C * H * W);
    float final_sum = global_sum_tensor.item<float>();
    float final_mean = final_sum / total_elements;

    // Return the mean as a scalar tensor
    return torch::tensor(final_mean, options);
}

PYBIND11_MODULE(TORCH_EXTENSION_NAME, m) {
    m.def("forward", &module_fn, "Fused elementwise multiplication and reduction");
}