```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 pack_channels_last_cuda(torch::Tensor x);
"""

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

template <typename scalar_t>
__global__ void pack_channels_last_kernel(
    const scalar_t* __restrict__ input,
    scalar_t* __restrict__ output,
    const int64_t channels,
    const int64_t spatial,
    const int64_t spatial_tiles,
    const int64_t channel_tiles,
    const int64_t total_tiles)
{
    __shared__ scalar_t tile[32][33];

    const int tx = threadIdx.x;
    const int ty = threadIdx.y;

    for (int64_t tile_id = blockIdx.x;
         tile_id < total_tiles;
         tile_id += gridDim.x) {
        const int64_t spatial_tile = tile_id % spatial_tiles;
        const int64_t remaining = tile_id / spatial_tiles;
        const int64_t channel_tile = remaining % channel_tiles;
        const int64_t batch = remaining / channel_tiles;

        const int64_t spatial_base = spatial_tile * 32;
        const int64_t channel_base = channel_tile * 32;
        const int64_t batch_base = batch * channels * spatial;

        #pragma unroll
        for (int j = 0; j < 32; j += 8) {
            const int64_t c = channel_base + ty + j;
            const int64_t s = spatial_base + tx;
            if (c < channels && s < spatial) {
                tile[ty + j][tx] =
                    input[batch_base + c * spatial + s];
            }
        }

        __syncthreads();

        #pragma unroll
        for (int j = 0; j < 32; j += 8) {
            const int64_t c = channel_base + tx;
            const int64_t s = spatial_base + ty + j;
            if (c < channels && s < spatial) {
                output[batch_base + s * channels + c] =
                    tile[tx][ty + j];
            }
        }

        __syncthreads();
    }
}

torch::Tensor pack_channels_last_cuda(torch::Tensor x)
{
    TORCH_CHECK(x.is_cuda(), "Expected a CUDA tensor");
    TORCH_CHECK(x.dim() == 4, "Expected a four-dimensional tensor");
    TORCH_CHECK(x.is_contiguous(), "Expected contiguous NCHW storage");

    const c10::cuda::CUDAGuard device_guard(x.device());

    auto output = torch::empty(
        x.sizes(),
        x.options().memory_format(at::MemoryFormat::ChannelsLast));

    if (x.numel() == 0) {
        return output;
    }

    const int64_t channels = x.size(1);
    const int64_t spatial = x.size(2) * x.size(3);
    const int64_t spatial_tiles = (spatial + 31) / 32;
    const int64_t channel_tiles = (channels + 31) / 32;
    const int64_t total_tiles =
        x.size(0) * spatial_tiles * channel_tiles;

    const int blocks =
        static_cast<int>(total_tiles < 65535 ? total_tiles : 65535);
    const dim3 threads(32, 8);
    const auto stream =
        at::cuda::getCurrentCUDAStream(x.get_device());

    AT_DISPATCH_FLOATING_TYPES_AND2(
        at::ScalarType::Half,
        at::ScalarType::BFloat16,
        x.scalar_type(),
        "pack_channels_last_cuda",
        [&] {
            pack_channels_last_kernel<scalar_t>
                <<<blocks, threads, 0, stream>>>(
                    x.data_ptr<scalar_t>(),
                    output.data_ptr<scalar_t>(),
                    channels,
                    spatial,
                    spatial_tiles,
                    channel_tiles,
                    total_tiles);
        });

    C10_CUDA_KERNEL_LAUNCH_CHECK();
    return output;
}
"""

_layout_ops = load_inline(
    name="asymmetric_conv_channels_last_v1",
    cpp_sources=_cpp_source,
    cuda_sources=_cuda_source,
    functions=["pack_channels_last_cuda"],
    extra_cflags=["-O3"],
    extra_cuda_cflags=["-O3"],
    verbose=False,
)


class _PackChannelsLast(torch.autograd.Function):
    @staticmethod
    def forward(ctx, x):
        return _layout_ops.pack_channels_last_cuda(x)

    @staticmethod
    def backward(ctx, grad_output):
        return grad_output


class ModelNew(nn.Module):
    def __init__(
        self,
        in_channels: int,
        out_channels: int,
        kernel_size: tuple,
        stride: int = 1,
        padding: int = 0,
        dilation: int = 1,
        groups: int = 1,
        bias: bool = False,
    ):
        super().__init__()
        self.conv2d = nn.Conv2d(
            in_channels,
            out_channels,
            kernel_size,
            stride=stride,
            padding=padding,
            dilation=dilation,
            groups=groups,
            bias=bias,
        )
        self.conv2d.to(memory_format=torch.channels_last)

    def forward(self, x: torch.Tensor) -> torch.Tensor:
        if x.is_cuda and x.dim() == 4:
            if not x.is_contiguous(memory_format=torch.channels_last):
                if x.is_contiguous() and x.dtype in (
                    torch.float16,
                    torch.bfloat16,
                    torch.float32,
                    torch.float64,
                ):
                    x = _PackChannelsLast.apply(x)
                else:
                    x = x.contiguous(memory_format=torch.channels_last)

        return self.conv2d(x)
```