MATLAB: Linprogr in Matlab not finding a solution when a solution exists

linprogrsolution

I am solving a linear programming in Matlab using linprogr. All the matrices are attached. 
clear
rng default
load matrices
f=zeros(size(x,1),1);
possible_solution = linprog(f,Aineq,bineq,Aeq,beq, lb, ub);
The output is 
No feasible solution found.
Linprog stopped because no point satisfies the constraints.
From some theory behind the problem, I know that the vector x (uploaded together with the other matrices) should be a solution of linprog(f,Aineq,bineq,Aeq,beq, lb, ub).
In fact, x seems to satisfy all the constraints, except the equality constraints. However, regarding the equality constraints, the maximum absolute difference between Aeq*x and beq is around 3*10^(-15). Hence, also the equality constraints are almost satisfied by x.
load x
check_lb=sum(x>=0)==size(x,1);
check_ub=sum(x<=1)==size(x,1);
check_ineq=(sum(Aineq*x<=bineq)==size(Aineq,1));
check_eq=max(abs((Aeq*x-beq)));
Therefore, my question is: why linprog(f,Aineq,bineq,Aeq,beq, lb, ub) does not find x as solution? It does not seem to be a problem of the equality constraints, because if I remove the inequality constraints and run 
rng default
possible_solution_no_inequalities = linprog(f,[],[],Aeq,beq, lb, ub);
the program finds a solution. Is it a matter of numerical precision? How can I control for that?

Best Answer

It looks to me LINPROG is flawed or buggy in your case. I try to remove suspected constraints and change beq so that the equality is satisfied by x, and LINPROG still fails.
load matrices
f=zeros(size(x,1),1);
Aeq = round(Aeq);
beq = Aeq*x;
keep = max(abs(Aineq),[],2)>1e-10;
all(x>=lb)
all(x<=ub)
all(Aineq*x<=bineq)
all(Aeq*x==beq)
% This will fail to return solution
options = optimoptions('linprog','Algorithm','interior-point','ConstraintTolerance',1e-3);
xlinprog = linprog(f,Aineq(keep,:),bineq(keep,:),Aeq,beq, lb, ub, options);
But if I call QUADPROG it will returns a solution
% This will returns a solution
H = sparse(size(x,1),size(x,1));
xquadprog = quadprog(H,f,Aineq(keep,:),bineq(keep,:),Aeq,beq, lb, ub);
Related Question