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


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

std::vector<torch::Tensor> pack_polyphase_cuda(
    torch::Tensor x, torch::Tensor weight);

torch::Tensor polyphase_pool_sum_cuda(
    torch::Tensor x, torch::Tensor bias);
"""

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

__global__ void pack_input_kernel(
    const float* __restrict__ x,
    float* __restrict__ out,
    int64_t count,
    int channels,
    int depth, int height, int width,
    int pd, int ph, int pw,
    int64_t sn, int64_t sc,
    int64_t sd, int64_t sh, int64_t sw
) {
    for (int64_t index = int64_t(blockIdx.x) * blockDim.x + threadIdx.x;
         index < count;
         index += int64_t(blockDim.x) * gridDim.x) {
        int64_t t = index;
        const int c = t % channels;
        t /= channels;
        const int w = t % pw - 1;
        t /= pw;
        const int h = t % ph - 1;
        t /= ph;
        const int d = t % pd - 1;
        const int n = t / pd;

        float value = 0.0f;
        if (d >= 0 && d < depth &&
            h >= 0 && h < height &&
            w >= 0 && w < width) {
            value = x[int64_t(n) * sn + int64_t(c) * sc +
                      int64_t(d) * sd + int64_t(h) * sh +
                      int64_t(w) * sw];
        }
        out[index] = value;
    }
}

__global__ void pack_weight_kernel(
    const float* __restrict__ weight,
    float* __restrict__ packed,
    int64_t count,
    int in_channels,
    int out_channels,
    int64_t si, int64_t so,
    int64_t sd, int64_t sh, int64_t sw
) {
    for (int64_t index = int64_t(blockIdx.x) * blockDim.x + threadIdx.x;
         index < count;
         index += int64_t(blockDim.x) * gridDim.x) {
        int64_t t = index;
        const int ic = t % in_channels;
        t /= in_channels;
        const int rw = t % 3;
        t /= 3;
        const int rh = t % 3;
        t /= 3;
        const int rd = t % 3;
        t /= 3;
        const int oc = t % out_channels;
        const int phase = t / out_channels;

        const int kd = ((phase >> 2) & 1) + 4 - 2 * rd;
        const int kh = ((phase >> 1) & 1) + 4 - 2 * rh;
        const int kw = (phase & 1) + 4 - 2 * rw;

        float value = 0.0f;
        if (kd < 5 && kh < 5 && kw < 5) {
            value = weight[int64_t(ic) * si + int64_t(oc) * so +
                           int64_t(kd) * sd + int64_t(kh) * sh +
                           int64_t(kw) * sw];
        }
        packed[index] = value;
    }
}

std::vector<torch::Tensor> pack_polyphase_cuda(
    torch::Tensor x, torch::Tensor weight
) {
    const c10::cuda::CUDAGuard guard(x.device());

    TORCH_CHECK(x.scalar_type() == at::kFloat &&
                weight.scalar_type() == at::kFloat,
                "Packing requires float32 tensors");
    TORCH_CHECK(x.device() == weight.device(), "Device mismatch");

    const int n = x.size(0);
    const int ic = x.size(1);
    const int oc = weight.size(1);
    const int d = x.size(2);
    const int h = x.size(3);
    const int w = x.size(4);

    const int pd = ((2 * d - 1) / 6) * 3 + 2;
    const int ph = ((2 * h - 1) / 6) * 3 + 2;
    const int pw = ((2 * w - 1) / 6) * 3 + 2;

    auto packed_x = torch::empty({n, pd, ph, pw, ic}, x.options())
                        .permute({0, 4, 1, 2, 3});
    auto packed_weight = torch::empty(
        {8 * oc, 3, 3, 3, ic}, weight.options()
    ).permute({0, 4, 1, 2, 3});

    const auto stream = at::cuda::getCurrentCUDAStream();
    const int threads = 256;
    const int input_blocks = static_cast<int>(
        std::min<int64_t>((packed_x.numel() + threads - 1) / threads, 65535)
    );
    const int weight_blocks = static_cast<int>(
        std::min<int64_t>((packed_weight.numel() + threads - 1) / threads, 65535)
    );

    if (packed_x.numel() != 0) {
        pack_input_kernel<<<input_blocks, threads, 0, stream>>>(
            x.data_ptr<float>(), packed_x.data_ptr<float>(),
            packed_x.numel(), ic, d, h, w, pd, ph, pw,
            x.stride(0), x.stride(1), x.stride(2),
            x.stride(3), x.stride(4)
        );
    }
    pack_weight_kernel<<<weight_blocks, threads, 0, stream>>>(
        weight.data_ptr<float>(), packed_weight.data_ptr<float>(),
        packed_weight.numel(), ic, oc,
        weight.stride(0), weight.stride(1), weight.stride(2),
        weight.stride(3), weight.stride(4)
    );
    C10_CUDA_KERNEL_LAUNCH_CHECK();
    return {packed_x, packed_weight};
}

template <typename T>
__device__ __forceinline__ T max_nan(T a, T b) {
    return isnan(a) ? a : (isnan(b) ? b : (a > b ? a : b));
}

template <typename scalar_t>
__global__ void polyphase_pool_sum_kernel(
    const scalar_t* __restrict__ x,
    const scalar_t* __restrict__ bias,
    scalar_t* __restrict__ out,
    int channels,
    int od, int oh, int ow,
    int64_t sn, int64_t sc,
    int64_t sd, int64_t sh, int64_t sw
) {
    using acc_t = at::acc_type<scalar_t, true>;
    __shared__ acc_t maxima[4][64];
    __shared__ acc_t warp_sums[2];

    const int tid = threadIdx.x;
    const int lane = tid & 63;
    const int group = tid >> 6;
    const int index = blockIdx.x;
    const int spatial_size = od * oh * ow;
    const int n = index / spatial_size;

    int spatial = index % spatial_size;
    const int w = spatial % ow;
    spatial /= ow;
    const int h = spatial % oh;
    const int d = spatial / oh;

    const int64_t origin =
        int64_t(n) * sn + int64_t(3 * d) * sd +
        int64_t(3 * h) * sh + int64_t(3 * w) * sw;

    acc_t sum = acc_t(0);

    for (int first = 0; first < channels; first += 64) {
        const int c = first + lane;
        acc_t maximum = -static_cast<acc_t>(INFINITY);

        if (c < channels) {
            const scalar_t* base =
                x + origin + int64_t(c + group * channels) * sc;
            const int64_t other_phase = int64_t(4 * channels) * sc;

            #pragma unroll
            for (int z = 0; z < 3; ++z) {
                #pragma unroll
                for (int y = 0; y < 3; ++y) {
                    #pragma unroll
                    for (int xx = 0; xx < 3; ++xx) {
                        const int64_t offset =
                            int64_t(z) * sd + int64_t(y) * sh +
                            int64_t(xx) * sw;
                        maximum = max_nan(
                            maximum, static_cast<acc_t>(base[offset])
                        );
                        maximum = max_nan(
                            maximum,
                            static_cast<acc_t>(base[offset + other_phase])
                        );
                    }
                }
            }
        }

        maxima[group][lane] = maximum;
        __syncthreads();

        if (group == 0 && c < channels) {
            const acc_t pooled = max_nan(
                max_nan(maxima[0][lane], maxima[1][lane]),
                max_nan(maxima[2][lane], maxima[3][lane])
            );
            const scalar_t biased = static_cast<scalar_t>(
                pooled + static_cast<acc_t>(bias[c])
            );
            sum += static_cast<acc_t>(biased);
        }
        __syncthreads();
    }

    if (tid < 64) {
        #pragma unroll
        for (int offset = 16; offset > 0; offset >>= 1) {
            sum += __shfl_down_sync(0xffffffff, sum, offset);
        }
        if ((tid & 31) == 0) {
            warp_sums[tid >> 5] = sum;
        }
    }
    __syncthreads();

    if (tid == 0) {
        out[index] = static_cast<scalar_t>(warp_sums[0] + warp_sums[1]);
    }
}

torch::Tensor polyphase_pool_sum_cuda(
    torch::Tensor x, torch::Tensor bias
) {
    const c10::cuda::CUDAGuard guard(x.device());
    const int channels = bias.numel();
    const int od = x.size(2) / 3;
    const int oh = x.size(3) / 3;
    const int ow = x.size(4) / 3;

    TORCH_CHECK(x.size(1) == 8 * channels, "Invalid phase channels");
    TORCH_CHECK(x.device() == bias.device(), "Device mismatch");

    bias = bias.to(x.scalar_type()).contiguous();
    auto output = torch::empty(
        {x.size(0), 1, od, oh, ow}, x.options()
    );
    if (output.numel() == 0) {
        return output;
    }

    const auto stream = at::cuda::getCurrentCUDAStream();
    AT_DISPATCH_FLOATING_TYPES_AND2(
        at::ScalarType::Half, at::ScalarType::BFloat16,
        x.scalar_type(), "polyphase_pool_sum_cuda", [&] {
            polyphase_pool_sum_kernel<scalar_t>
                <<<static_cast<unsigned int>(output.numel()), 256, 0, stream>>>(
                    x.data_ptr<scalar_t>(),
                    bias.data_ptr<scalar_t>(),
                    output.data_ptr<scalar_t>(),
                    channels, od, oh, ow,
                    x.stride(0), x.stride(1), x.stride(2),
                    x.stride(3), x.stride(4)
                );
        }
    );
    C10_CUDA_KERNEL_LAUNCH_CHECK();
    return output;
}
"""

_ops = load_inline(
    name="polyphase_conv3d_pool_sum_v3",
    cpp_sources=_cpp_source,
    cuda_sources=_cuda_source,
    functions=["pack_polyphase_cuda", "polyphase_pool_sum_cuda"],
    extra_cflags=["-O3"],
    extra_cuda_cflags=["-O3"],
    verbose=False,
)


class ModelNew(nn.Module):
    def __init__(self, in_channels, out_channels, kernel_size, stride, padding):
        super().__init__()
        self.conv_transpose = nn.ConvTranspose3d(
            in_channels,
            out_channels,
            kernel_size,
            stride=stride,
            padding=padding,
        )
        self.max_pool1 = nn.MaxPool3d(kernel_size=2)
        self.max_pool2 = nn.MaxPool3d(kernel_size=3)

    def forward(self, x):
        conv = self.conv_transpose
        use_polyphase = (
            x.is_cuda
            and x.dtype == torch.float32
            and conv.weight.dtype == torch.float32
            and x.ndim == 5
            and x.size(0) > 0
            and min(x.shape[2:]) >= 4
            and conv.kernel_size == (5, 5, 5)
            and conv.stride == (2, 2, 2)
            and conv.padding == (2, 2, 2)
            and conv.output_padding == (0, 0, 0)
            and conv.dilation == (1, 1, 1)
            and conv.groups == 1
            and conv.bias is not None
        )

        if not use_polyphase:
            y = conv(x)
            return self.max_pool2(self.max_pool1(y)).sum(dim=1, keepdim=True)

        packed_x, packed_weight = _ops.pack_polyphase_cuda(x, conv.weight)
        phases = F.conv3d(packed_x, packed_weight, bias=None)
        return _ops.polyphase_pool_sum_cuda(phases, conv.bias)
```