You write custom CUDA kernels to replace the pytorch operators in the given architecture to get speedups. 

    You have complete freedom to choose the set of operators you want to replace. You may make the decision to replace some operators with custom CUDA kernels and leave others unchanged. You may replace multiple operators with custom implementations, consider operator fusion opportunities (combining multiple operators into a single kernel, for example, combining matmul+relu), or algorithmic changes (such as online softmax). You are only limited by your imagination.


        Here's an example to show you the syntax of inline embedding custom CUDA operators in torch: The example given architecture is: 

        ``` 

        # Copyright (c) Meta Platforms, Inc. and affiliates.
# All rights reserved.
#
# This source code is licensed under the license found in the
# LICENSE file in the root directory of this source tree.

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


class Model(nn.Module):
    def __init__(self) -> None:
        super().__init__()

    def forward(self, a, b):
        return a + b


def get_inputs():
    # randomly generate input tensors based on the model architecture
    a = torch.randn(1, 128).cuda()
    b = torch.randn(1, 128).cuda()
    return [a, b]


def get_init_inputs():
    # randomly generate tensors required for initialization based on the model architecture
    return []

        ``` 

        The example new arch with custom CUDA kernels looks like this:
        ```
        # Copyright (c) Meta Platforms, Inc. and affiliates.
# All rights reserved.
#
# This source code is licensed under the license found in the
# LICENSE file in the root directory of this source tree.

import torch
import torch.nn as nn
import torch.nn.functional as F
from torch.utils.cpp_extension import load_inline

# Define the custom CUDA kernel for element-wise addition
elementwise_add_source = """
#include <torch/extension.h>
#include <cuda_runtime.h>

__global__ void elementwise_add_kernel(const float* a, const float* b, float* out, int size) {
    int idx = blockIdx.x * blockDim.x + threadIdx.x;
    if (idx < size) {
        out[idx] = a[idx] + b[idx];
    }
}

torch::Tensor elementwise_add_cuda(torch::Tensor a, torch::Tensor b) {
    auto size = a.numel();
    auto out = torch::zeros_like(a);

    const int block_size = 256;
    const int num_blocks = (size + block_size - 1) / block_size;

    elementwise_add_kernel<<<num_blocks, block_size>>>(a.data_ptr<float>(), b.data_ptr<float>(), out.data_ptr<float>(), size);

    return out;
}
"""

elementwise_add_cpp_source = (
    "torch::Tensor elementwise_add_cuda(torch::Tensor a, torch::Tensor b);"
)

# Compile the inline CUDA code for element-wise addition
elementwise_add = load_inline(
    name="elementwise_add",
    cpp_sources=elementwise_add_cpp_source,
    cuda_sources=elementwise_add_source,
    functions=["elementwise_add_cuda"],
    verbose=True,
    extra_cflags=[""],
    extra_ldflags=[""],
)


class ModelNew(nn.Module):
    def __init__(self) -> None:
        super().__init__()
        self.elementwise_add = elementwise_add

    def forward(self, a, b):
        return self.elementwise_add.elementwise_add_cuda(a, b)

        ``` 

        
    You are given the following architecture: 

    ```
    # Copyright (c) Meta Platforms, Inc. and affiliates.
# All rights reserved.
#
# This source code is licensed under the license found in the
# LICENSE file in the root directory of this source tree.

import torch
import torch.nn as nn

class Model(nn.Module):
    """
    Model that performs a 3D transposed convolution, LogSumExp, HardSwish, subtraction, clamp operations.
    """
    def __init__(self, in_channels, out_channels, kernel_size, stride, padding, bias_shape):
        super(Model, self).__init__()
        self.conv_transpose = nn.ConvTranspose3d(in_channels, out_channels, kernel_size, stride=stride, padding=padding)
        self.bias = nn.Parameter(torch.randn(1, 1, 1, 1)) 

    def forward(self, x):
        x = self.conv_transpose(x)
        x = torch.logsumexp(x, dim=1, keepdim=True)
        x = x * torch.sigmoid(x + 3) / 6
        x = x - self.bias
        x = torch.clamp(x, min=-1, max=1)
        return x

batch_size = 128
in_channels = 3
out_channels = 16
depth, height, width = 16, 32, 32
kernel_size = 3
stride = 2
padding = 1
bias_shape = (1, 1, 1, 1)  

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

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

    ```
    
Optimize the architecture named Model with custom CUDA operators! Name your optimized output architecture ModelNew. Output the new code in codeblocks. Please generate real code, NOT pseudocode, make sure the code compiles and is fully functional. Just output the new model code, no other text, and NO testing code! 



This is turn 4 of 6 for this task. Your previous solution (turn 3) was:
```python
# Copyright (c) Meta Platforms, Inc. and affiliates.
# All rights reserved.
#
# This source code is licensed under the license found in the
# LICENSE file in the root directory of this source tree.

import torch
import torch.nn as nn
from torch.utils.cpp_extension import load_inline

# Fused CUDA kernel (v3):
#   ConvTranspose3d (stride=2, pad=1, k=3, CI=3, CO=16)
#   + LogSumExp over channel dim
#   + HardSwish (x * sigmoid(x + 3) / 6)
#   + scalar bias subtraction
#   + clamp([-1, 1])
#
# Key ideas:
#   * The (N, 16, D, H, W) transposed-conv result is never materialized; each
#     thread accumulates all 16 channels for two horizontally adjacent output
#     voxels (even/odd parity) directly in registers.
#   * All convolution weights / biases live in __constant__ memory. Because a
#     warp follows a uniform kernel-tap schedule, the weights are fetched via
#     constant-cache broadcasts (often folded directly into FMA operands),
#     completely eliminating the shared-memory load traffic that dominated the
#     previous versions.
#   * The tap schedule is specialized on the (od, oh) parity with warp-uniform
#     branches, so every constant-memory index is a compile-time immediate.
#   * A template specialization for ID=16, IH=32, IW=32 makes all address
#     arithmetic and thread decomposition compile-time (shifts / magic mults).
#   * Fast SFU transcendentals (__expf/__logf) in the epilogue.
_cuda_src = r"""
#include <torch/extension.h>
#include <ATen/cuda/CUDAContext.h>
#include <c10/cuda/CUDAException.h>
#include <cuda_runtime.h>

namespace {

constexpr int CI = 3;
constexpr int CO = 16;

__constant__ float CW[CI * CO * 27];  // conv weight layout [ci][co][kd][kh][kw]
__constant__ float CB[CO];            // conv bias
__constant__ float SB[1];             // scalar bias

__device__ __forceinline__ float epilogue(const float (&v)[CO], const float sb) {
    // numerically stable logsumexp over channels
    float m = v[0];
    #pragma unroll
    for (int c = 1; c < CO; ++c) m = fmaxf(m, v[c]);
    float s = 0.f;
    #pragma unroll
    for (int c = 0; c < CO; ++c) s += __expf(v[c] - m);
    const float lse = m + __logf(s);
    // hardswish as in the reference: x * sigmoid(x + 3) / 6
    const float sig = __frcp_rn(1.f + __expf(-(lse + 3.f)));
    const float y = fmaf(lse, sig * 0.16666666666666666f, -sb);
    return fminf(fmaxf(y, -1.f), 1.f);
}

// Accumulate the transposed-conv taps for a pair of adjacent voxels.
//   even voxel: (od, oh, 2*wp)     gathers (id, ih, wp)   with kw = 1
//   odd voxel:  (od, oh, 2*wp + 1) gathers (id, ih, wp+1) with kw = 0
//                                    and (id, ih, wp)     with kw = 2
// PD/PH encode the parity of od/oh: parity 0 -> tap index 1 only,
// parity 1 -> tap indices {0, 2}. All loops are fully unrolled so every
// constant-memory offset is an immediate.
template <int PD, int PH>
__device__ __forceinline__ void compute_pair(
    const float* __restrict__ xn,
    const int ID, const int IH, const int IW,
    const int od, const int oh, const int wp, const bool has1,
    float (&v0)[CO], float (&v1)[CO])
{
    constexpr int KDA = PD ? 0 : 1;
    constexpr int KDB = PD ? 2 : 1;
    constexpr int NKD = PD ? 2 : 1;
    constexpr int KHA = PH ? 0 : 1;
    constexpr int KHB = PH ? 2 : 1;
    constexpr int NKH = PH ? 2 : 1;

    const int HW  = IH * IW;
    const int DHW = ID * HW;

    #pragma unroll
    for (int ci = 0; ci < CI; ++ci) {
        const float* xc = xn + ci * DHW;
        #pragma unroll
        for (int kdi = 0; kdi < NKD; ++kdi) {
            const int kd = (kdi == 0) ? KDA : KDB;
            const int id = (od + 1 - kd) >> 1;   // always in [0, ID)
            const float* xd = xc + id * HW;
            #pragma unroll
            for (int khi = 0; khi < NKH; ++khi) {
                const int kh = (khi == 0) ? KHA : KHB;
                const int ih = (oh + 1 - kh) >> 1;   // always in [0, IH)
                const float* p = xd + ih * IW + wp;
                const float x0 = __ldg(p);
                const float* wpv = CW + ci * (CO * 27) + kd * 9 + kh * 3;
                #pragma unroll
                for (int co = 0; co < CO; ++co)
                    v0[co] = fmaf(x0, wpv[co * 27 + 1], v0[co]);
                if (has1) {
                    const float x1 = __ldg(p + 1);
                    #pragma unroll
                    for (int co = 0; co < CO; ++co)
                        v1[co] = fmaf(x1, wpv[co * 27], v1[co]);
                    #pragma unroll
                    for (int co = 0; co < CO; ++co)
                        v1[co] = fmaf(x0, wpv[co * 27 + 2], v1[co]);
                }
            }
        }
    }
}

__device__ __forceinline__ void body(const float* __restrict__ x,
                                     float* __restrict__ out,
                                     const int N,
                                     const int ID, const int IH, const int IW)
{
    const int OD = 2 * ID - 1;
    const int OH = 2 * IH - 1;
    const int OW = 2 * IW - 1;
    const int OW2 = IW;   // one thread covers voxels (2*wp, 2*wp+1)

    const int idx = blockIdx.x * blockDim.x + threadIdx.x;
    if (idx >= N * OD * OH * OW2) return;

    const int wp = idx % OW2;
    int t = idx / OW2;
    const int oh = t % OH; t /= OH;
    const int od = t % OD;
    const int n  = t / OD;

    const bool has1 = (wp + 1) < IW;

    float v0[CO], v1[CO];
    #pragma unroll
    for (int c = 0; c < CO; ++c) {
        const float b = CB[c];
        v0[c] = b;
        v1[c] = b;
    }

    const float* xn = x + static_cast<long long>(n) * (CI * ID * IH * IW);

    const int pd = od & 1;
    const int ph = oh & 1;
    if (pd == 0) {
        if (ph == 0) compute_pair<0, 0>(xn, ID, IH, IW, od, oh, wp, has1, v0, v1);
        else         compute_pair<0, 1>(xn, ID, IH, IW, od, oh, wp, has1, v0, v1);
    } else {
        if (ph == 0) compute_pair<1, 0>(xn, ID, IH, IW, od, oh, wp, has1, v0, v1);
        else         compute_pair<1, 1>(xn, ID, IH, IW, od, oh, wp, has1, v0, v1);
    }

    const float sb = SB[0];
    const long long oi =
        (static_cast<long long>(n * OD + od) * OH + oh) * OW + 2 * wp;
    const float y0 = epilogue(v0, sb);
    if (has1) {
        const float y1 = epilogue(v1, sb);
        *reinterpret_cast<float2*>(out + oi) = make_float2(y0, y1);
    } else {
        out[oi] = y0;
    }
}

template <int ID, int IH, int IW>
__global__ __launch_bounds__(256)
void fused_kernel_spec(const float* __restrict__ x, float* __restrict__ out,
                       const int N) {
    body(x, out, N, ID, IH, IW);
}

__global__ __launch_bounds__(256)
void fused_kernel_gen(const float* __restrict__ x, float* __restrict__ out,
                      const int N, const int ID, const int IH, const int IW) {
    body(x, out, N, ID, IH, IW);
}

}  // namespace

torch::Tensor fused_forward(torch::Tensor x, torch::Tensor weight,
                            torch::Tensor conv_bias, torch::Tensor scalar_bias) {
    TORCH_CHECK(x.is_cuda() && weight.is_cuda() && conv_bias.is_cuda() && scalar_bias.is_cuda(),
                "all tensors must be CUDA tensors");
    TORCH_CHECK(x.dim() == 5, "input must be 5D (N, C, D, H, W)");
    TORCH_CHECK(x.scalar_type() == at::kFloat && weight.scalar_type() == at::kFloat &&
                conv_bias.scalar_type() == at::kFloat && scalar_bias.scalar_type() == at::kFloat,
                "float32 tensors required");

    auto xc = x.contiguous();
    auto wc = weight.contiguous();
    auto bc = conv_bias.contiguous();
    auto bs = scalar_bias.contiguous();

    const int N  = static_cast<int>(xc.size(0));
    const int CIv = static_cast<int>(xc.size(1));
    const int ID = static_cast<int>(xc.size(2));
    const int IH = static_cast<int>(xc.size(3));
    const int IW = static_cast<int>(xc.size(4));
    const int COv = static_cast<int>(wc.size(1));

    TORCH_CHECK(CIv == CI && COv == CO, "unsupported channel configuration");
    TORCH_CHECK(wc.numel() == CI * CO * 27, "unsupported kernel size");
    TORCH_CHECK(bs.numel() == 1, "scalar bias must contain exactly one element");

    const int OD = 2 * ID - 1;
    const int OH = 2 * IH - 1;
    const int OW = 2 * IW - 1;

    auto out = torch::empty({N, 1, OD, OH, OW}, xc.options());

    auto stream = at::cuda::getCurrentCUDAStream();

    // stage weights / biases in constant memory (device-to-device, async)
    cudaMemcpyToSymbolAsync(CW, wc.data_ptr<float>(), sizeof(float) * CI * CO * 27,
                            0, cudaMemcpyDeviceToDevice, stream.stream());
    cudaMemcpyToSymbolAsync(CB, bc.data_ptr<float>(), sizeof(float) * CO,
                            0, cudaMemcpyDeviceToDevice, stream.stream());
    cudaMemcpyToSymbolAsync(SB, bs.data_ptr<float>(), sizeof(float),
                            0, cudaMemcpyDeviceToDevice, stream.stream());

    const int total = N * OD * OH * IW;
    const int grid = (total + 255) / 256;

    if (ID == 16 && IH == 32 && IW == 32) {
        fused_kernel_spec<16, 32, 32><<<grid, 256, 0, stream.stream()>>>(
            xc.data_ptr<float>(), out.data_ptr<float>(), N);
    } else {
        fused_kernel_gen<<<grid, 256, 0, stream.stream()>>>(
            xc.data_ptr<float>(), out.data_ptr<float>(), N, ID, IH, IW);
    }
    C10_CUDA_KERNEL_LAUNCH_CHECK();

    return out;
}
"""

_cpp_src = ("torch::Tensor fused_forward(torch::Tensor x, torch::Tensor weight, "
            "torch::Tensor conv_bias, torch::Tensor scalar_bias);")

_fused = load_inline(
    name="fused_convt3d_lse_hswish_clamp_v3",
    cpp_sources=_cpp_src,
    cuda_sources=_cuda_src,
    functions=["fused_forward"],
    verbose=False,
    extra_cuda_cflags=["-O3"],
)


class ModelNew(nn.Module):
    """
    Optimized model: a single fused CUDA kernel performs the 3D transposed
    convolution, LogSumExp, HardSwish, bias subtraction and clamp without
    materializing the large (N, C_out, D, H, W) intermediate tensor. Weights
    are staged in constant memory so they are consumed as broadcast FMA
    operands with no shared/global load traffic in the inner loop.
    """
    def __init__(self, in_channels, out_channels, kernel_size, stride, padding, bias_shape):
        super(ModelNew, self).__init__()
        self.conv_transpose = nn.ConvTranspose3d(in_channels, out_channels, kernel_size,
                                                 stride=stride, padding=padding)
        self.bias = nn.Parameter(torch.randn(1, 1, 1, 1))

    def forward(self, x):
        ct = self.conv_transpose
        if (
            x.is_cuda
            and x.dtype == torch.float32
            and ct.bias is not None
            and ct.groups == 1
            and ct.in_channels == 3
            and ct.out_channels == 16
            and tuple(ct.kernel_size) == (3, 3, 3)
            and tuple(ct.stride) == (2, 2, 2)
            and tuple(ct.padding) == (1, 1, 1)
            and tuple(ct.dilation) == (1, 1, 1)
            and tuple(ct.output_padding) == (0, 0, 0)
        ):
            return _fused.fused_forward(x, ct.weight, ct.bias, self.bias)

        # generic fallback (reference implementation)
        x = ct(x)
        x = torch.logsumexp(x, dim=1, keepdim=True)
        x = x * torch.sigmoid(x + 3) / 6
        x = x - self.bias
        x = torch.clamp(x, min=-1, max=1)
        return x
```
Result of evaluating it:
The solution compiled but failed at run time. Error output:
```
torch.AcceleratorError: CUDA error: misaligned address
Search for `cudaErrorMisalignedAddress' in https://docs.nvidia.com/cuda/cuda-runtime-api/group__CUDART__TYPES.html for more information.
CUDA kernel errors might be asynchronously reported at some other API call, so the stacktrace below might be incorrect.
For debugging consider passing CUDA_LAUNCH_BLOCKING=1
The CUDA driver logged these messages, which may provide useful details:
Sticky error detected
Returning 716 (CUDA_ERROR_MISALIGNED_ADDRESS) from cuMemFree_v2
```

The fastest correct solution so far (28.027x speedup) is:
```python
# Copyright (c) Meta Platforms, Inc. and affiliates.
# All rights reserved.
#
# This source code is licensed under the license found in the
# LICENSE file in the root directory of this source tree.

import torch
import torch.nn as nn
from torch.utils.cpp_extension import load_inline

# Fused CUDA kernel:
#   ConvTranspose3d (stride=2, pad=1, k=3, CI=3, CO=16)
#   + LogSumExp over channel dim
#   + HardSwish (x * sigmoid(x + 3) / 6)
#   + scalar bias subtraction
#   + clamp([-1, 1])
#
# The transposed convolution result for all 16 output channels of a given
# spatial voxel is accumulated in registers, so the huge (N, 16, D, H, W)
# intermediate tensor never touches global memory. Each thread produces two
# horizontally adjacent output voxels (even/odd parity), which makes the
# kernel-tap schedule uniform across a warp and keeps input loads coalesced.
_cuda_src = r"""
#include <torch/extension.h>
#include <ATen/cuda/CUDAContext.h>
#include <cuda_runtime.h>

namespace {

template <int CO>
__device__ __forceinline__ float lse_hswish_bias_clamp(float (&vals)[CO], float sb) {
    // numerically stable logsumexp over channels
    float m = vals[0];
    #pragma unroll
    for (int co = 1; co < CO; ++co) m = fmaxf(m, vals[co]);
    float s = 0.f;
    #pragma unroll
    for (int co = 0; co < CO; ++co) s += expf(vals[co] - m);
    const float lse = m + logf(s);
    // hardswish as written in the reference: x * sigmoid(x + 3) / 6
    const float sig = 1.f / (1.f + expf(-(lse + 3.f)));
    float y = (lse * sig) / 6.f - sb;
    return fminf(fmaxf(y, -1.f), 1.f);
}

template <int CI, int CO>
__global__ __launch_bounds__(256)
void fused_convt3d_lse_hswish_clamp_kernel(
    const float* __restrict__ x,          // (N, CI, ID, IH, IW)
    const float* __restrict__ w,          // (CI, CO, 3, 3, 3)
    const float* __restrict__ conv_bias,  // (CO)
    const float* __restrict__ scalar_bias,// (1)
    float* __restrict__ out,              // (N, 1, OD, OH, OW)
    const int N, const int ID, const int IH, const int IW,
    const int OD, const int OH, const int OW)
{
    constexpr int K = 3;
    constexpr int WSIZE = CI * CO * K * K * K;

    __shared__ float ws[WSIZE];
    __shared__ float cb[CO];

    for (int i = threadIdx.x; i < WSIZE; i += blockDim.x) ws[i] = w[i];
    if (threadIdx.x < CO) cb[threadIdx.x] = conv_bias[threadIdx.x];
    __syncthreads();

    const int OW2 = (OW + 1) >> 1;   // one thread handles a pair of adjacent voxels
    const int idx = blockIdx.x * blockDim.x + threadIdx.x;
    const int total = N * OD * OH * OW2;
    if (idx >= total) return;

    const int wp = idx % OW2;
    int t = idx / OW2;
    const int oh = t % OH; t /= OH;
    const int od = t % OD;
    const int n  = t / OD;

    const bool has1 = (wp * 2 + 1) < OW;   // second (odd) voxel exists

    float vals0[CO];
    float vals1[CO];
    #pragma unroll
    for (int co = 0; co < CO; ++co) {
        vals0[co] = cb[co];
        vals1[co] = cb[co];
    }

    // out[o] = sum_{i,k} in[i] * w[k], with o = 2*i - 1 + k  (stride=2, pad=1)
    // => i = (o + 1 - k) / 2, valid when (o + 1 - k) is even and in range.
    const int kd0 = (od + 1) & 1;
    const int kh0 = (oh + 1) & 1;

    const float* xn = x + static_cast<long long>(n) * (CI * ID * IH * IW);

    #pragma unroll
    for (int ci = 0; ci < CI; ++ci) {
        const float* xc = xn + ci * (ID * IH * IW);
        const float* wc = ws + ci * (CO * 27);
        for (int kd = kd0; kd < K; kd += 2) {
            const int id = (od + 1 - kd) >> 1;
            if (static_cast<unsigned>(id) >= static_cast<unsigned>(ID)) continue;
            const float* xd = xc + id * (IH * IW);
            const float* wd = wc + kd * 9;
            for (int kh = kh0; kh < K; kh += 2) {
                const int ih = (oh + 1 - kh) >> 1;
                if (static_cast<unsigned>(ih) >= static_cast<unsigned>(IH)) continue;
                const float* xh = xd + ih * IW;
                const float* wh = wd + kh * 3;

                // kw = 1 contributes to the even voxel (ow = 2*wp), iw = wp
                {
                    const float xv = __ldg(xh + wp);
                    const float* wv = wh + 1;
                    #pragma unroll
                    for (int co = 0; co < CO; ++co)
                        vals0[co] = fmaf(xv, wv[co * 27], vals0[co]);
                }
                if (has1) {
                    // kw = 0 contributes to the odd voxel, iw = wp + 1
                    {
                        const float xv = __ldg(xh + wp + 1);
                        const float* wv = wh;
                        #pragma unroll
                        for (int co = 0; co < CO; ++co)
                            vals1[co] = fmaf(xv, wv[co * 27], vals1[co]);
                    }
                    // kw = 2 contributes to the odd voxel, iw = wp
                    {
                        const float xv = __ldg(xh + wp);
                        const float* wv = wh + 2;
                        #pragma unroll
                        for (int co = 0; co < CO; ++co)
                            vals1[co] = fmaf(xv, wv[co * 27], vals1[co]);
                    }
                }
            }
        }
    }

    const float sb = __ldg(scalar_bias);
    const float y0 = lse_hswish_bias_clamp<CO>(vals0, sb);

    const int out_idx = ((n * OD + od) * OH + oh) * OW + wp * 2;
    out[out_idx] = y0;
    if (has1) {
        const float y1 = lse_hswish_bias_clamp<CO>(vals1, sb);
        out[out_idx + 1] = y1;
    }
}

}  // namespace

torch::Tensor fused_forward(torch::Tensor x, torch::Tensor weight,
                            torch::Tensor conv_bias, torch::Tensor scalar_bias) {
    TORCH_CHECK(x.is_cuda() && weight.is_cuda() && conv_bias.is_cuda() && scalar_bias.is_cuda(),
                "all tensors must be CUDA tensors");
    TORCH_CHECK(x.dim() == 5, "input must be 5D (N, C, D, H, W)");
    TORCH_CHECK(x.scalar_type() == at::kFloat && weight.scalar_type() == at::kFloat &&
                conv_bias.scalar_type() == at::kFloat && scalar_bias.scalar_type() == at::kFloat,
                "float32 tensors required");

    auto xc = x.contiguous();
    auto wc = weight.contiguous();
    auto bc = conv_bias.contiguous();
    auto bs = scalar_bias.contiguous();

    const int N  = static_cast<int>(xc.size(0));
    const int CI = static_cast<int>(xc.size(1));
    const int ID = static_cast<int>(xc.size(2));
    const int IH = static_cast<int>(xc.size(3));
    const int IW = static_cast<int>(xc.size(4));
    const int CO = static_cast<int>(wc.size(1));

    TORCH_CHECK(CI == 3 && CO == 16, "unsupported channel configuration");
    TORCH_CHECK(bs.numel() == 1, "scalar bias must contain exactly one element");

    // output size of ConvTranspose3d with stride=2, padding=1, kernel=3
    const int OD = (ID - 1) * 2 - 2 + 3;
    const int OH = (IH - 1) * 2 - 2 + 3;
    const int OW = (IW - 1) * 2 - 2 + 3;

    auto out = torch::empty({N, 1, OD, OH, OW}, xc.options());

    const int OW2 = (OW + 1) / 2;
    const int total = N * OD * OH * OW2;
    const int block = 256;
    const int grid = (total + block - 1) / block;

    auto stream = at::cuda::getCurrentCUDAStream();
    fused_convt3d_lse_hswish_clamp_kernel<3, 16><<<grid, block, 0, stream>>>(
        xc.data_ptr<float>(), wc.data_ptr<float>(), bc.data_ptr<float>(), bs.data_ptr<float>(),
        out.data_ptr<float>(), N, ID, IH, IW, OD, OH, OW);

    return out;
}
"""

_cpp_src = ("torch::Tensor fused_forward(torch::Tensor x, torch::Tensor weight, "
            "torch::Tensor conv_bias, torch::Tensor scalar_bias);")

_fused = load_inline(
    name="fused_convt3d_lse_hswish_clamp",
    cpp_sources=_cpp_src,
    cuda_sources=_cuda_src,
    functions=["fused_forward"],
    verbose=False,
    extra_cuda_cflags=["-O3"],
)


class ModelNew(nn.Module):
    """
    Optimized model: a single fused CUDA kernel performs the 3D transposed
    convolution, LogSumExp, HardSwish, bias subtraction and clamp without
    materializing the large (N, C_out, D, H, W) intermediate tensor.
    """
    def __init__(self, in_channels, out_channels, kernel_size, stride, padding, bias_shape):
        super(ModelNew, self).__init__()
        self.conv_transpose = nn.ConvTranspose3d(in_channels, out_channels, kernel_size,
                                                 stride=stride, padding=padding)
        self.bias = nn.Parameter(torch.randn(1, 1, 1, 1))

    def forward(self, x):
        ct = self.conv_transpose
        if (
            x.is_cuda
            and x.dtype == torch.float32
            and ct.bias is not None
            and ct.groups == 1
            and ct.in_channels == 3
            and ct.out_channels == 16
            and tuple(ct.kernel_size) == (3, 3, 3)
            and tuple(ct.stride) == (2, 2, 2)
            and tuple(ct.padding) == (1, 1, 1)
            and tuple(ct.dilation) == (1, 1, 1)
            and tuple(ct.output_padding) == (0, 0, 0)
        ):
            return _fused.fused_forward(x, ct.weight, ct.bias, self.bias)

        # generic fallback (reference implementation)
        x = ct(x)
        x = torch.logsumexp(x, dim=1, keepdim=True)
        x = x * torch.sigmoid(x + 3) / 6
        x = x - self.bias
        x = torch.clamp(x, min=-1, max=1)
        return x
```

Write an improved complete ModelNew: fix any error, and if the solution is already correct, make it faster. Output the new code in one ```python codeblock. Just the model code, no other text, and NO testing code.


Precision rule: the timed forward must compute in FP32 (TF32 tensor-core math is allowed, as in the PyTorch reference). Do not cast inputs, weights or intermediate results to float16 or bfloat16 and do not call half-precision kernels in forward().
