From b001d7a15676ead5d9b38d57609fd390c3749f81 Mon Sep 17 00:00:00 2001 From: Jay Stanley Date: Fri, 13 May 2022 11:52:55 -0500 Subject: remove calls to cerr (#377) typo in releases --- ot/lp/network_simplex_simple_omp.h | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) (limited to 'ot/lp/network_simplex_simple_omp.h') diff --git a/ot/lp/network_simplex_simple_omp.h b/ot/lp/network_simplex_simple_omp.h index 87e4c05..dde84fd 100644 --- a/ot/lp/network_simplex_simple_omp.h +++ b/ot/lp/network_simplex_simple_omp.h @@ -1610,9 +1610,7 @@ namespace lemon_omp { } else { - char errMess[1000]; - sprintf( errMess, "RESULT MIGHT BE INACURATE\nMax number of iteration reached, currently \%d. Sometimes iterations go on in cycle even though the solution has been reached, to check if it's the case here have a look at the minimal reduced cost. If it is very close to machine precision, you might actually have the correct solution, if not try setting the maximum number of iterations a bit higher\n",iter_number ); - std::cerr << errMess; + // max iters retVal = MAX_ITER_REACHED; break; } -- cgit v1.2.3 From 1f307594244dd4c274b64d028823cbcfff302f37 Mon Sep 17 00:00:00 2001 From: Nathan Cassereau <84033440+ncassereau-idris@users.noreply.github.com> Date: Wed, 1 Jun 2022 08:52:47 +0200 Subject: [MRG] numItermax in 64 bits in EMD solver (#380) * Correct test_mm_convergence for cupy * Fix bug where number of iterations is limited to 2^31 * Update RELEASES.md * Replace size_t with long long * Use uint64_t instead of long long --- RELEASES.md | 5 ++++- ot/lp/EMD.h | 5 +++-- ot/lp/EMD_wrapper.cpp | 4 ++-- ot/lp/emd_wrap.pyx | 9 +++++---- ot/lp/network_simplex_simple.h | 8 ++++---- ot/lp/network_simplex_simple_omp.h | 6 +++--- test/test_unbalanced.py | 38 ++++++++++++++++++++------------------ 7 files changed, 41 insertions(+), 34 deletions(-) (limited to 'ot/lp/network_simplex_simple_omp.h') diff --git a/RELEASES.md b/RELEASES.md index e761e64..fdaff59 100644 --- a/RELEASES.md +++ b/RELEASES.md @@ -10,7 +10,10 @@ - Fixed an issue where we could not ask TorchBackend to place a random tensor on GPU (Issue #371, PR #373) -- Fixed an issue where hitting iteration limits would be reported to stderr by std::cerr regardless of Python's stderr stream status. +- Fixed an issue where Sinkhorn solver assumed a symmetric cost matrix (Issue #374, PR #375) +- Fixed an issue where hitting iteration limits would be reported to stderr by std::cerr regardless of Python's stderr stream status (PR #377) +- Fixed an issue where the metric argument in ot.dist did not allow a callable parameter (Issue #378, PR #379) +- Fixed an issue where the max number of iterations in ot.emd was not allow to go beyond 2^31 (PR #380) ## 0.8.2 diff --git a/ot/lp/EMD.h b/ot/lp/EMD.h index 8a1f9ac..b56f060 100644 --- a/ot/lp/EMD.h +++ b/ot/lp/EMD.h @@ -18,6 +18,7 @@ #include #include +#include typedef unsigned int node_id_type; @@ -28,8 +29,8 @@ enum ProblemType { MAX_ITER_REACHED }; -int EMD_wrap(int n1,int n2, double *X, double *Y,double *D, double *G, double* alpha, double* beta, double *cost, int maxIter); -int EMD_wrap_omp(int n1,int n2, double *X, double *Y,double *D, double *G, double* alpha, double* beta, double *cost, int maxIter, int numThreads); +int EMD_wrap(int n1,int n2, double *X, double *Y,double *D, double *G, double* alpha, double* beta, double *cost, uint64_t maxIter); +int EMD_wrap_omp(int n1,int n2, double *X, double *Y,double *D, double *G, double* alpha, double* beta, double *cost, uint64_t maxIter, int numThreads); diff --git a/ot/lp/EMD_wrapper.cpp b/ot/lp/EMD_wrapper.cpp index 2bdc172..457216b 100644 --- a/ot/lp/EMD_wrapper.cpp +++ b/ot/lp/EMD_wrapper.cpp @@ -20,7 +20,7 @@ int EMD_wrap(int n1, int n2, double *X, double *Y, double *D, double *G, - double* alpha, double* beta, double *cost, int maxIter) { + double* alpha, double* beta, double *cost, uint64_t maxIter) { // beware M and C are stored in row major C style!!! using namespace lemon; @@ -122,7 +122,7 @@ int EMD_wrap(int n1, int n2, double *X, double *Y, double *D, double *G, int EMD_wrap_omp(int n1, int n2, double *X, double *Y, double *D, double *G, - double* alpha, double* beta, double *cost, int maxIter, int numThreads) { + double* alpha, double* beta, double *cost, uint64_t maxIter, int numThreads) { // beware M and C are stored in row major C style!!! using namespace lemon_omp; diff --git a/ot/lp/emd_wrap.pyx b/ot/lp/emd_wrap.pyx index 42e08f4..e5cec89 100644 --- a/ot/lp/emd_wrap.pyx +++ b/ot/lp/emd_wrap.pyx @@ -14,13 +14,14 @@ from ..utils import dist cimport cython cimport libc.math as math +from libc.stdint cimport uint64_t import warnings cdef extern from "EMD.h": - int EMD_wrap(int n1,int n2, double *X, double *Y,double *D, double *G, double* alpha, double* beta, double *cost, int maxIter) nogil - int EMD_wrap_omp(int n1,int n2, double *X, double *Y,double *D, double *G, double* alpha, double* beta, double *cost, int maxIter, int numThreads) nogil + int EMD_wrap(int n1,int n2, double *X, double *Y,double *D, double *G, double* alpha, double* beta, double *cost, uint64_t maxIter) nogil + int EMD_wrap_omp(int n1,int n2, double *X, double *Y,double *D, double *G, double* alpha, double* beta, double *cost, uint64_t maxIter, int numThreads) nogil cdef enum ProblemType: INFEASIBLE, OPTIMAL, UNBOUNDED, MAX_ITER_REACHED @@ -39,7 +40,7 @@ def check_result(result_code): @cython.boundscheck(False) @cython.wraparound(False) -def emd_c(np.ndarray[double, ndim=1, mode="c"] a, np.ndarray[double, ndim=1, mode="c"] b, np.ndarray[double, ndim=2, mode="c"] M, int max_iter, int numThreads): +def emd_c(np.ndarray[double, ndim=1, mode="c"] a, np.ndarray[double, ndim=1, mode="c"] b, np.ndarray[double, ndim=2, mode="c"] M, uint64_t max_iter, int numThreads): """ Solves the Earth Movers distance problem and returns the optimal transport matrix @@ -75,7 +76,7 @@ def emd_c(np.ndarray[double, ndim=1, mode="c"] a, np.ndarray[double, ndim=1, mod target histogram M : (ns,nt) numpy.ndarray, float64 loss matrix - max_iter : int + max_iter : uint64_t The maximum number of iterations before stopping the optimization algorithm if it has not converged. diff --git a/ot/lp/network_simplex_simple.h b/ot/lp/network_simplex_simple.h index 5b1038f..9612a8a 100644 --- a/ot/lp/network_simplex_simple.h +++ b/ot/lp/network_simplex_simple.h @@ -233,7 +233,7 @@ namespace lemon { /// mixed order in the internal data structure. /// In special cases, it could lead to better overall performance, /// but it is usually slower. Therefore it is disabled by default. - NetworkSimplexSimple(const GR& graph, bool arc_mixing, int nbnodes, ArcsType nb_arcs, size_t maxiters) : + NetworkSimplexSimple(const GR& graph, bool arc_mixing, int nbnodes, ArcsType nb_arcs, uint64_t maxiters) : _graph(graph), //_arc_id(graph), _arc_mixing(arc_mixing), _init_nb_nodes(nbnodes), _init_nb_arcs(nb_arcs), MAX(std::numeric_limits::max()), @@ -242,7 +242,7 @@ namespace lemon { { // Reset data structures reset(); - max_iter=maxiters; + max_iter = maxiters; } /// The type of the flow amounts, capacity bounds and supply values @@ -293,7 +293,7 @@ namespace lemon { private: - size_t max_iter; + uint64_t max_iter; TEMPLATE_DIGRAPH_TYPEDEFS(GR); typedef std::vector IntVector; @@ -1427,7 +1427,7 @@ namespace lemon { // Perform heuristic initial pivots if (!initialPivots()) return UNBOUNDED; - size_t iter_number=0; + uint64_t iter_number = 0; //pivot.setDantzig(true); // Execute the Network Simplex algorithm while (pivot.findEnteringArc()) { diff --git a/ot/lp/network_simplex_simple_omp.h b/ot/lp/network_simplex_simple_omp.h index dde84fd..5f19d73 100644 --- a/ot/lp/network_simplex_simple_omp.h +++ b/ot/lp/network_simplex_simple_omp.h @@ -244,7 +244,7 @@ namespace lemon_omp { /// mixed order in the internal data structure. /// In special cases, it could lead to better overall performance, /// but it is usually slower. Therefore it is disabled by default. - NetworkSimplexSimple(const GR& graph, bool arc_mixing, int nbnodes, ArcsType nb_arcs, size_t maxiters = 0, int numThreads=-1) : + NetworkSimplexSimple(const GR& graph, bool arc_mixing, int nbnodes, ArcsType nb_arcs, uint64_t maxiters = 0, int numThreads=-1) : _graph(graph), //_arc_id(graph), _arc_mixing(arc_mixing), _init_nb_nodes(nbnodes), _init_nb_arcs(nb_arcs), MAX(std::numeric_limits::max()), @@ -317,7 +317,7 @@ namespace lemon_omp { private: - size_t max_iter; + uint64_t max_iter; int num_threads; TEMPLATE_DIGRAPH_TYPEDEFS(GR); @@ -1563,7 +1563,7 @@ namespace lemon_omp { // Perform heuristic initial pivots if (!initialPivots()) return UNBOUNDED; - size_t iter_number = 0; + uint64_t iter_number = 0; // Execute the Network Simplex algorithm while (pivot.findEnteringArc()) { if ((++iter_number <= max_iter&&max_iter > 0) || max_iter<=0) { diff --git a/test/test_unbalanced.py b/test/test_unbalanced.py index 02b3fc3..fc40df0 100644 --- a/test/test_unbalanced.py +++ b/test/test_unbalanced.py @@ -295,26 +295,27 @@ def test_mm_convergence(nx): x = rng.randn(n, 2) rng = np.random.RandomState(75) y = rng.randn(n, 2) - a = ot.utils.unif(n) - b = ot.utils.unif(n) + a_np = ot.utils.unif(n) + b_np = ot.utils.unif(n) M = ot.dist(x, y) M = M / M.max() reg_m = 100 - a, b, M = nx.from_numpy(a, b, M) + a, b, M = nx.from_numpy(a_np, b_np, M) G_kl, _ = ot.unbalanced.mm_unbalanced(a, b, M, reg_m=reg_m, div='kl', - verbose=True, log=True) - loss_kl = nx.to_numpy(ot.unbalanced.mm_unbalanced2( - a, b, M, reg_m, div='kl', verbose=True)) + verbose=False, log=True) + loss_kl = nx.to_numpy( + ot.unbalanced.mm_unbalanced2(a, b, M, reg_m, div='kl', verbose=True) + ) G_l2, _ = ot.unbalanced.mm_unbalanced(a, b, M, reg_m=reg_m, div='l2', verbose=False, log=True) # check if the marginals come close to the true ones when large reg - np.testing.assert_allclose(np.sum(nx.to_numpy(G_kl), 1), a, atol=1e-03) - np.testing.assert_allclose(np.sum(nx.to_numpy(G_kl), 0), b, atol=1e-03) - np.testing.assert_allclose(np.sum(nx.to_numpy(G_l2), 1), a, atol=1e-03) - np.testing.assert_allclose(np.sum(nx.to_numpy(G_l2), 0), b, atol=1e-03) + np.testing.assert_allclose(np.sum(nx.to_numpy(G_kl), 1), a_np, atol=1e-03) + np.testing.assert_allclose(np.sum(nx.to_numpy(G_kl), 0), b_np, atol=1e-03) + np.testing.assert_allclose(np.sum(nx.to_numpy(G_l2), 1), a_np, atol=1e-03) + np.testing.assert_allclose(np.sum(nx.to_numpy(G_l2), 0), b_np, atol=1e-03) # check if mm_unbalanced2 returns the correct loss np.testing.assert_allclose(nx.to_numpy(nx.sum(G_kl * M)), loss_kl, @@ -324,15 +325,16 @@ def test_mm_convergence(nx): a_np, b_np = np.array([]), np.array([]) a, b = nx.from_numpy(a_np, b_np) - G_kl_null = ot.unbalanced.mm_unbalanced(a, b, M, reg_m=reg_m, div='kl') - G_l2_null = ot.unbalanced.mm_unbalanced(a, b, M, reg_m=reg_m, div='l2') - np.testing.assert_allclose(G_kl_null, G_kl) - np.testing.assert_allclose(G_l2_null, G_l2) + G_kl_null = ot.unbalanced.mm_unbalanced(a, b, M, reg_m=reg_m, div='kl', verbose=False) + G_l2_null = ot.unbalanced.mm_unbalanced(a, b, M, reg_m=reg_m, div='l2', verbose=False) + np.testing.assert_allclose(nx.to_numpy(G_kl_null), nx.to_numpy(G_kl)) + np.testing.assert_allclose(nx.to_numpy(G_l2_null), nx.to_numpy(G_l2)) # test when G0 is given G0 = ot.emd(a, b, M) + G0_np = nx.to_numpy(G0) reg_m = 10000 - G_kl = ot.unbalanced.mm_unbalanced(a, b, M, reg_m=reg_m, div='kl', G0=G0) - G_l2 = ot.unbalanced.mm_unbalanced(a, b, M, reg_m=reg_m, div='l2', G0=G0) - np.testing.assert_allclose(G0, G_kl, atol=1e-05) - np.testing.assert_allclose(G0, G_l2, atol=1e-05) + G_kl = ot.unbalanced.mm_unbalanced(a, b, M, reg_m=reg_m, div='kl', G0=G0, verbose=False) + G_l2 = ot.unbalanced.mm_unbalanced(a, b, M, reg_m=reg_m, div='l2', G0=G0, verbose=False) + np.testing.assert_allclose(G0_np, nx.to_numpy(G_kl), atol=1e-05) + np.testing.assert_allclose(G0_np, nx.to_numpy(G_l2), atol=1e-05) -- cgit v1.2.3 From e547fe30c59be72ae93c9f017786477b2652776f Mon Sep 17 00:00:00 2001 From: Nathan Cassereau <84033440+ncassereau-idris@users.noreply.github.com> Date: Mon, 13 Jun 2022 14:49:55 +0200 Subject: [MRG] Correct pointer overflow in EMD (#381) * avoid overflow on openmp version of emd solver * monothread version updated * Fixed typo in readme * added PR in releases * typo in releases.md * added a precision to releases.md * added a precision to releases.md * correct readme * forgot to cast * lower error --- README.md | 4 ++-- RELEASES.md | 5 ++++- ot/lp/EMD_wrapper.cpp | 36 ++++++++++++++++++------------------ ot/lp/network_simplex_simple_omp.h | 4 ++-- 4 files changed, 26 insertions(+), 23 deletions(-) (limited to 'ot/lp/network_simplex_simple_omp.h') diff --git a/README.md b/README.md index 12340d5..7c9475b 100644 --- a/README.md +++ b/README.md @@ -26,8 +26,8 @@ POT provides the following generic OT solvers (links to examples): * Debiased Sinkhorn barycenters [Sinkhorn divergence barycenter](https://pythonot.github.io/auto_examples/barycenters/plot_debiased_barycenter.html) [37] * [Smooth optimal transport solvers](https://pythonot.github.io/auto_examples/plot_OT_1D_smooth.html) (dual and semi-dual) for KL and squared L2 regularizations [17]. * Weak OT solver between empirical distributions [39] -* Non regularized [Wasserstein barycenters [16] ](https://pythonot.github.io/auto_examples/barycenters/plot_barycenter_lp_vs_entropic.html)) with LP solver (only small scale). -* [Gromov-Wasserstein distances](https://pythonot.github.io/auto_examples/gromov/plot_gromov.html) and [GW barycenters](https://pythonot.github.io/auto_examples/gromov/plot_gromov_barycenter.html) (exact [13] and regularized [12]), differentiable using gradients from +* Non regularized [Wasserstein barycenters [16] ](https://pythonot.github.io/auto_examples/barycenters/plot_barycenter_lp_vs_entropic.html) with LP solver (only small scale). +* [Gromov-Wasserstein distances](https://pythonot.github.io/auto_examples/gromov/plot_gromov.html) and [GW barycenters](https://pythonot.github.io/auto_examples/gromov/plot_gromov_barycenter.html) (exact [13] and regularized [12]), differentiable using gradients from Graph Dictionary Learning [38] * [Fused-Gromov-Wasserstein distances solver](https://pythonot.github.io/auto_examples/gromov/plot_fgw.html#sphx-glr-auto-examples-plot-fgw-py) and [FGW barycenters](https://pythonot.github.io/auto_examples/gromov/plot_barycenter_fgw.html) [24] * [Stochastic solver](https://pythonot.github.io/auto_examples/others/plot_stochastic.html) and diff --git a/RELEASES.md b/RELEASES.md index fdaff59..b384617 100644 --- a/RELEASES.md +++ b/RELEASES.md @@ -13,7 +13,10 @@ - Fixed an issue where Sinkhorn solver assumed a symmetric cost matrix (Issue #374, PR #375) - Fixed an issue where hitting iteration limits would be reported to stderr by std::cerr regardless of Python's stderr stream status (PR #377) - Fixed an issue where the metric argument in ot.dist did not allow a callable parameter (Issue #378, PR #379) -- Fixed an issue where the max number of iterations in ot.emd was not allow to go beyond 2^31 (PR #380) +- Fixed an issue where the max number of iterations in ot.emd was not allowed to go beyond 2^31 (PR #380) +- Fixed an issue where pointers would overflow in the EMD solver, returning an +incomplete transport plan above a certain size (slightly above 46k, its square being +roughly 2^31) (PR #381) ## 0.8.2 diff --git a/ot/lp/EMD_wrapper.cpp b/ot/lp/EMD_wrapper.cpp index 457216b..4aa5a6e 100644 --- a/ot/lp/EMD_wrapper.cpp +++ b/ot/lp/EMD_wrapper.cpp @@ -24,7 +24,7 @@ int EMD_wrap(int n1, int n2, double *X, double *Y, double *D, double *G, // beware M and C are stored in row major C style!!! using namespace lemon; - int n, m, cur; + uint64_t n, m, cur; typedef FullBipartiteDigraph Digraph; DIGRAPH_TYPEDEFS(Digraph); @@ -51,15 +51,15 @@ int EMD_wrap(int n1, int n2, double *X, double *Y, double *D, double *G, // Define the graph - std::vector indI(n), indJ(m); + std::vector indI(n), indJ(m); std::vector weights1(n), weights2(m); Digraph di(n, m); - NetworkSimplexSimple net(di, true, n+m, ((int64_t)n)*((int64_t)m), maxIter); + NetworkSimplexSimple net(di, true, (int) (n + m), n * m, maxIter); // Set supply and demand, don't account for 0 values (faster) cur=0; - for (int i=0; i0) { weights1[ cur ] = val; @@ -70,7 +70,7 @@ int EMD_wrap(int n1, int n2, double *X, double *Y, double *D, double *G, // Demand is actually negative supply... cur=0; - for (int i=0; i0) { weights2[ cur ] = -val; @@ -79,12 +79,12 @@ int EMD_wrap(int n1, int n2, double *X, double *Y, double *D, double *G, } - net.supplyMap(&weights1[0], n, &weights2[0], m); + net.supplyMap(&weights1[0], (int) n, &weights2[0], (int) m); // Set the cost of each edge int64_t idarc = 0; - for (int i=0; i indI(n), indJ(m); + std::vector indI(n), indJ(m); std::vector weights1(n), weights2(m); Digraph di(n, m); - NetworkSimplexSimple net(di, true, n+m, ((int64_t)n)*((int64_t)m), maxIter, numThreads); + NetworkSimplexSimple net(di, true, (int) (n + m), n * m, maxIter, numThreads); // Set supply and demand, don't account for 0 values (faster) cur=0; - for (int i=0; i0) { weights1[ cur ] = val; @@ -172,7 +172,7 @@ int EMD_wrap_omp(int n1, int n2, double *X, double *Y, double *D, double *G, // Demand is actually negative supply... cur=0; - for (int i=0; i0) { weights2[ cur ] = -val; @@ -181,12 +181,12 @@ int EMD_wrap_omp(int n1, int n2, double *X, double *Y, double *D, double *G, } - net.supplyMap(&weights1[0], n, &weights2[0], m); + net.supplyMap(&weights1[0], (int) n, &weights2[0], (int) m); // Set the cost of each edge int64_t idarc = 0; - for (int i=0; i::epsilon()*10 -#define _EPSILON 1e-8 +#define EPSILON std::numeric_limits::epsilon() +#define _EPSILON 1e-14 #define MAX_DEBUG_ITER 100000 /// \ingroup min_cost_flow_algs -- cgit v1.2.3 From fa0d4f2afff73284f4b79bfebb085eed332c112f Mon Sep 17 00:00:00 2001 From: Nathan Cassereau <84033440+ncassereau-idris@users.noreply.github.com> Date: Mon, 21 Nov 2022 17:27:50 +0100 Subject: [MRG] Replaces numpy compiler with setuptools (#409) * Numpy ccompiler deprecation handled with setuptools ccompiler * Remove useless OMP Macro, already provides _OPENMP * RELEASES.md * Remove forgotten temporary bug added for logging purposes --- RELEASES.md | 1 + ot/helpers/pre_build_helpers.py | 24 ++---------------------- ot/lp/network_simplex_simple_omp.h | 6 +++--- setup.py | 2 +- test/test_ot.py | 16 ++++++++++++++++ 5 files changed, 23 insertions(+), 26 deletions(-) (limited to 'ot/lp/network_simplex_simple_omp.h') diff --git a/RELEASES.md b/RELEASES.md index 1c7b7da..564fd4a 100644 --- a/RELEASES.md +++ b/RELEASES.md @@ -25,6 +25,7 @@ roughly 2^31) (PR #381) - Fixed an issue where a pytorch example would throw an error if executed on a GPU (Issue #389, PR #391) - Added a work-around for scipy's bug, where you cannot compute the Hamming distance with a "None" weight attribute. (Issue #400, PR #402) - Fixed an issue where the doc could not be built due to some changes in matplotlib's API (Issue #403, PR #402) +- Replaced Numpy C Compiler with Setuptools C Compiler due to deprecation issues (Issue #408, PR #409) ## 0.8.2 diff --git a/ot/helpers/pre_build_helpers.py b/ot/helpers/pre_build_helpers.py index 93ecd6a..2930036 100644 --- a/ot/helpers/pre_build_helpers.py +++ b/ot/helpers/pre_build_helpers.py @@ -4,34 +4,14 @@ import os import sys import glob import tempfile -import setuptools # noqa import subprocess -from distutils.dist import Distribution -from distutils.sysconfig import customize_compiler -from numpy.distutils.ccompiler import new_compiler -from numpy.distutils.command.config_compiler import config_cc +from setuptools.command.build_ext import customize_compiler, new_compiler def _get_compiler(): - """Get a compiler equivalent to the one that will be used to build POT - Handles compiler specified as follows: - - python setup.py build_ext --compiler= - - CC= python setup.py build_ext - """ - dist = Distribution({'script_name': os.path.basename(sys.argv[0]), - 'script_args': sys.argv[1:], - 'cmdclass': {'config_cc': config_cc}}) - - cmd_opts = dist.command_options.get('build_ext') - if cmd_opts is not None and 'compiler' in cmd_opts: - compiler = cmd_opts['compiler'][1] - else: - compiler = None - - ccompiler = new_compiler(compiler=compiler) + ccompiler = new_compiler() customize_compiler(ccompiler) - return ccompiler diff --git a/ot/lp/network_simplex_simple_omp.h b/ot/lp/network_simplex_simple_omp.h index c324d4c..890b7ab 100644 --- a/ot/lp/network_simplex_simple_omp.h +++ b/ot/lp/network_simplex_simple_omp.h @@ -67,7 +67,7 @@ //#include "core.h" //#include "lmath.h" -#ifdef OMP +#ifdef _OPENMP #include #endif #include @@ -254,7 +254,7 @@ namespace lemon_omp { // Reset data structures reset(); max_iter = maxiters; -#ifdef OMP +#ifdef _OPENMP if (max_threads < 0) { max_threads = omp_get_max_threads(); } @@ -513,7 +513,7 @@ namespace lemon_omp { int j; #pragma omp parallel { -#ifdef OMP +#ifdef _OPENMP int t = omp_get_thread_num(); #else int t = 0; diff --git a/setup.py b/setup.py index c03191a..dc9066d 100644 --- a/setup.py +++ b/setup.py @@ -37,7 +37,7 @@ compile_args = ["/O2" if sys.platform == "win32" else "-O3"] link_args = [] if openmp_supported: - compile_args += flags + ["/DOMP" if sys.platform == 'win32' else "-DOMP"] + compile_args += flags link_args += flags if sys.platform.startswith('darwin'): diff --git a/test/test_ot.py b/test/test_ot.py index 9a4e175..f2338ac 100644 --- a/test/test_ot.py +++ b/test/test_ot.py @@ -204,6 +204,22 @@ def test_emd_emd2(): np.testing.assert_allclose(w, 0) +def test_omp_emd2(): + # test emd2 and emd2 with openmp for simple identity + n = 100 + rng = np.random.RandomState(0) + + x = rng.randn(n, 2) + u = ot.utils.unif(n) + + M = ot.dist(x, x) + + w = ot.emd2(u, u, M) + w2 = ot.emd2(u, u, M, numThreads=2) + + np.testing.assert_allclose(w, w2) + + def test_emd_empty(): # test emd and emd2 for simple identity n = 100 -- cgit v1.2.3