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


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

torch::Tensor maxpool3d_fast_cuda(torch::Tensor x);
"""

_cuda_source = r"""
#include <torch/extension.h>
#include <c10/cuda/CUDAGuard.h>
#include <c10/cuda/CUDAStream.h>
#include <c10/cuda/CUDAException.h>
#include <cuda_runtime.h>
#include <stdint.h>

__global__ void maxpool3d_128_kernel(
    const float* __restrict__ input,
    float* __restrict__ output,
    const int64_t rows)
{
    // Each pair of warps processes one complete output row.
    const int ow = threadIdx.x;
    const int64_t row = int64_t(blockIdx.x) * 4 + threadIdx.y;

    if (ow >= 62 || row >= rows)
        return;

    const int oh = row % 62;
    const int64_t depth_row = row / 62;
    const int od = depth_row % 62;
    const int64_t plane = depth_row / 62;

    const int z0 = od * 2 - 1;
    const int y0 = oh * 2 - 1;
    const int x0 = ow * 2 - 1;

    const float* base = input + plane * int64_t(128 * 128 * 128);
    float best = -CUDART_INF_F;

    // The largest sampled coordinate is 127, so only lower bounds
    // need checking for this specialization.
    #pragma unroll
    for (int kz = 0; kz < 3; ++kz) {
        const int z = z0 + kz * 3;
        if (z >= 0) {
            #pragma unroll
            for (int ky = 0; ky < 3; ++ky) {
                const int y = y0 + ky * 3;
                if (y >= 0) {
                    const int offset = (z * 128 + y) * 128;

                    #pragma unroll
                    for (int kx = 0; kx < 3; ++kx) {
                        const int w = x0 + kx * 3;
                        if (w >= 0) {
                            const float value = base[offset + w];
                            // Match PyTorch's NaN and equal-value behavior.
                            if (value > best || isnan(value))
                                best = value;
                        }
                    }
                }
            }
        }
    }

    output[row * 62 + ow] = best;
}

torch::Tensor maxpool3d_fast_cuda(torch::Tensor x)
{
    TORCH_CHECK(x.is_cuda(), "Expected a CUDA tensor");
    TORCH_CHECK(x.scalar_type() == at::kFloat,
                "Expected float32 input");
    TORCH_CHECK(x.is_contiguous(), "Expected contiguous input");
    TORCH_CHECK(x.dim() == 4 || x.dim() == 5,
                "Expected a 4D or 5D input");

    const int rank = x.dim();
    TORCH_CHECK(
        x.size(rank - 3) == 128 &&
        x.size(rank - 2) == 128 &&
        x.size(rank - 1) == 128,
        "Expected spatial dimensions 128 x 128 x 128");

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

    auto shape = x.sizes().vec();
    shape[rank - 3] = 62;
    shape[rank - 2] = 62;
    shape[rank - 1] = 62;
    auto output = torch::empty(shape, x.options());

    const int64_t rows = output.numel() / 62;
    if (rows == 0)
        return output;

    const dim3 threads(64, 4);
    const dim3 blocks(static_cast<unsigned int>((rows + 3) / 4));
    const cudaStream_t stream =
        c10::cuda::getCurrentCUDAStream(x.get_device());

    maxpool3d_128_kernel<<<blocks, threads, 0, stream>>>(
        x.data_ptr<float>(), output.data_ptr<float>(), rows);
    C10_CUDA_KERNEL_LAUNCH_CHECK();

    return output;
}
"""

_extension = load_inline(
    name="maxpool3d_dilated_128_v2",
    cpp_sources=_cpp_source,
    cuda_sources=_cuda_source,
    functions=["maxpool3d_fast_cuda"],
    extra_cflags=["-O3"],
    extra_cuda_cflags=["-O3"],
    verbose=False,
)


class ModelNew(nn.Module):
    def __init__(
        self,
        kernel_size: int,
        stride: int = None,
        padding: int = 0,
        dilation: int = 1,
        return_indices: bool = False,
        ceil_mode: bool = False,
    ):
        super().__init__()

        self.maxpool = nn.MaxPool3d(
            kernel_size=kernel_size,
            stride=stride,
            padding=padding,
            dilation=dilation,
            return_indices=return_indices,
            ceil_mode=ceil_mode,
        )

        self._specialized = (
            tuple(_triple(kernel_size)) == (3, 3, 3)
            and tuple(_triple(kernel_size if stride is None else stride))
            == (2, 2, 2)
            and tuple(_triple(padding)) == (1, 1, 1)
            and tuple(_triple(dilation)) == (3, 3, 3)
            and not return_indices
            and not ceil_mode
        )

    def forward(self, x: torch.Tensor):
        if (
            self._specialized
            and x.is_cuda
            and x.dtype == torch.float32
            and x.ndim in (4, 5)
            and x.shape[-3:] == (128, 128, 128)
            and x.shape[-4] > 0
            and x.is_contiguous()
            and not (torch.is_grad_enabled() and x.requires_grad)
        ):
            return _extension.maxpool3d_fast_cuda(x)

        return self.maxpool(x)
```