function [A,b,v,x0,zo,options] = genProb(params,options)
% function [A,b,v,x0] = genProb(params)
%
% Take provided parameters and output a randomly generated problem:
% min ||v-x|| s.t. Ax=b, x>=0
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% INPUT
% params - struct with option parameters
% .free - There exists free variables 
% .m - number of constraints of A
% .m1 - number of basic variables (when free is true)
% .n - number of variables of A
% .n2 - number of free variables
% .density - density of matrix A
% .simplex - m=1
% .degen - generate a problem with degenerate solution 
% .numdegen - number of degenerate variables in BFS
% .xvertex - generate problem with vertex closest solution
% .polarconev - v is in polarcone of vertex
% .randxpos - ?
% .lengthv - norm of v

if isfield(params, 'free')
    free=params.free;
else
    free=false;
end

if free
    %%% Total number of variables
    if isfield(params, 'n')
        n=params.n;
    else
        n=1000;
    end
    
    %%% Number of free variables
    if isfield(params, 'n2')
        n2=params.n2;
    else
        n2=100;
    end

    %%% Number of nonfree variables
    n1=n-n2;

    %%% Number of constraints
    if isfield(params, 'm')
        m=params.m;
    else
        m=300;
    end
    
    %%% Number of nonfree basic variables
    m1=m-n2;

else
    if isfield(params, 'm')
        m=params.m;
    else
        m=300;
    end
    
    if isfield(params, 'n')
        n=params.n;
    else
        n=1000;
    end
end

if isfield(params, 'density')
    density=params.density;
else
    density=0.08;
end

if isfield(params, 'simplex')
    simplex=params.simplex;
else
    simplex=false;
end

if isfield(params, 'degen')
    degen=params.degen;
else
    degen=false;
end

if degen 
    if isfield(params, 'numdegen')
        numdegen=params.numdegen; % number of deg. basic var.
    else
        numdegen=5;
    end
end

if isfield(params, 'xvertex')
    xvertex=params.xvertex;
else
    xvertex=true;
end

if isfield(params, 'polarconev')
    polarconev=params.polarconev;
else
    polarconev=true;
end

if isfield(params, 'randxpos')
    randxpos=params.randxpos;
else
    randxpos=0;
end

if isfield(params, 'lengthv')
    lengthv=params.lengthv;
else
    lengthv=1;
end

if simplex
	m = 1;
	density = 1;
    A = ones(1,n);     % to compare onto simplex with Ting Kei code???
    b = 10^sqrt(expon);
else
    if free             %%% If there are free variables
        if degen        %%% If the solution is degenerate
            tic
            rc = .01;
            xB = 10*(rand(m1,1)+.1);   % nondeg BFS also ensures strict feas.
            nperm = randperm(n);
            m1perm = nperm(1:m1);
            ndegenperm = randperm(length(m1perm),numdegen);
            xB(ndegenperm) = 0;       % makes nondegen BFS a degen BFS
            nonbas = setdiff(nperm(1:n1),m1perm);  % x nonbasic variables so z>0
            free = nperm(n1+1:end);
            notfree = nperm(1:n1);
	        options.free = free;

            A(:,nonbas) = sprandn(m,n-m,density,rc);
            A(:,[m1perm free]) = sprandn(m,m,density,rc); % ensure ONE proper basis matrix
            x0 = zeros(n,1);
            x0(m1perm) = xB;   % for a ***nondeg.*** BFS
            x0(free) = randn(n2,1);
            x0 = x0/norm(x0);
            nA = norm(A,'fro');
            A = sparse(A/nA);
            b = A*x0;
            fprintf('time to generate A,b,xo, random problem %g\n',toc)
            options.deginds=m1perm(ndegenperm);
            options.m1perm=m1perm;
        else        
            tic
            rc = .01;
            xB = 10*(rand(m1,1)+.1);   % nondeg BFS also ensures strict feas.
            nperm = randperm(n);
            m1perm = nperm(1:m1);
            nonbas = setdiff(nperm(1:n1),m1perm);  % x nonbasic variables so z1>0
            free = nperm(n1+1:end);
            notfree = nperm(1:n1);
	        options.free = free;

            A(:,nonbas) = sprandn(m,n-m,density,rc);
            A(:,[m1perm free]) = sprandn(m,m,density,rc); % ensure ONE proper basis matrix
            x0 = zeros(n,1);
            x0(m1perm) = xB;   % for a ***nondeg.*** BFS
            x0(free) = randn(n2,1);
            x0 = x0/norm(x0);
            nA = norm(A,'fro');
            A = sparse(A/nA);
            b = A*x0;
            fprintf('time to generate random problem %g\n',toc)
            deginds=0;
            options.m1perm=m1perm;
        end % End of if degen
    else            %%% Else there are no free variables
        if degen
            tic
            rc = .01;
            xB = 10*(rand(m,1)+.1);   % would be a BFS also ensures strict feas.
            nperm = randperm(n,m);
            nonbas = setdiff(1:n,nperm);  % x nonbasic variables so z>0
            ndegenperm = randperm(length(nperm),numdegen);
            xB(ndegenperm) = 0;       % makes nondegen BFS a degen BFS
            A(:,nonbas) = sprandn(m,n-m,density,rc);
            A(:,nperm) = sprandn(m,m,density,rc); % ensure ONE proper basis matrix
            x0 = zeros(n,1);
            x0(nperm) = xB;   % for a ***deg.*** BFS
            x0 = x0/norm(x0); % Need this still?
            nA = norm(A,'fro');
            A = sparse(A/nA); % Need this still?
            b = A*x0;
            fprintf('time to generate random problem %g\n',toc)
            options.deginds=nperm(ndegenperm);
            options.nperm=nperm;
        else
            tic
            rc = .01;
            xB = 10*(rand(m,1)+.1);   % nondeg BFS also ensures strict feas.
            nperm = randperm(n,m);
            nonbas = setdiff(1:n,nperm);  % x nonbasic variables so z>0
            A(:,nonbas) = sprandn(m,n-m,density,rc);
            A(:,nperm) = sprandn(m,m,density,rc); % ensure ONE proper basis matrix
            x0 = zeros(n,1);
            x0(nperm) = xB;   % for a ***nondeg.*** BFS
            x0 = x0/norm(x0);
            nA = norm(A,'fro');
            A = sparse(A/nA);
            b = A*x0;
            fprintf('time to generate random problem %g\n',toc)
            options.deginds=0;
            options.nperm=nperm;
        end     % End of if degen
    end     % End of if free
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% Start Generating v
if polarconev   % choose obj so x*-v is in polar cone
    if xvertex % ensure a vertex is opt.
        fprintf('Choosing v to ensure a NICE vertex is optimal\n')
        options.nperm = nperm;
        xo = x0;  % choose the constructed nice vertex
        zo = zeros(n,1);
        zo(nonbas) = lengthv*(rand(length(nonbas),1)+.1);
        yo = randn(m,1);
        v = xo -A'*yo - zo;
    elseif randxpos == 0 % choose a random number of pos
        xpnum = randi(n-m,1,1); % m+xpnum els of x > 0
        fprintf('Choosing v to ensure opt fixed x with ')
        fprintf('m+%i = %i positive xi''s \n',xpnum,m+xpnum)
        xo = x0; 
        randxN = randperm(length(nonbas));
        %% add pos els to BFS xo
        xo(nonbas(randxN(1:xpnum))) = 10*(rand(xpnum,1)+.1);
        xo = xo/norm(xo);
        b = A*xo;
        xzdeg = randi(min(m,n-m)-1,1,1); % # overlapping zeros
        zo = zeros(n,1);
        zo(nonbas(randxN(xpnum+1:end))) = ...
            lengthv*(rand(length(nonbas)-xpnum,1)+.1);
        indszpos = find(zo>0);
        randzpos = randperm(length(indszpos));
        xzdeg = min(20,length(randzpos));  %????
        zo(indszpos(randzpos(1:xzdeg))) = 0; % set degeneracy
        fprintf('xzdeg: setting %i els of zo to zero \n',xzdeg)
        yo = randn(m,1);
        v = xo -A'*yo - zo;
    end
else
    v = randn(n,1);
    v = lengthv*v/norm(v);   % how large should v be?????
    zo=0;
end
%%% Finish Generating v
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

end  % end of function

