n = 3;

Q = randn(n); Q = Q'+Q; 
Q1 = randn(n); 
Q1 = Q1'+Q1;  
Q1 = Q1'*Q1;  % posdef
Q2 = randn(n); 
Q2 = Q2'+Q2;  
% Q2 = Q2'*Q2; 
q = randn(n,1); q1 = randn(n,1); q2 = randn(n,1);
gam = randn; 
gam1 = -10+randn; gam2 = -10+randn; 


[x,lamsol,objvalsol] = qcqp2(Q,Q1,Q2,q,q1,q2,gam1,gam2);


% compare with random solution
objval = inf;
dobreak = 0;
for ii = 1:300
    d = randn(n,1);        
%    d = d/norm(d);    d = d*Delta*rand;
    it = 0;
    while d'*Q1*d/2+q1'*d+gam1>0 | d'*Q2*d/2+q2'*d+gam2>0
    d = randn(n,1);   %   d = d*Delta*rand;    
    it = it+1;
    if it>1000, dobreak = 1; break,     end
    end
    if dobreak, 'no feasible point seems to be present',  break,end
    objval = min(objval,d'*Q*d/2+q'*d);
    if objval ==d'*Q*d/2+q'*d,  dbestrand = d; end
end
objval

dsol = x;
H = (Q+lamsol(1)*Q1+lamsol(2)*Q2); 
%result
if exist('lamsol','var'), 
    if exist('dbestrand','var')==0, dbestrand = nan(n,1);end
randvecinfo = [dbestrand'*Q1*dbestrand/2+q1'*dbestrand+gam1 dbestrand'*Q2*dbestrand/2+q2'*dbestrand+gam2]

res1=dsol'*Q1*dsol/2+q1'*dsol+gam1;
res2=dsol'*Q2*dsol/2+q2'*dsol+gam2;
disp('   objdiff      objvalsol     objvalrand       lam         mu           res1        res2         cond(H)  ')
disp([objvalsol-objval objvalsol objval lamsol(1) lamsol(2) dsol'*Q1*dsol/2+q1'*dsol+gam1 dsol'*Q2*dsol/2+q2'*dsol+gam2  cond(H)])
disp('    check neg,          * <= *       ,       nonneg,     nonneg.      and res1,res2 are nonpos, 0 if solution on boundary')
else
'seems infeasible'    
objval    
if exist('dbestrand','var'), 
randvecinfo = [dbestrand'*Q1*dbestrand/2+q1'*dbestrand+gam1 dbestrand'*Q2*dbestrand/2+q2'*dbestrand+gam2]
end
end

%

if objvalsol>objval
    'suspicious solution. possibly the solution is -infty. check last value, is it e.g. -1e10?'
% infinity check
dbestrand = 1e5*dbestrand;
disp([dbestrand'*Q1*dbestrand/2+q1'*dbestrand+gam1 dbestrand'*Q2*dbestrand/2+q2'*dbestrand+gam2 dbestrand'*Q*dbestrand/2+q'*dbestrand+gam])
[V,D] = eig(Q);[dd,ix] = min(diag(D));dbestrand = 1e5*V(:,ix);
disp([dbestrand'*Q1*dbestrand/2+q1'*dbestrand+gam1 dbestrand'*Q2*dbestrand/2+q2'*dbestrand+gam2 dbestrand'*Q*dbestrand/2+q'*dbestrand+gam])
end


