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

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

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

__device__ __forceinline__ float4 streaming_load(const float4* ptr) {
    float4 v;
    asm volatile(
        "ld.global.cs.v4.f32 {%0, %1, %2, %3}, [%4];"
        : "=f"(v.x), "=f"(v.y), "=f"(v.z), "=f"(v.w)
        : "l"(ptr)
    );
    return v;
}

__global__ __launch_bounds__(256, 4)
void triplet_rows(
    const float4* __restrict__ anchor,
    const float4* __restrict__ positive,
    const float4* __restrict__ negative,
    float* __restrict__ partials,
    int64_t rows,
    float margin
) {
    constexpr int ThreadsPerRow = 64;
    constexpr int RowsPerBlock = 4;
    constexpr int VectorsPerRow = 2048;

    const int tid = threadIdx.x;
    const int lane = tid & 31;
    const int warp = tid >> 5;
    const int row_thread = tid & 63;
    const int64_t row = int64_t(blockIdx.x) * RowsPerBlock + tid / ThreadsPerRow;

    float p0 = 0.f, p1 = 0.f, p2 = 0.f, p3 = 0.f;
    float n0 = 0.f, n1 = 0.f, n2 = 0.f, n3 = 0.f;

    if (row < rows) {
        const int64_t base = row * VectorsPerRow;

        #pragma unroll 4
        for (int j = row_thread; j < VectorsPerRow; j += ThreadsPerRow) {
            const float4 av = streaming_load(anchor + base + j);
            const float4 pv = streaming_load(positive + base + j);
            const float4 nv = streaming_load(negative + base + j);

            const float dp0 = (av.x - pv.x) + 1.e-6f;
            const float dp1 = (av.y - pv.y) + 1.e-6f;
            const float dp2 = (av.z - pv.z) + 1.e-6f;
            const float dp3 = (av.w - pv.w) + 1.e-6f;

            const float dn0 = (av.x - nv.x) + 1.e-6f;
            const float dn1 = (av.y - nv.y) + 1.e-6f;
            const float dn2 = (av.z - nv.z) + 1.e-6f;
            const float dn3 = (av.w - nv.w) + 1.e-6f;

            p0 = fmaf(dp0, dp0, p0);
            p1 = fmaf(dp1, dp1, p1);
            p2 = fmaf(dp2, dp2, p2);
            p3 = fmaf(dp3, dp3, p3);

            n0 = fmaf(dn0, dn0, n0);
            n1 = fmaf(dn1, dn1, n1);
            n2 = fmaf(dn2, dn2, n2);
            n3 = fmaf(dn3, dn3, n3);
        }
    }

    float ps = warp_sum((p0 + p1) + (p2 + p3));
    float ns = warp_sum((n0 + n1) + (n2 + n3));

    __shared__ float positive_sums[8];
    __shared__ float negative_sums[8];

    if (lane == 0) {
        positive_sums[warp] = ps;
        negative_sums[warp] = ns;
    }
    __syncthreads();

    if (warp == 0) {
        ps = lane < 8 ? positive_sums[lane] : 0.f;
        ns = lane < 8 ? negative_sums[lane] : 0.f;

        ps += __shfl_down_sync(0xffffffffu, ps, 1, 2);
        ns += __shfl_down_sync(0xffffffffu, ns, 1, 2);

        float loss = 0.f;
        if (lane < 8 && (lane & 1) == 0 &&
            int64_t(blockIdx.x) * RowsPerBlock + lane / 2 < rows) {
            loss = (margin + sqrtf(ps)) - sqrtf(ns);
            loss = loss < 0.f ? 0.f : loss;
        }

        loss = warp_sum(loss);
        if (lane == 0)
            partials[blockIdx.x] = loss;
    }
}

__global__ __launch_bounds__(256)
void triplet_mean(
    const float* __restrict__ partials,
    float* __restrict__ output,
    int64_t count,
    int64_t rows
) {
    float s0 = 0.f, s1 = 0.f, s2 = 0.f, s3 = 0.f;
    const float4* vectors = reinterpret_cast<const float4*>(partials);

    for (int64_t i = threadIdx.x; i < count / 4; i += 256) {
        const float4 v = vectors[i];
        s0 += v.x;
        s1 += v.y;
        s2 += v.z;
        s3 += v.w;
    }
    for (int64_t i = (count / 4) * 4 + threadIdx.x; i < count; i += 256)
        s0 += partials[i];

    float sum = warp_sum((s0 + s1) + (s2 + s3));
    __shared__ float warp_sums[8];

    const int lane = threadIdx.x & 31;
    const int warp = threadIdx.x >> 5;

    if (lane == 0)
        warp_sums[warp] = sum;
    __syncthreads();

    if (warp == 0) {
        sum = lane < 8 ? warp_sums[lane] : 0.f;
        sum = warp_sum(sum);
        if (lane == 0)
            output[0] = sum / static_cast<float>(rows);
    }
}

torch::Tensor triplet_forward(
    torch::Tensor anchor,
    torch::Tensor positive,
    torch::Tensor negative,
    double margin
) {
    const c10::cuda::CUDAGuard guard(anchor.device());
    const int64_t rows = anchor.size(0);
    const int64_t groups = (rows + 3) / 4;

    auto partials = torch::empty({groups}, anchor.options());
    auto output = torch::empty({}, anchor.options());
    cudaStream_t stream = at::cuda::getCurrentCUDAStream();

    triplet_rows<<<static_cast<unsigned int>(groups), 256, 0, stream>>>(
        reinterpret_cast<const float4*>(anchor.data_ptr<float>()),
        reinterpret_cast<const float4*>(positive.data_ptr<float>()),
        reinterpret_cast<const float4*>(negative.data_ptr<float>()),
        partials.data_ptr<float>(),
        rows,
        static_cast<float>(margin)
    );
    C10_CUDA_KERNEL_LAUNCH_CHECK();

    triplet_mean<<<1, 256, 0, stream>>>(
        partials.data_ptr<float>(),
        output.data_ptr<float>(),
        groups,
        rows
    );
    C10_CUDA_KERNEL_LAUNCH_CHECK();

    return output;
}
"""

_triplet_cuda = load_inline(
    name="fused_triplet_margin_streaming_v5",
    cpp_sources=r"""
        torch::Tensor triplet_forward(
            torch::Tensor anchor,
            torch::Tensor positive,
            torch::Tensor negative,
            double margin);
    """,
    cuda_sources=_cuda_source,
    functions=["triplet_forward"],
    extra_cflags=["-O3"],
    extra_cuda_cflags=["-O3"],
    verbose=False,
)


class ModelNew(nn.Module):
    def __init__(self, margin=1.0):
        super().__init__()
        if margin <= 0:
            raise ValueError("margin must be greater than zero")
        self.margin = float(margin)

    def forward(self, anchor, positive, negative):
        tensors = (anchor, positive, negative)
        use_custom = (
            anchor.ndim == 2
            and anchor.shape[1] == 8192
            and 0 < anchor.shape[0] <= 8589934588
            and anchor.shape == positive.shape == negative.shape
            and anchor.device == positive.device == negative.device
            and all(
                t.is_cuda
                and t.dtype == torch.float32
                and t.is_contiguous()
                and t.data_ptr() % 16 == 0
                for t in tensors
            )
            and not (
                torch.is_grad_enabled()
                and any(t.requires_grad for t in tensors)
            )
        )
        if use_custom:
            return _triplet_cuda.triplet_forward(
                anchor, positive, negative, self.margin
            )
        return F.triplet_margin_loss(
            anchor, positive, negative, margin=self.margin
        )
```