ecos crash on a very sparse problem #48

RuntimeError: Internal problem occurred in ECOS while setting up the problem.
Please send a bug report with data to Alexander Domahidi.
Email: domahidi@control.ee.ethz.ch

import numpy as np
import time
import math
import scipy
import qpsolvers

np.random.seed(1)

def relu_inv(x): d=np.empty_like(x); d[:]=x[:]; d[x[:]<=0]=np.nan; return d
def relu(x): return (x>0)*x
def invH(x): return math.sqrt(2)*scipy.special.erfcinv(2*x)

solver = "ecos"
def optimize_fixed_points(X, H, threshold=0):
    (N, P) = X.shape
    assert np.all(H.shape == (N, P))
    
    # Solve quadratic problem for each row
    W = np.empty((N, N)); W[:] = np.nan
    Q = np.eye(N-1)
    found = 0
    for ii in range(N):
        I = np.concatenate((np.array(range(ii),dtype=int),np.array(range(ii+1,N),dtype=int)))
        I = I.reshape(I.shape[0],1)
        J = np.nonzero(np.isfinite(H[ii,:])*1)[0]
        coJ = np.nonzero(np.isnan(H[ii,:])*1)[0]
        assert len(J)+len(coJ) == P
        c = np.zeros(I.shape[0])
        lb = None
        A = X[I,J].T
        b = H[ii,J]
        G = X[I,coJ].T
        z = np.ones(G.shape[0])*threshold
        w = qpsolvers.solve_qp(Q, c, G, z, A, b, lb, solver=solver, abstol=1e-8, reltol=1e-6)
        if w is not None:
            assert np.all(np.isfinite(w))
            assert np.linalg.norm(A@w-b)<1e-6
            assert np.all(G@w<z+1e-6)
            W[ii,ii] = 0
            W[ii,I.T] = w
            found += 1
    return (W, found)

# Problem parameters
N = 128; f = 1/64; K = np.round(f*N); P = 1
print("N=%d K=%d P=%d" % (N,K,P))

for i in range(100):
    tictocT = time.time()
    # Create patterns
    X = relu(np.random.randn(N, P) - invH(f))
    H = relu_inv(X)
    (W, nfound) = optimize_fixed_points(X, H)
    d = X - relu(W@X)
    assert np.linalg.norm(d[np.isfinite(d[:])]) < 1e-8, "Fixed-points violation: %1.1e"%np.linalg.norm(d[np.isfinite(d[:])])
    assert nfound<N or np.all(np.isfinite(W))
    print("Done #%d (found %d/%d, took %1.1fsec)."%(i, nfound, N, time.time()-tictocT))