Another issue on Windows with ot.da.OTDA()

[toto6]

The following script:

import numpy as np
import ot
a = np.random.rand(1500, 95)
b = np.random.rand(50000, 95)
opt = ot.da.OTDA()
opt.fit(a,b)
print(np.sum(opt.G))

returns me on Windows
0

Windows-10-10.0.15063-SP0
Python 3.6.1 |Anaconda custom (64-bit)| (default, May 11 2017, 13:25:24) [MSC v.1900 64 bit (AMD64)]
NumPy 1.12.1
SciPy 0.19.0
POT 0.3.1

while it returns me on debian
0.716

Linux-3.2.0-4-amd64-x86_64-with-debian-7.11
('Python', '2.7.3 (default, Jun 21 2016, 18:38:19) \n[GCC 4.7.2]')
('NumPy', '1.13.1')
('SciPy', '0.19.1')
('POT', '0.3.1')

If I'm not mistaken, it should always return 1

[rflamary]

The difference between the Windows and Debian result is indeed a serious problem. I personally cannot investigate this bug because I don't have any Windows available. Could you confirm that it happens only on very large problems?

The fact that it does not sum to 1 is because the LP solver has not converged.
It comes from the fact that we limit the number of iteration in the C code and it should be a parameters of emd function.

If you set a larger value manually, in the C code and compile manually it should work better. We really need to do it because your are the second to ask this question.

[toto6]

with

import numpy as np
import ot
a = np.random.rand(1500, 95)
b = np.random.rand(4000, 95)
opt = ot.da.OTDA()
opt.fit(a,b)
print(np.sum(opt.G))

it returns 1 on both Windows and Debian.
By increasing by 1000:

import numpy as np
import ot
a = np.random.rand(1500, 95)
b = np.random.rand(5000, 95)
opt = ot.da.OTDA()
opt.fit(a,b)
print(np.sum(opt.G))

it returns me 0.9992 on debian and 0 on Windows.

[toto6]

I increased as you said in file EMD_wrapper.cpp the variable max_iter and it worked!

True that it could be useful to modify the max number of iteration, but also to display in the python console if the problem is not feasable.
I saw that there is

std::cout << "Infeasible problem";

in the source code but it is not displayed in the python console

[rflamary]

Ok this clearly comes from this line:
https://github.com/rflamary/POT/blob/master/ot/lp/EMD_wrapper.cpp#L22

10000 iterations max is very small for large matrices as you are investigating.

@aje if you want to implement the parameter maxiter (should be numItermax in the emd function) you are welcome to submit a PR because I don't know if @ncourty or me have time to do that onn short notice.

I agree that the emd_c/EMD_wrap function should also return some kind of log which should be a log dictionnary for emd (as is done with sinkhorn). the cout in EMD_wrap does not print in the console indeed so it should happen in emd_c.

[toto6]

I try to make this change