%% Works
options = optimset('Disp','none','Algorithm','active-set');
x0 = [0 1.5];
lb = [-3 -3];
ub = - lb;
tic
[x,fval,exitflag,output] = fmincon(@peaksObj,x0,[],[],[],[],lb,ub,@peaksCon,options);
fminconTot = output.funcCount;
tf1 = toc;
% rerun to remove JIT influence
tic
[x,fval,exitflag,output] = fmincon(@peaksObj,x0,[],[],[],[],lb,ub,@peaksCon,options);
fminconTot = output.funcCount;
tf1 = toc;

%% Try a random starting point (uniform distribution grid of 4 points)
a = -3;
b = 3;
x1 = a + (b-a) * rand(2);
x2 = a + (b-a) * rand(2);
x0r = [x1; x2];
fminconCount = 0;
tic
for i = 1:length(x0r)
    [xopt(i,:),fval,exitflag,output] = fmincon(@peaksObj,x0r(i,:),[],[],[],[],lb,ub,@peaksCon,options);
    fminconCount = output.funcCount + fminconCount;
end
tf2 = toc;

%% GA
tic
options = gaoptimset('Disp','none');
[x,fval,exitflag,output] = ga(@peaksObj,2,[],[],[],[],lb,ub,@peaksCon);%,options)
gaTot = output.funccount;
tga = toc;

%% SA
tic
options = saoptimset('Display','none');
[x,fval,exitflag,output] = simulannealbnd(@peaksObj,x0,lb,ub,options);
saTot = output.funccount;
tsa = toc;

%% PS
tic
options = psoptimset('Display','none','CompleteSearch','on','SearchMethod',@searchlhs);
[x,fval,exitflag,output] = patternsearch(@peaksObj,x0,[],[],[],[],lb,ub,@peaksCon,options);
psTot = output.funccount;
tps = toc;

%% Tabulate Results
fprintf('\n\nResults from different optimization solvers (default settings):\n');
fprintf('\nSolver\t\t\t\t\tFcn Calls\t\tTime(s)\n');
fprintf('--------------------------------------------------------------\n')
fprintf('fmincon (lucky guess)   \t%d \t\t\t%f\n',fminconTot,tf1);
fprintf('fmincon (random start)  \t%d \t\t%f\n',fminconCount,tf2);
fprintf('Genetic Algorithm       \t%d \t\t%f\n',gaTot,tga);
fprintf('Simulated Annealing     \t%d \t\t%f\n',saTot,tsa);
fprintf('Pattern Search          \t%d \t\t%f\n',psTot,tps);

%% Plot result
close all

