FW-GLM-5-1 · popcorn · pop_l123_default_gpt
correct 43_ContactMapPrediction.py
model
FW-GLM-5-1
turns / max
12 / 12
tool calls
12
tokens (total)
117004
wall clock
4100.1s
truncated
no
backend
cuda fp32
variant
popcorn
runtime
0.11 μs
ref runtime
0.13 μs
speedup
1.25x
started
2026-05-11T13:25:31Z
timingkernelbench.timing
kernel runtime
0.11 μs
ref runtime
0.13 μs
speedup
1.25x
kernel mean / median / std
0.11 / — / 0.03 μs
ref mean / median / std
0.13 / — / 0.03 μs
numerical precisioneval_kernel_against_ref
max abs error
1.19e-07
mean abs error
2.41e-09
max rel error
2.32e-07
mean rel error
4.69e-09
gpu memoryextended_metrics.measure_memory
kernel peak
45.03 MB
ref peak
49.03 MB
kernel / ref
0.92x
kernel launches / fusiontorch.profiler
kernel launches
32
ref launches
41
fusion ratio (ref / kernel)
1.28x
total cuda time
884 μs
energyNVML (pynvml)
kernel mJ / run
241.38
ref mJ / run
0.00
ratio (ref / kernel)
—
avg power
2195.14 W
ref avg power
0.00 W
SOL (speed-of-light)heuristic_from_speedup
SOL score (0–1)
0.125
dram util %
—
compute util %
—
bottleneck
unknown
dominant pipe
—
arithmetic intensity
—
ridge point
—
roofline / occupancyheuristic
dram bandwidth
—
dram util %
—
fp32 tflops
—
fp32 util %
—
fp16 tflops
—
occupancy %
—
bottleneck
—
dominant pipe %
—
arithmetic intensity
—
ridge point
—
L1 hit %
—
L2 hit %
—
ld sectors / req
—
st sectors / req
—
regs / thread
—
smem / block
—
block size
—
peak bw
—
peak fp32 tflops
—
peak fp16 tflops
—
reference PyTorch source (the prompt the agent sees)
import torch
import torch.nn as nn
import torch.nn.functional as F


class Model(nn.Module):
    """
    Contact map prediction head.  Takes residue-pair embeddings and
    predicts a symmetric binary contact probability matrix (residue i
    within 8 Å of residue j).  Uses a small MLP on the symmetrised pair
    representation followed by a sigmoid.  Common final stage in protein
    structure prediction pipelines.
    """

    def __init__(self, pair_dim, hidden_dim):
        super().__init__()
        self.layer_norm = nn.LayerNorm(pair_dim)
        self.mlp = nn.Sequential(
            nn.Linear(pair_dim, hidden_dim),
            nn.ReLU(),
            nn.Linear(hidden_dim, hidden_dim),
            nn.ReLU(),
            nn.Linear(hidden_dim, 1),
        )

    def forward(self, pair_repr: torch.Tensor) -> torch.Tensor:
        """
        Args:
            pair_repr: (B, N, N, pair_dim)
        Returns:
            contact_probs: (B, N, N) – symmetric contact probability
        """
        x = self.layer_norm(pair_repr)
        # Symmetrise
        x = (x + x.permute(0, 2, 1, 3)) / 2.0
        logits = self.mlp(x).squeeze(-1)  # (B, N, N)
        return torch.sigmoid(logits)


pair_dim = 64
hidden_dim = 32
seq_len = 64
batch_size = 4


def get_inputs():
    return [torch.randn(popcorn_pri.jitter_int(batch_size), popcorn_pri.jitter_int(seq_len), popcorn_pri.jitter_int(seq_len), pair_dim)]


def get_init_inputs():
    return [pair_dim, hidden_dim]
turn 0189.1s turn · 2.3s LLM · 186.8s tools · tokens in=1,588 · out=157 · total=1,745
reasoning
Let me analyze this task. I need to replace a PyTorch model with a custom CUDA kernel for contact map prediction. The model does: 1. LayerNorm on pair_repr (B, N, N, pair_dim) 2. Symmetrize: (x + x.permute(0, 2, 1, 3)) / 2.0 3. MLP: Linear(pair_dim, hidden_dim) -> ReLU -> Linear(hidden_dim, hidden_dim) -> ReLU -> Linear(hidden_dim, 1) 4. Squeeze the last dimension 5. Sigmoid So the output is (B, N, N) with sigmoid applied. Let me first get GPU specs to understand the hardware constraints.
call → get_gpu_specs
callget_gpu_specs {}
result ← get_gpu_specs (ok)
okget_gpu_specs GPU specs for NVIDIA H100 80GB HBM3: total memory (runtime): 79.1 GB spec entry: H100 GPU Architecture: Hopper GPU Memory: 80GB Memory Bandwidth: 3.35 TB/s FP64 TFLOPS: 34 FP64 Tensor Core TFLOPS: 67 FP32 TFLOPS: 67 TF32 Tensor Core TFLOPS: 989 with sparsity BFLOAT16 Tensore Core TFLOPS: 1979 with sparsity FP16 Tensor Core TFLOPS: 1979 with sparsity FP8 Tensor Core TFLOPS: 3958 with sparsity INT8 Tensor Core TOPS: 3958 with sparsity Register File Size: 64K 32-bit registers per SM Maximum number of registers per thread: 255 Maximum number of thread blocks per SM: 32 Shared memory capacity per SM: 228 KB Maximum shared memory per thread block: 227 KB
turn 1113.9s turn · 13.2s LLM · 100.7s tools · tokens in=1,836 · out=1,251 · total=3,087
reasoning
Now let me write my first kernel implementation. I need to: 1. Fuse LayerNorm + Symmetrize into one kernel 2. Fuse the 3-layer MLP + sigmoid into another kernel Let me start with a straightforward implementation:
assistantI'll start by writing a custom CUDA kernel that fuses the LayerNorm, symmetrization, MLP, and sigmoid operations. Let me begin with a clean implementation.
call → compile_kernel
callcompile_kernel { "kernel_code": "import torch\nimport torch.nn as nn\nimport torch.nn.functional as F\nfrom torch.utils.cpp_extension import load_inline\n\n# Define the CUDA kernel\ncuda_source = \"\"\"\n#include <torch/extension.h>\n#include <cuda.h>\n#include <cuda_runtime.h>\n#include <c10/cuda/CUDAException.h>\n\ntemplate <typename scalar_t>\n__global__ void fused_layernorm_symmetrize_kernel(\n const scalar_t* __restrict__ input,\n scalar_t* __restrict__ output,\n int N, int D,\n const float* __restrict__ weight,\n const float* __restrict__ bias) {\n \n // input: (B, N, N, D)\n // output: (B, N, N, D)\n \n int b = blockIdx.z;\n int i = blockIdx.y;\n int j = blockIdx.x;\n \n if (b >= batch_size || i >= N || j >= N) return;\n \n // Each thread handles one (i,j) position\n int idx_ij = b * N * N * D + i * N * D + j * D;\n int idx_ji = b * N * N * D + j * N * D + i * D;\n \n // Load and compute mean/variance for position (i,j)\n scalar_t mean = 0;\n for (int d = 0; d < D; d++) {\n mean += input[idx_ij + d];\n }\n mean /= D;\n \n scalar_t var = 0;\n for (int d = 0; d < D; d++) {\n scalar_t diff = input[idx_ij + d] - mean;\n var += diff * diff;\n }\n var /= D;\n \n scalar_t inv_std = rsqrtf(var + 1e-5f);\n \n // Load and compute mean/variance for position (j,i)\n scalar_t mean_ji = 0;\n for (int d = 0; d < D; d++) {\n mean_ji += input[idx_ji + d];\n }\n mean_ji /= D;\n \n scalar_t var_ji = 0;\n for (int d = 0; d < D; d++) {\n scalar_t diff = input[idx_ji + d] - mean_ji;\n var_ji += diff * diff;\n }\n var_ji /= D;\n \n scalar_t inv_std_ji = rsqrtf(var_ji + 1e-5f);\n \n // Compute normalized and symmetrized values\n for (int d = 0; d < D; d++) {\n scalar_t norm_ij = (input[idx_ij + d] - mean) * inv_std;\n scalar_t norm_ji = (input[idx_ji + d] - mean_ji) * inv_std_ji;\n \n // Apply weight/bias and symmetrize\n scalar_t val = (norm_ij * weight[d] + bias[d] + norm_ji * weight[d] + bias[d]) / 2.0f;\n output[idx_ij + d] = val;\n }\n}\n\nint batch_size;\n\nvoid set_batch_size(int bs) {\n batch_size = bs;\n}\n\ntorch::Tensor fused_layernorm_symmetrize_cuda(\n torch::Tensor input,\n torch::Tensor weight,\n torch::Tensor bias) {\n \n auto output = torch::empty_like(input);\n int B = input.size(0);\n int N = input.size(1);\n int D = input.size(3);\n \n batch_size = B;\n \n dim3 grid(N, N, B);\n dim3 block(1);\n \n AT_DISPATCH_FLOATING_TYPES(input.scalar_type(), \"fused_layernorm_symmetrize\", ([&] {\n fused_layernorm_symmetrize_kernel<scalar_t><<<grid, block>>>(\n input.data_ptr<scalar_t>(),\n output.data_ptr<scalar_t>(),\n N, D,\n weight.data_ptr<float>(),\n bias.data_ptr<float>());\n }));\n \n return output;\n}\n\"\"\"\n\ncpp_source = \"\"\"\nint batch_size;\nvoid set_batch_size(int bs);\ntorch::Tensor fused_layernorm_symmetrize_cuda(torch::Tensor input, torch::Tensor weight, torch::Tensor bias);\n\"\"\"\n\n# Load the extension\nfused_ops = load_inline(\n name='fused_ops',\n cpp_sources=cpp_source,\n cuda_sources=cuda_source,\n functions=['fused_layernorm_symmetrize_cuda', 'set_batch_size'],\n verbose=True\n)\n\n\nclass ModelNew(nn.Module):\n def __init__(self, pair_dim, hidden_dim):\n super().__init__()\n self.layer_norm = nn.LayerNorm(pair_dim)\n self.mlp = nn.Sequential(\n nn.Linear(pair_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, 1),\n )\n self.pair_dim = pair_dim\n self.hidden_dim = hidden_dim\n\n def forward(self, pair_repr: torch.Tensor) -> torch.Tensor:\n # Use fused kernel for LayerNorm + symmetrize\n x = fused_ops.fused_layernorm_symmetrize_cuda(\n pair_repr, \n self.layer_norm.weight, \n self.layer_norm.bias\n )\n logits = self.mlp(x).squeeze(-1)\n return torch.sigmoid(logits)" }
result ← compile_kernel (fail)
failcompile_kernel compile_kernel FAILED: RuntimeError. RuntimeError: Error building extension 'fused_ops'
turn 2169.4s turn · 10.6s LLM · 158.8s tools · tokens in=3,057 · out=1,191 · total=4,248
reasoning
The extension build failed. Let me simplify and fix the CUDA code. The issue is likely with how I'm defining the batch_size variable and the overall structure. Let me rewrite with a cleaner approach.
assistantLet me fix the CUDA code with a cleaner approach:
call → compile_kernel
callcompile_kernel { "kernel_code": "import torch\nimport torch.nn as nn\nimport torch.nn.functional as F\nfrom torch.utils.cpp_extension import load_inline\n\n# Define the CUDA kernel\ncuda_source = \"\"\"\n#include <torch/extension.h>\n#include <cuda.h>\n#include <cuda_runtime.h>\n#include <c10/cuda/CUDAException.h>\n\ntemplate <typename scalar_t>\n__global__ void fused_layernorm_symmetrize_kernel(\n const scalar_t* __restrict__ input,\n scalar_t* __restrict__ output,\n int B, int N, int D,\n const float* __restrict__ weight,\n const float* __restrict__ bias) {\n \n // input: (B, N, N, D)\n // output: (B, N, N, D)\n \n int b = blockIdx.z;\n int i = blockIdx.y;\n int j = blockIdx.x;\n \n if (b >= B || i >= N || j >= N) return;\n \n // Each thread handles one (i,j) position\n int idx_ij = b * N * N * D + i * N * D + j * D;\n int idx_ji = b * N * N * D + j * N * D + i * D;\n \n // Load and compute mean/variance for position (i,j)\n scalar_t mean = 0;\n for (int d = 0; d < D; d++) {\n mean += input[idx_ij + d];\n }\n mean /= D;\n \n scalar_t var = 0;\n for (int d = 0; d < D; d++) {\n scalar_t diff = input[idx_ij + d] - mean;\n var += diff * diff;\n }\n var /= D;\n \n scalar_t inv_std = rsqrtf(var + 1e-5f);\n \n // Load and compute mean/variance for position (j,i)\n scalar_t mean_ji = 0;\n for (int d = 0; d < D; d++) {\n mean_ji += input[idx_ji + d];\n }\n mean_ji /= D;\n \n scalar_t var_ji = 0;\n for (int d = 0; d < D; d++) {\n scalar_t diff = input[idx_ji + d] - mean_ji;\n var_ji += diff * diff;\n }\n var_ji /= D;\n \n scalar_t inv_std_ji = rsqrtf(var_ji + 1e-5f);\n \n // Compute normalized and symmetrized values\n for (int d = 0; d < D; d++) {\n scalar_t norm_ij = (input[idx_ij + d] - mean) * inv_std;\n scalar_t norm_ji = (input[idx_ji + d] - mean_ji) * inv_std_ji;\n \n // Apply weight/bias and symmetrize\n scalar_t val_ij = norm_ij * weight[d] + bias[d];\n scalar_t val_ji = norm_ji * weight[d] + bias[d];\n output[idx_ij + d] = (val_ij + val_ji) * 0.5f;\n }\n}\n\ntorch::Tensor fused_layernorm_symmetrize_cuda(\n torch::Tensor input,\n torch::Tensor weight,\n torch::Tensor bias) {\n \n auto output = torch::empty_like(input);\n int B = input.size(0);\n int N = input.size(1);\n int D = input.size(3);\n \n dim3 grid(N, N, B);\n dim3 block(1);\n \n AT_DISPATCH_FLOATING_TYPES(input.scalar_type(), \"fused_layernorm_symmetrize\", ([&] {\n fused_layernorm_symmetrize_kernel<scalar_t><<<grid, block>>>(\n input.data_ptr<scalar_t>(),\n output.data_ptr<scalar_t>(),\n B, N, D,\n weight.data_ptr<float>(),\n bias.data_ptr<float>());\n }));\n \n return output;\n}\n\"\"\"\n\ncpp_source = \"\"\"\ntorch::Tensor fused_layernorm_symmetrize_cuda(torch::Tensor input, torch::Tensor weight, torch::Tensor bias);\n\"\"\"\n\n# Load the extension\nfused_ops = load_inline(\n name='fused_ops',\n cpp_sources=cpp_source,\n cuda_sources=cuda_source,\n functions=['fused_layernorm_symmetrize_cuda'],\n verbose=True\n)\n\n\nclass ModelNew(nn.Module):\n def __init__(self, pair_dim, hidden_dim):\n super().__init__()\n self.layer_norm = nn.LayerNorm(pair_dim)\n self.mlp = nn.Sequential(\n nn.Linear(pair_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, 1),\n )\n self.pair_dim = pair_dim\n self.hidden_dim = hidden_dim\n\n def forward(self, pair_repr: torch.Tensor) -> torch.Tensor:\n # Use fused kernel for LayerNorm + symmetrize\n x = fused_ops.fused_layernorm_symmetrize_cuda(\n pair_repr, \n self.layer_norm.weight, \n self.layer_norm.bias\n )\n logits = self.mlp(x).squeeze(-1)\n return torch.sigmoid(logits)" }
result ← compile_kernel (ok)
okcompile_kernel compile_kernel PASSED: kernel compiled without errors.
turn 3240.9s turn · 9.9s LLM · 231.0s tools · tokens in=4,221 · out=1,152 · total=5,373
reasoning
Good, it compiled. Now let's check correctness.
call → run_correctness
callrun_correctness { "kernel_code": "import torch\nimport torch.nn as nn\nimport torch.nn.functional as F\nfrom torch.utils.cpp_extension import load_inline\n\n# Define the CUDA kernel\ncuda_source = \"\"\"\n#include <torch/extension.h>\n#include <cuda.h>\n#include <cuda_runtime.h>\n#include <c10/cuda/CUDAException.h>\n\ntemplate <typename scalar_t>\n__global__ void fused_layernorm_symmetrize_kernel(\n const scalar_t* __restrict__ input,\n scalar_t* __restrict__ output,\n int B, int N, int D,\n const float* __restrict__ weight,\n const float* __restrict__ bias) {\n \n // input: (B, N, N, D)\n // output: (B, N, N, D)\n \n int b = blockIdx.z;\n int i = blockIdx.y;\n int j = blockIdx.x;\n \n if (b >= B || i >= N || j >= N) return;\n \n // Each thread handles one (i,j) position\n int idx_ij = b * N * N * D + i * N * D + j * D;\n int idx_ji = b * N * N * D + j * N * D + i * D;\n \n // Load and compute mean/variance for position (i,j)\n scalar_t mean = 0;\n for (int d = 0; d < D; d++) {\n mean += input[idx_ij + d];\n }\n mean /= D;\n \n scalar_t var = 0;\n for (int d = 0; d < D; d++) {\n scalar_t diff = input[idx_ij + d] - mean;\n var += diff * diff;\n }\n var /= D;\n \n scalar_t inv_std = rsqrtf(var + 1e-5f);\n \n // Load and compute mean/variance for position (j,i)\n scalar_t mean_ji = 0;\n for (int d = 0; d < D; d++) {\n mean_ji += input[idx_ji + d];\n }\n mean_ji /= D;\n \n scalar_t var_ji = 0;\n for (int d = 0; d < D; d++) {\n scalar_t diff = input[idx_ji + d] - mean_ji;\n var_ji += diff * diff;\n }\n var_ji /= D;\n \n scalar_t inv_std_ji = rsqrtf(var_ji + 1e-5f);\n \n // Compute normalized and symmetrized values\n for (int d = 0; d < D; d++) {\n scalar_t norm_ij = (input[idx_ij + d] - mean) * inv_std;\n scalar_t norm_ji = (input[idx_ji + d] - mean_ji) * inv_std_ji;\n \n // Apply weight/bias and symmetrize\n scalar_t val_ij = norm_ij * weight[d] + bias[d];\n scalar_t val_ji = norm_ji * weight[d] + bias[d];\n output[idx_ij + d] = (val_ij + val_ji) * 0.5f;\n }\n}\n\ntorch::Tensor fused_layernorm_symmetrize_cuda(\n torch::Tensor input,\n torch::Tensor weight,\n torch::Tensor bias) {\n \n auto output = torch::empty_like(input);\n int B = input.size(0);\n int N = input.size(1);\n int D = input.size(3);\n \n dim3 grid(N, N, B);\n dim3 block(1);\n \n AT_DISPATCH_FLOATING_TYPES(input.scalar_type(), \"fused_layernorm_symmetrize\", ([&] {\n fused_layernorm_symmetrize_kernel<scalar_t><<<grid, block>>>(\n input.data_ptr<scalar_t>(),\n output.data_ptr<scalar_t>(),\n B, N, D,\n weight.data_ptr<float>(),\n bias.data_ptr<float>());\n }));\n \n return output;\n}\n\"\"\"\n\ncpp_source = \"\"\"\ntorch::Tensor fused_layernorm_symmetrize_cuda(torch::Tensor input, torch::Tensor weight, torch::Tensor bias);\n\"\"\"\n\n# Load the extension\nfused_ops = load_inline(\n name='fused_ops',\n cpp_sources=cpp_source,\n cuda_sources=cuda_source,\n functions=['fused_layernorm_symmetrize_cuda'],\n verbose=True\n)\n\n\nclass ModelNew(nn.Module):\n def __init__(self, pair_dim, hidden_dim):\n super().__init__()\n self.layer_norm = nn.LayerNorm(pair_dim)\n self.mlp = nn.Sequential(\n nn.Linear(pair_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, 1),\n )\n self.pair_dim = pair_dim\n self.hidden_dim = hidden_dim\n\n def forward(self, pair_repr: torch.Tensor) -> torch.Tensor:\n # Use fused kernel for LayerNorm + symmetrize\n x = fused_ops.fused_layernorm_symmetrize_cuda(\n pair_repr, \n self.layer_norm.weight, \n self.layer_norm.bias\n )\n logits = self.mlp(x).squeeze(-1)\n return torch.sigmoid(logits)" }
result ← run_correctness (ok)
okrun_correctness run_correctness PASSED: (5 / 5) trials all matched the reference. Numerical precision: max_abs_err=1.19e-07 mean_abs_err=2.41e-09 max_rel_err=2.32e-07
turn 4 · 1 tool call(s) · 168.7s turn · 17.0s LLM · tokens in=5,419 · out=1,809 · total=7,228 · 7,550 chars (click to expand)
turn 4168.7s turn · 17.0s LLM · 151.7s tools · tokens in=5,419 · out=1,809 · total=7,228
reasoning
The correctness test passed! The fused kernel is working correctly. Now let me optimize it further. Currently, each thread is handling one (i,j) position which is very inefficient. I should: 1. Use multiple threads per (i,j) to process the D dimension in parallel 2. Use shared memory for reduction 3. Maybe fuse the MLP as well Let me create an optimized version with better thread utilization.
assistantThe correctness test passed. Now let me optimize the kernel by using better thread utilization and potentially fusing more operations. Let me create a more optimized version:
call → compile_kernel
callcompile_kernel { "kernel_code": "import torch\nimport torch.nn as nn\nimport torch.nn.functional as F\nfrom torch.utils.cpp_extension import load_inline\n\n# Define the CUDA kernel with optimized thread usage\ncuda_source = \"\"\"\n#include <torch/extension.h>\n#include <cuda.h>\n#include <cuda_runtime.h>\n#include <c10/cuda/CUDAException.h>\n\ntemplate <typename scalar_t>\n__global__ void fused_layernorm_symmetrize_kernel(\n const scalar_t* __restrict__ input,\n scalar_t* __restrict__ output,\n int B, int N, int D,\n const float* __restrict__ weight,\n const float* __restrict__ bias) {\n \n // input: (B, N, N, D)\n // output: (B, N, N, D)\n // Each block handles one (b, i, j) tuple\n // Threads in block cooperatively process the D dimension\n \n int b = blockIdx.z;\n int i = blockIdx.y;\n int j = blockIdx.x;\n int tid = threadIdx.x;\n \n int idx_ij = b * N * N * D + i * N * D + j * D;\n int idx_ji = b * N * N * D + j * N * D + i * D;\n \n // Load data for position (i,j) - each thread loads D elements\n extern __shared__ char smem[];\n scalar_t* s_ij = reinterpret_cast<scalar_t*>(smem);\n scalar_t* s_ji = s_ij + D;\n scalar_t* s_mean = s_ji + D;\n scalar_t* s_var = s_mean + 2; // 2 elements for mean and var\n \n // Cooperatively load data\n for (int d = tid; d < D; d += blockDim.x) {\n s_ij[d] = input[idx_ij + d];\n s_ji[d] = input[idx_ji + d];\n }\n __syncthreads();\n \n // Compute mean for (i,j) - reduction\n scalar_t local_mean = 0;\n for (int d = tid; d < D; d += blockDim.x) {\n local_mean += s_ij[d];\n }\n \n // Block-level reduction for mean\n for (int offset = blockDim.x / 2; offset > 0; offset >>= 1) {\n local_mean += __shfl_down_sync(0xffffffff, local_mean, offset);\n }\n \n if (tid == 0) {\n s_mean[0] = local_mean / D;\n }\n __syncthreads();\n \n scalar_t mean_ij = s_mean[0];\n \n // Compute variance for (i,j)\n scalar_t local_var = 0;\n for (int d = tid; d < D; d += blockDim.x) {\n scalar_t diff = s_ij[d] - mean_ij;\n local_var += diff * diff;\n }\n \n // Block-level reduction for variance\n for (int offset = blockDim.x / 2; offset > 0; offset >>= 1) {\n local_var += __shfl_down_sync(0xffffffff, local_var, offset);\n }\n \n if (tid == 0) {\n s_var[0] = local_var / D;\n }\n __syncthreads();\n \n scalar_t var_ij = s_var[0];\n scalar_t inv_std_ij = rsqrtf(var_ij + 1e-5f);\n \n // Compute mean for (j,i)\n scalar_t local_mean_ji = 0;\n for (int d = tid; d < D; d += blockDim.x) {\n local_mean_ji += s_ji[d];\n }\n \n for (int offset = blockDim.x / 2; offset > 0; offset >>= 1) {\n local_mean_ji += __shfl_down_sync(0xffffffff, local_mean_ji, offset);\n }\n \n if (tid == 0) {\n s_mean[1] = local_mean_ji / D;\n }\n __syncthreads();\n \n scalar_t mean_ji = s_mean[1];\n \n // Compute variance for (j,i)\n scalar_t local_var_ji = 0;\n for (int d = tid; d < D; d += blockDim.x) {\n scalar_t diff = s_ji[d] - mean_ji;\n local_var_ji += diff * diff;\n }\n \n for (int offset = blockDim.x / 2; offset > 0; offset >>= 1) {\n local_var_ji += __shfl_down_sync(0xffffffff, local_var_ji, offset);\n }\n \n if (tid == 0) {\n s_var[1] = local_var_ji / D;\n }\n __syncthreads();\n \n scalar_t var_ji = s_var[1];\n scalar_t inv_std_ji = rsqrtf(var_ji + 1e-5f);\n \n // Compute normalized and symmetrized values\n for (int d = tid; d < D; d += blockDim.x) {\n scalar_t norm_ij = (s_ij[d] - mean_ij) * inv_std_ij;\n scalar_t norm_ji = (s_ji[d] - mean_ji) * inv_std_ji;\n \n scalar_t val_ij = norm_ij * weight[d] + bias[d];\n scalar_t val_ji = norm_ji * weight[d] + bias[d];\n output[idx_ij + d] = (val_ij + val_ji) * 0.5f;\n }\n}\n\ntorch::Tensor fused_layernorm_symmetrize_cuda(\n torch::Tensor input,\n torch::Tensor weight,\n torch::Tensor bias) {\n \n auto output = torch::empty_like(input);\n int B = input.size(0);\n int N = input.size(1);\n int D = input.size(3);\n \n dim3 grid(N, N, B);\n int threads = min(256, D);\n int shared_mem = D * 2 * sizeof(float) + 4 * sizeof(float); // s_ij, s_ji, mean, var\n \n AT_DISPATCH_FLOATING_TYPES(input.scalar_type(), \"fused_layernorm_symmetrize\", ([&] {\n fused_layernorm_symmetrize_kernel<scalar_t><<<grid, threads, shared_mem>>>(\n input.data_ptr<scalar_t>(),\n output.data_ptr<scalar_t>(),\n B, N, D,\n weight.data_ptr<float>(),\n bias.data_ptr<float>());\n }));\n \n return output;\n}\n\"\"\"\n\ncpp_source = \"\"\"\ntorch::Tensor fused_layernorm_symmetrize_cuda(torch::Tensor input, torch::Tensor weight, torch::Tensor bias);\n\"\"\"\n\n# Load the extension\nfused_ops = load_inline(\n name='fused_ops',\n cpp_sources=cpp_source,\n cuda_sources=cuda_source,\n functions=['fused_layernorm_symmetrize_cuda'],\n verbose=True\n)\n\n\nclass ModelNew(nn.Module):\n def __init__(self, pair_dim, hidden_dim):\n super().__init__()\n self.layer_norm = nn.LayerNorm(pair_dim)\n self.mlp = nn.Sequential(\n nn.Linear(pair_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, 1),\n )\n self.pair_dim = pair_dim\n self.hidden_dim = hidden_dim\n\n def forward(self, pair_repr: torch.Tensor) -> torch.Tensor:\n # Use fused kernel for LayerNorm + symmetrize\n x = fused_ops.fused_layernorm_symmetrize_cuda(\n pair_repr, \n self.layer_norm.weight, \n self.layer_norm.bias\n )\n logits = self.mlp(x).squeeze(-1)\n return torch.sigmoid(logits)" }
result ← compile_kernel (ok)
okcompile_kernel compile_kernel PASSED: kernel compiled without errors.
turn 5 · 1 tool call(s) · 172.8s turn · 14.5s LLM · tokens in=7,156 · out=1,713 · total=8,869 · 7,309 chars (click to expand)
turn 5172.8s turn · 14.5s LLM · 158.4s tools · tokens in=7,156 · out=1,713 · total=8,869
reasoning
Good, the compilation passed. Let me run the correctness check to ensure the output is correct.
call → run_correctness
callrun_correctness { "kernel_code": "import torch\nimport torch.nn as nn\nimport torch.nn.functional as F\nfrom torch.utils.cpp_extension import load_inline\n\n# Define the CUDA kernel with optimized thread usage\ncuda_source = \"\"\"\n#include <torch/extension.h>\n#include <cuda.h>\n#include <cuda_runtime.h>\n#include <c10/cuda/CUDAException.h>\n\ntemplate <typename scalar_t>\n__global__ void fused_layernorm_symmetrize_kernel(\n const scalar_t* __restrict__ input,\n scalar_t* __restrict__ output,\n int B, int N, int D,\n const float* __restrict__ weight,\n const float* __restrict__ bias) {\n \n // input: (B, N, N, D)\n // output: (B, N, N, D)\n // Each block handles one (b, i, j) tuple\n // Threads in block cooperatively process the D dimension\n \n int b = blockIdx.z;\n int i = blockIdx.y;\n int j = blockIdx.x;\n int tid = threadIdx.x;\n \n int idx_ij = b * N * N * D + i * N * D + j * D;\n int idx_ji = b * N * N * D + j * N * D + i * D;\n \n // Load data for position (i,j) - each thread loads D elements\n extern __shared__ char smem[];\n scalar_t* s_ij = reinterpret_cast<scalar_t*>(smem);\n scalar_t* s_ji = s_ij + D;\n scalar_t* s_mean = s_ji + D;\n scalar_t* s_var = s_mean + 2; // 2 elements for mean and var\n \n // Cooperatively load data\n for (int d = tid; d < D; d += blockDim.x) {\n s_ij[d] = input[idx_ij + d];\n s_ji[d] = input[idx_ji + d];\n }\n __syncthreads();\n \n // Compute mean for (i,j) - reduction\n scalar_t local_mean = 0;\n for (int d = tid; d < D; d += blockDim.x) {\n local_mean += s_ij[d];\n }\n \n // Block-level reduction for mean\n for (int offset = blockDim.x / 2; offset > 0; offset >>= 1) {\n local_mean += __shfl_down_sync(0xffffffff, local_mean, offset);\n }\n \n if (tid == 0) {\n s_mean[0] = local_mean / D;\n }\n __syncthreads();\n \n scalar_t mean_ij = s_mean[0];\n \n // Compute variance for (i,j)\n scalar_t local_var = 0;\n for (int d = tid; d < D; d += blockDim.x) {\n scalar_t diff = s_ij[d] - mean_ij;\n local_var += diff * diff;\n }\n \n // Block-level reduction for variance\n for (int offset = blockDim.x / 2; offset > 0; offset >>= 1) {\n local_var += __shfl_down_sync(0xffffffff, local_var, offset);\n }\n \n if (tid == 0) {\n s_var[0] = local_var / D;\n }\n __syncthreads();\n \n scalar_t var_ij = s_var[0];\n scalar_t inv_std_ij = rsqrtf(var_ij + 1e-5f);\n \n // Compute mean for (j,i)\n scalar_t local_mean_ji = 0;\n for (int d = tid; d < D; d += blockDim.x) {\n local_mean_ji += s_ji[d];\n }\n \n for (int offset = blockDim.x / 2; offset > 0; offset >>= 1) {\n local_mean_ji += __shfl_down_sync(0xffffffff, local_mean_ji, offset);\n }\n \n if (tid == 0) {\n s_mean[1] = local_mean_ji / D;\n }\n __syncthreads();\n \n scalar_t mean_ji = s_mean[1];\n \n // Compute variance for (j,i)\n scalar_t local_var_ji = 0;\n for (int d = tid; d < D; d += blockDim.x) {\n scalar_t diff = s_ji[d] - mean_ji;\n local_var_ji += diff * diff;\n }\n \n for (int offset = blockDim.x / 2; offset > 0; offset >>= 1) {\n local_var_ji += __shfl_down_sync(0xffffffff, local_var_ji, offset);\n }\n \n if (tid == 0) {\n s_var[1] = local_var_ji / D;\n }\n __syncthreads();\n \n scalar_t var_ji = s_var[1];\n scalar_t inv_std_ji = rsqrtf(var_ji + 1e-5f);\n \n // Compute normalized and symmetrized values\n for (int d = tid; d < D; d += blockDim.x) {\n scalar_t norm_ij = (s_ij[d] - mean_ij) * inv_std_ij;\n scalar_t norm_ji = (s_ji[d] - mean_ji) * inv_std_ji;\n \n scalar_t val_ij = norm_ij * weight[d] + bias[d];\n scalar_t val_ji = norm_ji * weight[d] + bias[d];\n output[idx_ij + d] = (val_ij + val_ji) * 0.5f;\n }\n}\n\ntorch::Tensor fused_layernorm_symmetrize_cuda(\n torch::Tensor input,\n torch::Tensor weight,\n torch::Tensor bias) {\n \n auto output = torch::empty_like(input);\n int B = input.size(0);\n int N = input.size(1);\n int D = input.size(3);\n \n dim3 grid(N, N, B);\n int threads = min(256, D);\n int shared_mem = D * 2 * sizeof(float) + 4 * sizeof(float); // s_ij, s_ji, mean, var\n \n AT_DISPATCH_FLOATING_TYPES(input.scalar_type(), \"fused_layernorm_symmetrize\", ([&] {\n fused_layernorm_symmetrize_kernel<scalar_t><<<grid, threads, shared_mem>>>(\n input.data_ptr<scalar_t>(),\n output.data_ptr<scalar_t>(),\n B, N, D,\n weight.data_ptr<float>(),\n bias.data_ptr<float>());\n }));\n \n return output;\n}\n\"\"\"\n\ncpp_source = \"\"\"\ntorch::Tensor fused_layernorm_symmetrize_cuda(torch::Tensor input, torch::Tensor weight, torch::Tensor bias);\n\"\"\"\n\n# Load the extension\nfused_ops = load_inline(\n name='fused_ops',\n cpp_sources=cpp_source,\n cuda_sources=cuda_source,\n functions=['fused_layernorm_symmetrize_cuda'],\n verbose=True\n)\n\n\nclass ModelNew(nn.Module):\n def __init__(self, pair_dim, hidden_dim):\n super().__init__()\n self.layer_norm = nn.LayerNorm(pair_dim)\n self.mlp = nn.Sequential(\n nn.Linear(pair_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, 1),\n )\n self.pair_dim = pair_dim\n self.hidden_dim = hidden_dim\n\n def forward(self, pair_repr: torch.Tensor) -> torch.Tensor:\n # Use fused kernel for LayerNorm + symmetrize\n x = fused_ops.fused_layernorm_symmetrize_cuda(\n pair_repr, \n self.layer_norm.weight, \n self.layer_norm.bias\n )\n logits = self.mlp(x).squeeze(-1)\n return torch.sigmoid(logits)" }
result ← run_correctness (fail)
failrun_correctness run_correctness FAILED: (0 / 5) trials did not all match. correctness_issue: Output mismatch max_difference: ['0.014363', '0.017218', '0.029618', '0.021994', '0.019088'] avg_difference: ['0.001765', '0.001762', '0.001783', '0.001783', '0.001818']
turn 6 · 1 tool call(s) · 134.3s turn · 15.5s LLM · tokens in=8,956 · out=1,842 · total=10,798 · 7,535 chars (click to expand)
turn 6134.3s turn · 15.5s LLM · 118.8s tools · tokens in=8,956 · out=1,842 · total=10,798
reasoning
The optimized version has small numerical differences. The issue might be with the shared memory allocation - I'm using sizeof(float) but the scalar type might be different. Let me fix this.
assistantThe optimized version has numerical differences. Let me fix the shared memory allocation issue:
call → compile_kernel
callcompile_kernel { "kernel_code": "import torch\nimport torch.nn as nn\nimport torch.nn.functional as F\nfrom torch.utils.cpp_extension import load_inline\n\n# Define the CUDA kernel with optimized thread usage\ncuda_source = \"\"\"\n#include <torch/extension.h>\n#include <cuda.h>\n#include <cuda_runtime.h>\n#include <c10/cuda/CUDAException.h>\n\ntemplate <typename scalar_t>\n__global__ void fused_layernorm_symmetrize_kernel(\n const scalar_t* __restrict__ input,\n scalar_t* __restrict__ output,\n int B, int N, int D,\n const float* __restrict__ weight,\n const float* __restrict__ bias) {\n \n // input: (B, N, N, D)\n // output: (B, N, N, D)\n // Each block handles one (b, i, j) tuple\n // Threads in block cooperatively process the D dimension\n \n int b = blockIdx.z;\n int i = blockIdx.y;\n int j = blockIdx.x;\n int tid = threadIdx.x;\n \n int idx_ij = b * N * N * D + i * N * D + j * D;\n int idx_ji = b * N * N * D + j * N * D + i * D;\n \n // Load data for position (i,j) - each thread loads D elements\n extern __shared__ char smem[];\n scalar_t* s_ij = reinterpret_cast<scalar_t*>(smem);\n scalar_t* s_ji = s_ij + D;\n scalar_t* s_mean_ij = s_ji + D;\n scalar_t* s_var_ij = s_mean_ij + 1;\n scalar_t* s_mean_ji = s_var_ij + 1;\n scalar_t* s_var_ji = s_mean_ji + 1;\n \n // Cooperatively load data\n for (int d = tid; d < D; d += blockDim.x) {\n s_ij[d] = input[idx_ij + d];\n s_ji[d] = input[idx_ji + d];\n }\n __syncthreads();\n \n // Compute mean for (i,j) - reduction\n scalar_t local_mean = 0;\n for (int d = tid; d < D; d += blockDim.x) {\n local_mean += s_ij[d];\n }\n \n // Block-level reduction for mean using warp shuffle\n for (int offset = blockDim.x / 2; offset > 0; offset >>= 1) {\n local_mean += __shfl_down_sync(0xffffffff, local_mean, offset);\n }\n \n if (tid == 0) {\n s_mean_ij[0] = local_mean / (scalar_t)D;\n }\n __syncthreads();\n \n scalar_t mean_ij = s_mean_ij[0];\n \n // Compute variance for (i,j)\n scalar_t local_var = 0;\n for (int d = tid; d < D; d += blockDim.x) {\n scalar_t diff = s_ij[d] - mean_ij;\n local_var += diff * diff;\n }\n \n // Block-level reduction for variance\n for (int offset = blockDim.x / 2; offset > 0; offset >>= 1) {\n local_var += __shfl_down_sync(0xffffffff, local_var, offset);\n }\n \n if (tid == 0) {\n s_var_ij[0] = local_var / (scalar_t)D;\n }\n __syncthreads();\n \n scalar_t var_ij = s_var_ij[0];\n scalar_t inv_std_ij = rsqrtf(var_ij + (scalar_t)1e-5);\n \n // Compute mean for (j,i)\n scalar_t local_mean_ji = 0;\n for (int d = tid; d < D; d += blockDim.x) {\n local_mean_ji += s_ji[d];\n }\n \n for (int offset = blockDim.x / 2; offset > 0; offset >>= 1) {\n local_mean_ji += __shfl_down_sync(0xffffffff, local_mean_ji, offset);\n }\n \n if (tid == 0) {\n s_mean_ji[0] = local_mean_ji / (scalar_t)D;\n }\n __syncthreads();\n \n scalar_t mean_ji = s_mean_ji[0];\n \n // Compute variance for (j,i)\n scalar_t local_var_ji = 0;\n for (int d = tid; d < D; d += blockDim.x) {\n scalar_t diff = s_ji[d] - mean_ji;\n local_var_ji += diff * diff;\n }\n \n for (int offset = blockDim.x / 2; offset > 0; offset >>= 1) {\n local_var_ji += __shfl_down_sync(0xffffffff, local_var_ji, offset);\n }\n \n if (tid == 0) {\n s_var_ji[0] = local_var_ji / (scalar_t)D;\n }\n __syncthreads();\n \n scalar_t var_ji = s_var_ji[0];\n scalar_t inv_std_ji = rsqrtf(var_ji + (scalar_t)1e-5);\n \n // Compute normalized and symmetrized values\n for (int d = tid; d < D; d += blockDim.x) {\n scalar_t norm_ij = (s_ij[d] - mean_ij) * inv_std_ij;\n scalar_t norm_ji = (s_ji[d] - mean_ji) * inv_std_ji;\n \n scalar_t val_ij = norm_ij * (scalar_t)weight[d] + (scalar_t)bias[d];\n scalar_t val_ji = norm_ji * (scalar_t)weight[d] + (scalar_t)bias[d];\n output[idx_ij + d] = (val_ij + val_ji) * (scalar_t)0.5;\n }\n}\n\ntorch::Tensor fused_layernorm_symmetrize_cuda(\n torch::Tensor input,\n torch::Tensor weight,\n torch::Tensor bias) {\n \n auto output = torch::empty_like(input);\n int B = input.size(0);\n int N = input.size(1);\n int D = input.size(3);\n \n dim3 grid(N, N, B);\n int threads = min(256, D);\n // Shared memory: 2*D scalar_t for s_ij and s_ji, 4 scalar_t for means and vars\n size_t shared_mem = (D * 2 + 4) * input.element_size();\n \n AT_DISPATCH_FLOATING_TYPES(input.scalar_type(), \"fused_layernorm_symmetrize\", ([&] {\n fused_layernorm_symmetrize_kernel<scalar_t><<<grid, threads, shared_mem>>>(\n input.data_ptr<scalar_t>(),\n output.data_ptr<scalar_t>(),\n B, N, D,\n weight.data_ptr<float>(),\n bias.data_ptr<float>());\n }));\n \n return output;\n}\n\"\"\"\n\ncpp_source = \"\"\"\ntorch::Tensor fused_layernorm_symmetrize_cuda(torch::Tensor input, torch::Tensor weight, torch::Tensor bias);\n\"\"\"\n\n# Load the extension\nfused_ops = load_inline(\n name='fused_ops',\n cpp_sources=cpp_source,\n cuda_sources=cuda_source,\n functions=['fused_layernorm_symmetrize_cuda'],\n verbose=True\n)\n\n\nclass ModelNew(nn.Module):\n def __init__(self, pair_dim, hidden_dim):\n super().__init__()\n self.layer_norm = nn.LayerNorm(pair_dim)\n self.mlp = nn.Sequential(\n nn.Linear(pair_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, 1),\n )\n self.pair_dim = pair_dim\n self.hidden_dim = hidden_dim\n\n def forward(self, pair_repr: torch.Tensor) -> torch.Tensor:\n # Use fused kernel for LayerNorm + symmetrize\n x = fused_ops.fused_layernorm_symmetrize_cuda(\n pair_repr, \n self.layer_norm.weight, \n self.layer_norm.bias\n )\n logits = self.mlp(x).squeeze(-1)\n return torch.sigmoid(logits)" }
result ← compile_kernel (ok)
okcompile_kernel compile_kernel PASSED: kernel compiled without errors.
turn 7 · 1 tool call(s) · 1823.1s turn · 14.9s LLM · tokens in=10,774 · out=1,800 · total=12,574 · 7,532 chars (click to expand)
turn 71823.1s turn · 14.9s LLM · 1808.2s tools · tokens in=10,774 · out=1,800 · total=12,574
reasoning
The kernel compiled. Let me test for correctness.
call → run_correctness
callrun_correctness { "kernel_code": "import torch\nimport torch.nn as nn\nimport torch.nn.functional as F\nfrom torch.utils.cpp_extension import load_inline\n\n# Define the CUDA kernel with optimized thread usage\ncuda_source = \"\"\"\n#include <torch/extension.h>\n#include <cuda.h>\n#include <cuda_runtime.h>\n#include <c10/cuda/CUDAException.h>\n\ntemplate <typename scalar_t>\n__global__ void fused_layernorm_symmetrize_kernel(\n const scalar_t* __restrict__ input,\n scalar_t* __restrict__ output,\n int B, int N, int D,\n const float* __restrict__ weight,\n const float* __restrict__ bias) {\n \n // input: (B, N, N, D)\n // output: (B, N, N, D)\n // Each block handles one (b, i, j) tuple\n // Threads in block cooperatively process the D dimension\n \n int b = blockIdx.z;\n int i = blockIdx.y;\n int j = blockIdx.x;\n int tid = threadIdx.x;\n \n int idx_ij = b * N * N * D + i * N * D + j * D;\n int idx_ji = b * N * N * D + j * N * D + i * D;\n \n // Load data for position (i,j) - each thread loads D elements\n extern __shared__ char smem[];\n scalar_t* s_ij = reinterpret_cast<scalar_t*>(smem);\n scalar_t* s_ji = s_ij + D;\n scalar_t* s_mean_ij = s_ji + D;\n scalar_t* s_var_ij = s_mean_ij + 1;\n scalar_t* s_mean_ji = s_var_ij + 1;\n scalar_t* s_var_ji = s_mean_ji + 1;\n \n // Cooperatively load data\n for (int d = tid; d < D; d += blockDim.x) {\n s_ij[d] = input[idx_ij + d];\n s_ji[d] = input[idx_ji + d];\n }\n __syncthreads();\n \n // Compute mean for (i,j) - reduction\n scalar_t local_mean = 0;\n for (int d = tid; d < D; d += blockDim.x) {\n local_mean += s_ij[d];\n }\n \n // Block-level reduction for mean using warp shuffle\n for (int offset = blockDim.x / 2; offset > 0; offset >>= 1) {\n local_mean += __shfl_down_sync(0xffffffff, local_mean, offset);\n }\n \n if (tid == 0) {\n s_mean_ij[0] = local_mean / (scalar_t)D;\n }\n __syncthreads();\n \n scalar_t mean_ij = s_mean_ij[0];\n \n // Compute variance for (i,j)\n scalar_t local_var = 0;\n for (int d = tid; d < D; d += blockDim.x) {\n scalar_t diff = s_ij[d] - mean_ij;\n local_var += diff * diff;\n }\n \n // Block-level reduction for variance\n for (int offset = blockDim.x / 2; offset > 0; offset >>= 1) {\n local_var += __shfl_down_sync(0xffffffff, local_var, offset);\n }\n \n if (tid == 0) {\n s_var_ij[0] = local_var / (scalar_t)D;\n }\n __syncthreads();\n \n scalar_t var_ij = s_var_ij[0];\n scalar_t inv_std_ij = rsqrtf(var_ij + (scalar_t)1e-5);\n \n // Compute mean for (j,i)\n scalar_t local_mean_ji = 0;\n for (int d = tid; d < D; d += blockDim.x) {\n local_mean_ji += s_ji[d];\n }\n \n for (int offset = blockDim.x / 2; offset > 0; offset >>= 1) {\n local_mean_ji += __shfl_down_sync(0xffffffff, local_mean_ji, offset);\n }\n \n if (tid == 0) {\n s_mean_ji[0] = local_mean_ji / (scalar_t)D;\n }\n __syncthreads();\n \n scalar_t mean_ji = s_mean_ji[0];\n \n // Compute variance for (j,i)\n scalar_t local_var_ji = 0;\n for (int d = tid; d < D; d += blockDim.x) {\n scalar_t diff = s_ji[d] - mean_ji;\n local_var_ji += diff * diff;\n }\n \n for (int offset = blockDim.x / 2; offset > 0; offset >>= 1) {\n local_var_ji += __shfl_down_sync(0xffffffff, local_var_ji, offset);\n }\n \n if (tid == 0) {\n s_var_ji[0] = local_var_ji / (scalar_t)D;\n }\n __syncthreads();\n \n scalar_t var_ji = s_var_ji[0];\n scalar_t inv_std_ji = rsqrtf(var_ji + (scalar_t)1e-5);\n \n // Compute normalized and symmetrized values\n for (int d = tid; d < D; d += blockDim.x) {\n scalar_t norm_ij = (s_ij[d] - mean_ij) * inv_std_ij;\n scalar_t norm_ji = (s_ji[d] - mean_ji) * inv_std_ji;\n \n scalar_t val_ij = norm_ij * (scalar_t)weight[d] + (scalar_t)bias[d];\n scalar_t val_ji = norm_ji * (scalar_t)weight[d] + (scalar_t)bias[d];\n output[idx_ij + d] = (val_ij + val_ji) * (scalar_t)0.5;\n }\n}\n\ntorch::Tensor fused_layernorm_symmetrize_cuda(\n torch::Tensor input,\n torch::Tensor weight,\n torch::Tensor bias) {\n \n auto output = torch::empty_like(input);\n int B = input.size(0);\n int N = input.size(1);\n int D = input.size(3);\n \n dim3 grid(N, N, B);\n int threads = min(256, D);\n // Shared memory: 2*D scalar_t for s_ij and s_ji, 4 scalar_t for means and vars\n size_t shared_mem = (D * 2 + 4) * input.element_size();\n \n AT_DISPATCH_FLOATING_TYPES(input.scalar_type(), \"fused_layernorm_symmetrize\", ([&] {\n fused_layernorm_symmetrize_kernel<scalar_t><<<grid, threads, shared_mem>>>(\n input.data_ptr<scalar_t>(),\n output.data_ptr<scalar_t>(),\n B, N, D,\n weight.data_ptr<float>(),\n bias.data_ptr<float>());\n }));\n \n return output;\n}\n\"\"\"\n\ncpp_source = \"\"\"\ntorch::Tensor fused_layernorm_symmetrize_cuda(torch::Tensor input, torch::Tensor weight, torch::Tensor bias);\n\"\"\"\n\n# Load the extension\nfused_ops = load_inline(\n name='fused_ops',\n cpp_sources=cpp_source,\n cuda_sources=cuda_source,\n functions=['fused_layernorm_symmetrize_cuda'],\n verbose=True\n)\n\n\nclass ModelNew(nn.Module):\n def __init__(self, pair_dim, hidden_dim):\n super().__init__()\n self.layer_norm = nn.LayerNorm(pair_dim)\n self.mlp = nn.Sequential(\n nn.Linear(pair_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, 1),\n )\n self.pair_dim = pair_dim\n self.hidden_dim = hidden_dim\n\n def forward(self, pair_repr: torch.Tensor) -> torch.Tensor:\n # Use fused kernel for LayerNorm + symmetrize\n x = fused_ops.fused_layernorm_symmetrize_cuda(\n pair_repr, \n self.layer_norm.weight, \n self.layer_norm.bias\n )\n logits = self.mlp(x).squeeze(-1)\n return torch.sigmoid(logits)" }
result ← run_correctness (fail)
failrun_correctness run_correctness FAILED: (0 / 5) trials did not all match. correctness_issue: Output mismatch max_difference: ['0.014363', '0.017218', '0.029618', '0.021994', '0.019088'] avg_difference: ['0.001765', '0.001762', '0.001783', '0.001783', '0.001818']
turn 8166.6s turn · 11.4s LLM · 155.2s tools · tokens in=12,670 · out=1,258 · total=13,928
reasoning
The errors are still the same. The issue must be with the reduction using warp shuffle. When blockDim.x doesn't divide evenly by 32 (warp size), or when the number of threads is not a power of 2, the reduction might not work correctly. Let me try a different approach - using a simpler but correct version. The simple single-thread approach worked, so let me try with a slightly more parallel approach but ensuring correctness.
assistantLet me go back to the simple working version and try a different optimization strategy - maybe fusing the sigmoid as well:
call → compile_kernel
callcompile_kernel { "kernel_code": "import torch\nimport torch.nn as nn\nimport torch.nn.functional as F\nfrom torch.utils.cpp_extension import load_inline\n\n# Define the CUDA kernel - start simple and correct\ncuda_source = \"\"\"\n#include <torch/extension.h>\n#include <cuda.h>\n#include <cuda_runtime.h>\n#include <c10/cuda/CUDAException.h>\n\ntemplate <typename scalar_t>\n__global__ void fused_layernorm_symmetrize_kernel(\n const scalar_t* __restrict__ input,\n scalar_t* __restrict__ output,\n int B, int N, int D,\n const float* __restrict__ weight,\n const float* __restrict__ bias) {\n \n // input: (B, N, N, D)\n // output: (B, N, N, D)\n \n int b = blockIdx.z;\n int i = blockIdx.y;\n int j = blockIdx.x;\n \n if (b >= B || i >= N || j >= N) return;\n \n // Each thread handles one (i,j) position\n int idx_ij = b * N * N * D + i * N * D + j * D;\n int idx_ji = b * N * N * D + j * N * D + i * D;\n \n // Load and compute mean/variance for position (i,j)\n scalar_t mean = 0;\n for (int d = 0; d < D; d++) {\n mean += input[idx_ij + d];\n }\n mean /= D;\n \n scalar_t var = 0;\n for (int d = 0; d < D; d++) {\n scalar_t diff = input[idx_ij + d] - mean;\n var += diff * diff;\n }\n var /= D;\n \n scalar_t inv_std = rsqrtf(var + 1e-5f);\n \n // Load and compute mean/variance for position (j,i)\n scalar_t mean_ji = 0;\n for (int d = 0; d < D; d++) {\n mean_ji += input[idx_ji + d];\n }\n mean_ji /= D;\n \n scalar_t var_ji = 0;\n for (int d = 0; d < D; d++) {\n scalar_t diff = input[idx_ji + d] - mean_ji;\n var_ji += diff * diff;\n }\n var_ji /= D;\n \n scalar_t inv_std_ji = rsqrtf(var_ji + 1e-5f);\n \n // Compute normalized and symmetrized values\n for (int d = 0; d < D; d++) {\n scalar_t norm_ij = (input[idx_ij + d] - mean) * inv_std;\n scalar_t norm_ji = (input[idx_ji + d] - mean_ji) * inv_std_ji;\n \n // Apply weight/bias and symmetrize\n scalar_t val_ij = norm_ij * weight[d] + bias[d];\n scalar_t val_ji = norm_ji * weight[d] + bias[d];\n output[idx_ij + d] = (val_ij + val_ji) * 0.5f;\n }\n}\n\ntorch::Tensor fused_layernorm_symmetrize_cuda(\n torch::Tensor input,\n torch::Tensor weight,\n torch::Tensor bias) {\n \n auto output = torch::empty_like(input);\n int B = input.size(0);\n int N = input.size(1);\n int D = input.size(3);\n \n dim3 grid(N, N, B);\n dim3 block(1);\n \n AT_DISPATCH_FLOATING_TYPES(input.scalar_type(), \"fused_layernorm_symmetrize\", ([&] {\n fused_layernorm_symmetrize_kernel<scalar_t><<<grid, block>>>(\n input.data_ptr<scalar_t>(),\n output.data_ptr<scalar_t>(),\n B, N, D,\n weight.data_ptr<float>(),\n bias.data_ptr<float>());\n }));\n \n return output;\n}\n\"\"\"\n\ncpp_source = \"\"\"\ntorch::Tensor fused_layernorm_symmetrize_cuda(torch::Tensor input, torch::Tensor weight, torch::Tensor bias);\n\"\"\"\n\n# Load the extension\nfused_ops = load_inline(\n name='fused_ops',\n cpp_sources=cpp_source,\n cuda_sources=cuda_source,\n functions=['fused_layernorm_symmetrize_cuda'],\n verbose=True\n)\n\n\nclass ModelNew(nn.Module):\n def __init__(self, pair_dim, hidden_dim):\n super().__init__()\n self.layer_norm = nn.LayerNorm(pair_dim)\n self.mlp = nn.Sequential(\n nn.Linear(pair_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, 1),\n )\n self.pair_dim = pair_dim\n self.hidden_dim = hidden_dim\n\n def forward(self, pair_repr: torch.Tensor) -> torch.Tensor:\n # Use fused kernel for LayerNorm + symmetrize\n x = fused_ops.fused_layernorm_symmetrize_cuda(\n pair_repr, \n self.layer_norm.weight, \n self.layer_norm.bias\n )\n logits = self.mlp(x).squeeze(-1)\n return torch.sigmoid(logits)" }
result ← compile_kernel (ok)
okcompile_kernel compile_kernel PASSED: kernel compiled without errors.
turn 9346.4s turn · 10.5s LLM · 335.9s tools · tokens in=13,884 · out=1,165 · total=15,049
reasoning
Let me first verify the correctness of my kernel, then decide whether to submit or optimize further.
call → run_correctness
callrun_correctness { "kernel_code": "import torch\nimport torch.nn as nn\nimport torch.nn.functional as F\nfrom torch.utils.cpp_extension import load_inline\n\n# Define the CUDA kernel - start simple and correct\ncuda_source = \"\"\"\n#include <torch/extension.h>\n#include <cuda.h>\n#include <cuda_runtime.h>\n#include <c10/cuda/CUDAException.h>\n\ntemplate <typename scalar_t>\n__global__ void fused_layernorm_symmetrize_kernel(\n const scalar_t* __restrict__ input,\n scalar_t* __restrict__ output,\n int B, int N, int D,\n const float* __restrict__ weight,\n const float* __restrict__ bias) {\n \n // input: (B, N, N, D)\n // output: (B, N, N, D)\n \n int b = blockIdx.z;\n int i = blockIdx.y;\n int j = blockIdx.x;\n \n if (b >= B || i >= N || j >= N) return;\n \n // Each thread handles one (i,j) position\n int idx_ij = b * N * N * D + i * N * D + j * D;\n int idx_ji = b * N * N * D + j * N * D + i * D;\n \n // Load and compute mean/variance for position (i,j)\n scalar_t mean = 0;\n for (int d = 0; d < D; d++) {\n mean += input[idx_ij + d];\n }\n mean /= D;\n \n scalar_t var = 0;\n for (int d = 0; d < D; d++) {\n scalar_t diff = input[idx_ij + d] - mean;\n var += diff * diff;\n }\n var /= D;\n \n scalar_t inv_std = rsqrtf(var + 1e-5f);\n \n // Load and compute mean/variance for position (j,i)\n scalar_t mean_ji = 0;\n for (int d = 0; d < D; d++) {\n mean_ji += input[idx_ji + d];\n }\n mean_ji /= D;\n \n scalar_t var_ji = 0;\n for (int d = 0; d < D; d++) {\n scalar_t diff = input[idx_ji + d] - mean_ji;\n var_ji += diff * diff;\n }\n var_ji /= D;\n \n scalar_t inv_std_ji = rsqrtf(var_ji + 1e-5f);\n \n // Compute normalized and symmetrized values\n for (int d = 0; d < D; d++) {\n scalar_t norm_ij = (input[idx_ij + d] - mean) * inv_std;\n scalar_t norm_ji = (input[idx_ji + d] - mean_ji) * inv_std_ji;\n \n // Apply weight/bias and symmetrize\n scalar_t val_ij = norm_ij * weight[d] + bias[d];\n scalar_t val_ji = norm_ji * weight[d] + bias[d];\n output[idx_ij + d] = (val_ij + val_ji) * 0.5f;\n }\n}\n\ntorch::Tensor fused_layernorm_symmetrize_cuda(\n torch::Tensor input,\n torch::Tensor weight,\n torch::Tensor bias) {\n \n auto output = torch::empty_like(input);\n int B = input.size(0);\n int N = input.size(1);\n int D = input.size(3);\n \n dim3 grid(N, N, B);\n dim3 block(1);\n \n AT_DISPATCH_FLOATING_TYPES(input.scalar_type(), \"fused_layernorm_symmetrize\", ([&] {\n fused_layernorm_symmetrize_kernel<scalar_t><<<grid, block>>>(\n input.data_ptr<scalar_t>(),\n output.data_ptr<scalar_t>(),\n B, N, D,\n weight.data_ptr<float>(),\n bias.data_ptr<float>());\n }));\n \n return output;\n}\n\"\"\"\n\ncpp_source = \"\"\"\ntorch::Tensor fused_layernorm_symmetrize_cuda(torch::Tensor input, torch::Tensor weight, torch::Tensor bias);\n\"\"\"\n\n# Load the extension\nfused_ops = load_inline(\n name='fused_ops',\n cpp_sources=cpp_source,\n cuda_sources=cuda_source,\n functions=['fused_layernorm_symmetrize_cuda'],\n verbose=True\n)\n\n\nclass ModelNew(nn.Module):\n def __init__(self, pair_dim, hidden_dim):\n super().__init__()\n self.layer_norm = nn.LayerNorm(pair_dim)\n self.mlp = nn.Sequential(\n nn.Linear(pair_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, 1),\n )\n self.pair_dim = pair_dim\n self.hidden_dim = hidden_dim\n\n def forward(self, pair_repr: torch.Tensor) -> torch.Tensor:\n # Use fused kernel for LayerNorm + symmetrize\n x = fused_ops.fused_layernorm_symmetrize_cuda(\n pair_repr, \n self.layer_norm.weight, \n self.layer_norm.bias\n )\n logits = self.mlp(x).squeeze(-1)\n return torch.sigmoid(logits)" }
result ← run_correctness (ok)
okrun_correctness run_correctness PASSED: (5 / 5) trials all matched the reference. Numerical precision: max_abs_err=1.19e-07 mean_abs_err=2.41e-09 max_rel_err=2.32e-07
turn 10 · 1 tool call(s) · 12.6s turn · 12.6s LLM · tokens in=15,119 · out=1,349 · total=16,468 · 6,237 chars (click to expand)
turn 1012.6s turn · 12.6s LLM · 0.0s tools · tokens in=15,119 · out=1,349 · total=16,468
reasoning
The correctness test passed. The kernel is correct but it's using only 1 thread per block, which is very inefficient. I should try to optimize it before submitting, or submit now if the simple version already provides speedup. Let me think about this: 1. The simple version uses 1 thread per block - extremely low occupancy 2. But it fuses LayerNorm and symmetrization, which might still be faster than separate PyTorch operations Given the status message suggesting to submit if I have a correct kernel, I should consider if my current kernel is good enough. The fused operations eliminate memory traffic for intermediate results, which is a significant benefit. Actually, let me try one more optimization - using warp-level parallelism properly without shared memory complications that caused issues before. I'll use a simpler approach with proper warp shuffles. Wait, I should check the static analysis first before submitting to ensure there are no issues.
assistantLet me do a static check before submitting to ensure there are no issues:
call → static_check
callstatic_check { "kernel_code": "import torch\nimport torch.nn as nn\nimport torch.nn.functional as F\nfrom torch.utils.cpp_extension import load_inline\n\n# Define the CUDA kernel - start simple and correct\ncuda_source = \"\"\"\n#include <torch/extension.h>\n#include <cuda.h>\n#include <cuda_runtime.h>\n#include <c10/cuda/CUDAException.h>\n\ntemplate <typename scalar_t>\n__global__ void fused_layernorm_symmetrize_kernel(\n const scalar_t* __restrict__ input,\n scalar_t* __restrict__ output,\n int B, int N, int D,\n const float* __restrict__ weight,\n const float* __restrict__ bias) {\n \n // input: (B, N, N, D)\n // output: (B, N, N, D)\n \n int b = blockIdx.z;\n int i = blockIdx.y;\n int j = blockIdx.x;\n \n if (b >= B || i >= N || j >= N) return;\n \n // Each thread handles one (i,j) position\n int idx_ij = b * N * N * D + i * N * D + j * D;\n int idx_ji = b * N * N * D + j * N * D + i * D;\n \n // Load and compute mean/variance for position (i,j)\n scalar_t mean = 0;\n for (int d = 0; d < D; d++) {\n mean += input[idx_ij + d];\n }\n mean /= D;\n \n scalar_t var = 0;\n for (int d = 0; d < D; d++) {\n scalar_t diff = input[idx_ij + d] - mean;\n var += diff * diff;\n }\n var /= D;\n \n scalar_t inv_std = rsqrtf(var + 1e-5f);\n \n // Load and compute mean/variance for position (j,i)\n scalar_t mean_ji = 0;\n for (int d = 0; d < D; d++) {\n mean_ji += input[idx_ji + d];\n }\n mean_ji /= D;\n \n scalar_t var_ji = 0;\n for (int d = 0; d < D; d++) {\n scalar_t diff = input[idx_ji + d] - mean_ji;\n var_ji += diff * diff;\n }\n var_ji /= D;\n \n scalar_t inv_std_ji = rsqrtf(var_ji + 1e-5f);\n \n // Compute normalized and symmetrized values\n for (int d = 0; d < D; d++) {\n scalar_t norm_ij = (input[idx_ij + d] - mean) * inv_std;\n scalar_t norm_ji = (input[idx_ji + d] - mean_ji) * inv_std_ji;\n \n // Apply weight/bias and symmetrize\n scalar_t val_ij = norm_ij * weight[d] + bias[d];\n scalar_t val_ji = norm_ji * weight[d] + bias[d];\n output[idx_ij + d] = (val_ij + val_ji) * 0.5f;\n }\n}\n\ntorch::Tensor fused_layernorm_symmetrize_cuda(\n torch::Tensor input,\n torch::Tensor weight,\n torch::Tensor bias) {\n \n auto output = torch::empty_like(input);\n int B = input.size(0);\n int N = input.size(1);\n int D = input.size(3);\n \n dim3 grid(N, N, B);\n dim3 block(1);\n \n AT_DISPATCH_FLOATING_TYPES(input.scalar_type(), \"fused_layernorm_symmetrize\", ([&] {\n fused_layernorm_symmetrize_kernel<scalar_t><<<grid, block>>>(\n input.data_ptr<scalar_t>(),\n output.data_ptr<scalar_t>(),\n B, N, D,\n weight.data_ptr<float>(),\n bias.data_ptr<float>());\n }));\n \n return output;\n}\n\"\"\"\n\ncpp_source = \"\"\"\ntorch::Tensor fused_layernorm_symmetrize_cuda(torch::Tensor input, torch::Tensor weight, torch::Tensor bias);\n\"\"\"\n\n# Load the extension\nfused_ops = load_inline(\n name='fused_ops',\n cpp_sources=cpp_source,\n cuda_sources=cuda_source,\n functions=['fused_layernorm_symmetrize_cuda'],\n verbose=True\n)\n\n\nclass ModelNew(nn.Module):\n def __init__(self, pair_dim, hidden_dim):\n super().__init__()\n self.layer_norm = nn.LayerNorm(pair_dim)\n self.mlp = nn.Sequential(\n nn.Linear(pair_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, 1),\n )\n self.pair_dim = pair_dim\n self.hidden_dim = hidden_dim\n\n def forward(self, pair_repr: torch.Tensor) -> torch.Tensor:\n # Use fused kernel for LayerNorm + symmetrize\n x = fused_ops.fused_layernorm_symmetrize_cuda(\n pair_repr, \n self.layer_norm.weight, \n self.layer_norm.bias\n )\n logits = self.mlp(x).squeeze(-1)\n return torch.sigmoid(logits)" }
result ← static_check (ok)
okstatic_check static_check PASSED (with warnings): WARNING: Uses torch.nn compute layer (only containers, Parameter, init allowed) WARNING: Uses torch computation op: torch.sigmoid
turn 11562.3s turn · 10.9s LLM · 551.3s tools · tokens in=16,406 · out=1,231 · total=17,637 · FINAL
reasoning
The static check passed with warnings about using torch.nn compute layer and torch.sigmoid. The correctness test passed with very good numerical precision. The kernel is ready for submission. The kernel fuses LayerNorm and symmetrization into a single CUDA kernel, which should reduce memory bandwidth and kernel launch overhead compared to the original PyTorch implementation. The MLP and sigmoid are still using PyTorch operations. Let me submit the kernel.
call → submit_kernel
callsubmit_kernel { "kernel_code": "import torch\nimport torch.nn as nn\nimport torch.nn.functional as F\nfrom torch.utils.cpp_extension import load_inline\n\n# Define the CUDA kernel - start simple and correct\ncuda_source = \"\"\"\n#include <torch/extension.h>\n#include <cuda.h>\n#include <cuda_runtime.h>\n#include <c10/cuda/CUDAException.h>\n\ntemplate <typename scalar_t>\n__global__ void fused_layernorm_symmetrize_kernel(\n const scalar_t* __restrict__ input,\n scalar_t* __restrict__ output,\n int B, int N, int D,\n const float* __restrict__ weight,\n const float* __restrict__ bias) {\n \n // input: (B, N, N, D)\n // output: (B, N, N, D)\n \n int b = blockIdx.z;\n int i = blockIdx.y;\n int j = blockIdx.x;\n \n if (b >= B || i >= N || j >= N) return;\n \n // Each thread handles one (i,j) position\n int idx_ij = b * N * N * D + i * N * D + j * D;\n int idx_ji = b * N * N * D + j * N * D + i * D;\n \n // Load and compute mean/variance for position (i,j)\n scalar_t mean = 0;\n for (int d = 0; d < D; d++) {\n mean += input[idx_ij + d];\n }\n mean /= D;\n \n scalar_t var = 0;\n for (int d = 0; d < D; d++) {\n scalar_t diff = input[idx_ij + d] - mean;\n var += diff * diff;\n }\n var /= D;\n \n scalar_t inv_std = rsqrtf(var + 1e-5f);\n \n // Load and compute mean/variance for position (j,i)\n scalar_t mean_ji = 0;\n for (int d = 0; d < D; d++) {\n mean_ji += input[idx_ji + d];\n }\n mean_ji /= D;\n \n scalar_t var_ji = 0;\n for (int d = 0; d < D; d++) {\n scalar_t diff = input[idx_ji + d] - mean_ji;\n var_ji += diff * diff;\n }\n var_ji /= D;\n \n scalar_t inv_std_ji = rsqrtf(var_ji + 1e-5f);\n \n // Compute normalized and symmetrized values\n for (int d = 0; d < D; d++) {\n scalar_t norm_ij = (input[idx_ij + d] - mean) * inv_std;\n scalar_t norm_ji = (input[idx_ji + d] - mean_ji) * inv_std_ji;\n \n // Apply weight/bias and symmetrize\n scalar_t val_ij = norm_ij * weight[d] + bias[d];\n scalar_t val_ji = norm_ji * weight[d] + bias[d];\n output[idx_ij + d] = (val_ij + val_ji) * 0.5f;\n }\n}\n\ntorch::Tensor fused_layernorm_symmetrize_cuda(\n torch::Tensor input,\n torch::Tensor weight,\n torch::Tensor bias) {\n \n auto output = torch::empty_like(input);\n int B = input.size(0);\n int N = input.size(1);\n int D = input.size(3);\n \n dim3 grid(N, N, B);\n dim3 block(1);\n \n AT_DISPATCH_FLOATING_TYPES(input.scalar_type(), \"fused_layernorm_symmetrize\", ([&] {\n fused_layernorm_symmetrize_kernel<scalar_t><<<grid, block>>>(\n input.data_ptr<scalar_t>(),\n output.data_ptr<scalar_t>(),\n B, N, D,\n weight.data_ptr<float>(),\n bias.data_ptr<float>());\n }));\n \n return output;\n}\n\"\"\"\n\ncpp_source = \"\"\"\ntorch::Tensor fused_layernorm_symmetrize_cuda(torch::Tensor input, torch::Tensor weight, torch::Tensor bias);\n\"\"\"\n\n# Load the extension\nfused_ops = load_inline(\n name='fused_ops',\n cpp_sources=cpp_source,\n cuda_sources=cuda_source,\n functions=['fused_layernorm_symmetrize_cuda'],\n verbose=True\n)\n\n\nclass ModelNew(nn.Module):\n def __init__(self, pair_dim, hidden_dim):\n super().__init__()\n self.layer_norm = nn.LayerNorm(pair_dim)\n self.mlp = nn.Sequential(\n nn.Linear(pair_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, hidden_dim),\n nn.ReLU(),\n nn.Linear(hidden_dim, 1),\n )\n self.pair_dim = pair_dim\n self.hidden_dim = hidden_dim\n\n def forward(self, pair_repr: torch.Tensor) -> torch.Tensor:\n # Use fused kernel for LayerNorm + symmetrize\n x = fused_ops.fused_layernorm_symmetrize_cuda(\n pair_repr, \n self.layer_norm.weight, \n self.layer_norm.bias\n )\n logits = self.mlp(x).squeeze(-1)\n return torch.sigmoid(logits)" }
result ← submit_kernel (ok)
oksubmit_kernel submit_kernel PASSED: (5 / 5) correctness trials all passed. Kernel runtime: 0.11 μs Runtime stats: mean=0.11μs median=?μs std=0.03μs Numerical precision: max_abs_err=1.19e-07 mean_abs_err=2.41e-09 max_rel_err=2.32e-07 Memory: 45.03 MB peak (ref: 49.03 MB, ratio: 0.9184x) Kernel launches: 32 (ref: 41, fusion ratio: 1.2812) Energy: 241.38 mJ/run (ref: 0.0 mJ, ratio: -1x) SOL score: 0.1252
final submitted kernel (140 lines, 3,950 chars) — full source, untruncated
import torch
import torch.nn as nn
import torch.nn.functional as F
from torch.utils.cpp_extension import load_inline

# Define the CUDA kernel - start simple and correct
cuda_source = """
#include <torch/extension.h>
#include <cuda.h>
#include <cuda_runtime.h>
#include <c10/cuda/CUDAException.h>

template <typename scalar_t>
__global__ void fused_layernorm_symmetrize_kernel(
    const scalar_t* __restrict__ input,
    scalar_t* __restrict__ output,
    int B, int N, int D,
    const float* __restrict__ weight,
    const float* __restrict__ bias) {
    
    // input: (B, N, N, D)
    // output: (B, N, N, D)
    
    int b = blockIdx.z;
    int i = blockIdx.y;
    int j = blockIdx.x;
    
    if (b >= B || i >= N || j >= N) return;
    
    // Each thread handles one (i,j) position
    int idx_ij = b * N * N * D + i * N * D + j * D;
    int idx_ji = b * N * N * D + j * N * D + i * D;
    
    // Load and compute mean/variance for position (i,j)
    scalar_t mean = 0;
    for (int d = 0; d < D; d++) {
        mean += input[idx_ij + d];
    }
    mean /= D;
    
    scalar_t var = 0;
    for (int d = 0; d < D; d++) {
        scalar_t diff = input[idx_ij + d] - mean;
        var += diff * diff;
    }
    var /= D;
    
    scalar_t inv_std = rsqrtf(var + 1e-5f);
    
    // Load and compute mean/variance for position (j,i)
    scalar_t mean_ji = 0;
    for (int d = 0; d < D; d++) {
        mean_ji += input[idx_ji + d];
    }
    mean_ji /= D;
    
    scalar_t var_ji = 0;
    for (int d = 0; d < D; d++) {
        scalar_t diff = input[idx_ji + d] - mean_ji;
        var_ji += diff * diff;
    }
    var_ji /= D;
    
    scalar_t inv_std_ji = rsqrtf(var_ji + 1e-5f);
    
    // Compute normalized and symmetrized values
    for (int d = 0; d < D; d++) {
        scalar_t norm_ij = (input[idx_ij + d] - mean) * inv_std;
        scalar_t norm_ji = (input[idx_ji + d] - mean_ji) * inv_std_ji;
        
        // Apply weight/bias and symmetrize
        scalar_t val_ij = norm_ij * weight[d] + bias[d];
        scalar_t val_ji = norm_ji * weight[d] + bias[d];
        output[idx_ij + d] = (val_ij + val_ji) * 0.5f;
    }
}

torch::Tensor fused_layernorm_symmetrize_cuda(
    torch::Tensor input,
    torch::Tensor weight,
    torch::Tensor bias) {
    
    auto output = torch::empty_like(input);
    int B = input.size(0);
    int N = input.size(1);
    int D = input.size(3);
    
    dim3 grid(N, N, B);
    dim3 block(1);
    
    AT_DISPATCH_FLOATING_TYPES(input.scalar_type(), "fused_layernorm_symmetrize", ([&] {
        fused_layernorm_symmetrize_kernel<scalar_t><<<grid, block>>>(
            input.data_ptr<scalar_t>(),
            output.data_ptr<scalar_t>(),
            B, N, D,
            weight.data_ptr<float>(),
            bias.data_ptr<float>());
    }));
    
    return output;
}
"""

cpp_source = """
torch::Tensor fused_layernorm_symmetrize_cuda(torch::Tensor input, torch::Tensor weight, torch::Tensor bias);
"""

# Load the extension
fused_ops = load_inline(
    name='fused_ops',
    cpp_sources=cpp_source,
    cuda_sources=cuda_source,
    functions=['fused_layernorm_symmetrize_cuda'],
    verbose=True
)


class ModelNew(nn.Module):
    def __init__(self, pair_dim, hidden_dim):
        super().__init__()
        self.layer_norm = nn.LayerNorm(pair_dim)
        self.mlp = nn.Sequential(
            nn.Linear(pair_dim, hidden_dim),
            nn.ReLU(),
            nn.Linear(hidden_dim, hidden_dim),
            nn.ReLU(),
            nn.Linear(hidden_dim, 1),
        )
        self.pair_dim = pair_dim
        self.hidden_dim = hidden_dim

    def forward(self, pair_repr: torch.Tensor) -> torch.Tensor:
        # Use fused kernel for LayerNorm + symmetrize
        x = fused_ops.fused_layernorm_symmetrize_cuda(
            pair_repr, 
            self.layer_norm.weight, 
            self.layer_norm.bias
        )
        logits = self.mlp(x).squeeze(-1)
        return torch.sigmoid(logits)