ot.emd / ot.emd2 crash the interpreter on empty inputs (out-of-bounds loop in setDenseCostMatrix) #901

Describe the bug

ot.emd and ot.emd2 crash the Python interpreter (SIGBUS on macOS) when called on an empty problem, i.e. with a cost matrix of shape (0, 0). No Python exception is raised, so the whole process dies.

The fault is in NetworkSimplexSimple::setDenseCostMatrix (ot/lp/network_simplex_simple.h#L926-L937), which EMD_wrap calls on the dense-cost fast path (ot/lp/EMD_wrapper.cpp#L273-L274):

_max_cost = D[0];
for (ArcsType i = 1; i != _arc_num; ++i) {
    if (D[i] > _max_cost) _max_cost = D[i];
}

With no arcs (_arc_num == 0), D[0] already reads out of bounds, and i != _arc_num never becomes false, so the loop keeps reading past the buffer until it faults. This code was added in #796 and is in the 0.9.7 / 0.9.7.post1 releases.

To Reproduce

Code sample

import numpy as np
import ot

M = np.zeros((0, 0))
ot.emd2([], [], M)  # interpreter crashes (SIGBUS on macOS), no Python exception

ot.emd([], [], M) crashes in the same way. Passing empty arrays for the weights explicitly (np.zeros(0)) does not change anything.

Stack trace of the faulting thread (from the macOS crash report):

emd_wrap.cpython-314-darwin.so  lemon::NetworkSimplexSimple<lemon::FullBipartiteDigraph, double, double, unsigned int, long long>::setDenseCostMatrix(double const*, int)
emd_wrap.cpython-314-darwin.so  EMD_wrap(int, int, double*, double*, double*, double*, double*, double*, double*, unsigned long long, double*, double*)
emd_wrap.cpython-314-darwin.so  __pyx_pf_2ot_2lp_8emd_wrap_2emd_c(...)
Python                          PyObject_Vectorcall
...
exception: EXC_BAD_ACCESS / SIGBUS, KERN_PROTECTION_FAILURE

Expected behavior

A Python-level error, or a well-defined result for the empty problem, but never an interpreter crash. Possible fixes include:

• an early return in EMD_wrap (and possibly in the other EMD_wrap_* entry points) when there are no source or target nodes,
• a guard in setDenseCostMatrix for _arc_num == 0,
• a check in emd_c or ot.emd / ot.emd2 that rejects empty inputs with a ValueError.

For comparison, other solvers on empty inputs currently fail with unrelated Python exceptions: ot.emd_1d raises ValueError: negative dimensions are not allowed, ot.wasserstein_1d and ot.sinkhorn2 raise ZeroDivisionError, and ot.lp.emd2_lazy raises an unrelated backend ValueError. None of them handles an empty problem explicitly.

Environment (please complete the following information):

• OS (e.g. MacOS, Windows, Linux): macOS 27.0.1 (arm64)
• Python version: 3.14.8
• How was POT installed (source, pip, conda): source, editable install of the [WIP] Add exact EMD-L1 solver for Cartesian grid histograms #863 branch, which includes master at 5cb82c9. ot/lp/network_simplex_simple.h and the dense EMD_wrap path are identical to master and to 0.9.7.post1.
• Build command you used (if compiling from source): uv pip install -e .

Output of the following code snippet:

macOS-27.0.1-arm64-arm-64bit-Mach-O
Python 3.14.8 (main, Sep 30 2026, 17:55:09) [Clang 21.0.0 (clang-2100.3.34.2)]
NumPy 2.5.3
SciPy 1.18.1
POT 0.9.8.dev

Additional context

Found while addressing a Copilot review comment on #863 about empty and invalid grid shapes reaching the native grid solver. While checking how other POT solvers handle an empty problem, ot.emd turned out to crash. #863 does not touch this code path.