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


_cpp_source = r"""
#include <torch/extension.h>

torch::Tensor fused_post_conv(
    torch::Tensor x,
    torch::Tensor weight,
    torch::Tensor bias,
    int64_t groups,
    double eps);
"""

_cuda_source = r"""
#include <torch/extension.h>
#include <ATen/cuda/CUDAContext.h>
#include <c10/cuda/CUDAGuard.h>
#include <c10/cuda/CUDAException.h>
#include <cuda.h>
#include <cuda_runtime.h>
#include <math.h>

constexpr int THREADS = 256;
constexpr int CHUNK = 2048;

__device__ __forceinline__ float warp_sum(float value) {
    #pragma unroll
    for (int offset = 16; offset > 0; offset >>= 1)
        value += __shfl_down_sync(0xffffffff, value, offset);
    return value;
}

__device__ __forceinline__ float block_sum(float value) {
    __shared__ float shared[8];
    const int lane = threadIdx.x & 31;
    const int warp = threadIdx.x >> 5;

    value = warp_sum(value);
    if (lane == 0)
        shared[warp] = value;
    __syncthreads();

    value = threadIdx.x < 8 ? shared[lane] : 0.0f;
    if (warp == 0)
        value = warp_sum(value);
    return value;
}

__global__ void partial_stats_kernel(
    const float* __restrict__ x,
    float2* __restrict__ partial,
    int group_elements,
    int parts)
{
    const int group = blockIdx.x;
    const int part = blockIdx.y;
    const int begin = part * CHUNK;
    const int end = min(begin + CHUNK, group_elements);
    const float* base = x + static_cast<int64_t>(group) * group_elements;

    // Welford accumulation avoids cancellation in the variance.
    float mean = 0.0f;
    float m2 = 0.0f;
    int count = 0;

    for (int i = begin + threadIdx.x; i < end; i += THREADS) {
        const float value = base[i];
        ++count;
        const float delta = value - mean;
        mean += delta / count;
        m2 += delta * (value - mean);
    }

    #pragma unroll
    for (int offset = 16; offset > 0; offset >>= 1) {
        const float other_mean =
            __shfl_down_sync(0xffffffff, mean, offset);
        const float other_m2 =
            __shfl_down_sync(0xffffffff, m2, offset);
        const int other_count =
            __shfl_down_sync(0xffffffff, count, offset);

        if ((threadIdx.x & 31) + offset < 32 && other_count > 0) {
            const int total = count + other_count;
            const float delta = other_mean - mean;
            m2 += other_m2 +
                  delta * delta * (float(count) * other_count / total);
            mean += delta * (float(other_count) / total);
            count = total;
        }
    }

    __shared__ float means[8];
    __shared__ float moments[8];
    __shared__ int counts[8];

    if ((threadIdx.x & 31) == 0) {
        const int warp = threadIdx.x >> 5;
        means[warp] = mean;
        moments[warp] = m2;
        counts[warp] = count;
    }
    __syncthreads();

    if (threadIdx.x == 0) {
        mean = 0.0f;
        m2 = 0.0f;
        count = 0;

        #pragma unroll
        for (int i = 0; i < 8; ++i) {
            const int other_count = counts[i];
            if (other_count > 0) {
                const int total = count + other_count;
                const float delta = means[i] - mean;
                m2 += moments[i] +
                      delta * delta * (float(count) * other_count / total);
                mean += delta * (float(other_count) / total);
                count = total;
            }
        }

        partial[static_cast<int64_t>(group) * parts + part] =
            make_float2(mean, m2);
    }
}

__global__ void finish_stats_kernel(
    const float2* __restrict__ partial,
    float2* __restrict__ stats,
    int group_elements,
    int parts,
    float eps)
{
    const int group = blockIdx.x;
    const float2* base =
        partial + static_cast<int64_t>(group) * parts;

    float weighted_mean = 0.0f;
    for (int p = threadIdx.x; p < parts; p += THREADS) {
        const int count = min(CHUNK, group_elements - p * CHUNK);
        weighted_mean += base[p].x * count;
    }

    const float total_mean = block_sum(weighted_mean);
    __shared__ float mean;
    if (threadIdx.x == 0)
        mean = total_mean / group_elements;
    __syncthreads();

    float moment = 0.0f;
    for (int p = threadIdx.x; p < parts; p += THREADS) {
        const int count = min(CHUNK, group_elements - p * CHUNK);
        const float2 value = base[p];
        const float delta = value.x - mean;
        moment += value.y + delta * delta * count;
    }

    const float total_moment = block_sum(moment);
    if (threadIdx.x == 0) {
        stats[group] = make_float2(
            mean,
            rsqrtf(fmaxf(total_moment / group_elements, 0.0f) + eps));
    }
}

__global__ void activation_logsumexp_kernel(
    const float* __restrict__ x,
    const float* __restrict__ weight,
    const float* __restrict__ bias,
    const float2* __restrict__ stats,
    float* __restrict__ output,
    int spatial,
    int groups)
{
    const int pixel = blockIdx.x * 32 + threadIdx.x;
    const int batch = blockIdx.y;
    const int lane = threadIdx.x;
    const int channel_lane = threadIdx.y;
    const int channels_per_group = 64 / groups;

    float values[8];
    float local_max = -CUDART_INF_F;

    #pragma unroll
    for (int i = 0; i < 8; ++i) {
        const int channel = channel_lane + i * 8;
        float result = -CUDART_INF_F;

        if (pixel < spatial) {
            const float value =
                x[(static_cast<int64_t>(batch) * 64 + channel) *
                  spatial + pixel];

            const float2 stat =
                stats[batch * groups + channel / channels_per_group];

            const float normalized =
                (value - stat.x) * stat.y;
            const float activated =
                tanhf(normalized * weight[channel] + bias[channel]);

            // tanh lies in [-1, 1], so the Hardswish clamp is redundant.
            result = value + activated * (activated + 3.0f) / 6.0f;
        }

        values[i] = result;
        local_max = fmaxf(local_max, result);
    }

    float local_sum = 0.0f;
    if (pixel < spatial) {
        #pragma unroll
        for (int i = 0; i < 8; ++i)
            local_sum += __expf(values[i] - local_max);
    }

    __shared__ float maxima[8][32];
    __shared__ float sums[8][32];
    maxima[channel_lane][lane] = local_max;
    sums[channel_lane][lane] = local_sum;
    __syncthreads();

    if (channel_lane == 0 && pixel < spatial) {
        float maximum = maxima[0][lane];
        #pragma unroll
        for (int i = 1; i < 8; ++i)
            maximum = fmaxf(maximum, maxima[i][lane]);

        float sum = 0.0f;
        #pragma unroll
        for (int i = 0; i < 8; ++i)
            sum += sums[i][lane] * __expf(maxima[i][lane] - maximum);

        output[static_cast<int64_t>(batch) * spatial + pixel] =
            maximum + logf(sum);
    }
}

torch::Tensor fused_post_conv(
    torch::Tensor x,
    torch::Tensor weight,
    torch::Tensor bias,
    int64_t groups,
    double eps)
{
    const c10::cuda::CUDAGuard device_guard(x.device());

    const int batch = static_cast<int>(x.size(0));
    const int spatial = static_cast<int>(x.size(2) * x.size(3));
    const int group_elements = (64 / groups) * spatial;
    const int parts = (group_elements + CHUNK - 1) / CHUNK;

    auto partial = torch::empty({batch, groups, parts, 2}, x.options());
    auto stats = torch::empty({batch, groups, 2}, x.options());
    auto output = torch::empty(
        {batch, 1, x.size(2), x.size(3)}, x.options());

    if (batch == 0 || spatial == 0)
        return output;

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

    partial_stats_kernel<<<dim3(batch * groups, parts), THREADS, 0, stream>>>(
        x.data_ptr<float>(),
        reinterpret_cast<float2*>(partial.data_ptr<float>()),
        group_elements,
        parts);
    C10_CUDA_KERNEL_LAUNCH_CHECK();

    finish_stats_kernel<<<batch * groups, THREADS, 0, stream>>>(
        reinterpret_cast<const float2*>(partial.data_ptr<float>()),
        reinterpret_cast<float2*>(stats.data_ptr<float>()),
        group_elements,
        parts,
        static_cast<float>(eps));
    C10_CUDA_KERNEL_LAUNCH_CHECK();

    activation_logsumexp_kernel<<<
        dim3((spatial + 31) / 32, batch), dim3(32, 8), 0, stream>>>(
        x.data_ptr<float>(),
        weight.data_ptr<float>(),
        bias.data_ptr<float>(),
        reinterpret_cast<const float2*>(stats.data_ptr<float>()),
        output.data_ptr<float>(),
        spatial,
        static_cast<int>(groups));
    C10_CUDA_KERNEL_LAUNCH_CHECK();

    return output;
}
"""

_fused_ops = load_inline(
    name="fused_gn_tanh_hardswish_residual_lse_v1",
    cpp_sources=_cpp_source,
    cuda_sources=_cuda_source,
    functions=["fused_post_conv"],
    extra_cflags=["-O3"],
    extra_cuda_cflags=["-O3", "--fmad=false"],
    verbose=False,
)


class ModelNew(nn.Module):
    def __init__(self, in_channels, out_channels, kernel_size, groups, eps=1e-5):
        super().__init__()
        self.conv = nn.Conv2d(in_channels, out_channels, kernel_size)
        self.group_norm = nn.GroupNorm(groups, out_channels, eps=eps)
        self.tanh = nn.Tanh()
        self.hard_swish = nn.Hardswish()

    def forward(self, x):
        x_conv = self.conv(x)

        if (
            x_conv.is_cuda
            and x_conv.dtype == torch.float32
            and x_conv.shape[1] == 64
            and x_conv.is_contiguous()
            and self.group_norm.weight.dtype == torch.float32
            and self.group_norm.bias.dtype == torch.float32
        ):
            return _fused_ops.fused_post_conv(
                x_conv,
                self.group_norm.weight,
                self.group_norm.bias,
                self.group_norm.num_groups,
                self.group_norm.eps,
            )

        x_norm = self.group_norm(x_conv)
        x_tanh = self.tanh(x_norm)
        x_hard_swish = self.hard_swish(x_tanh)
        return torch.logsumexp(x_conv + x_hard_swish, dim=1, keepdim=True)
```