clear;
clf;
d1=31; d2=22;
% Code from the exercise description
%%
n=8;
dummystep=31*d1+d2;
rng('default');
for i=1:dummystep
dummy=rand;
end
Corr=zeros(n,n);
for i=1:n
for j=1:n
Corr(i,j)=(-1)^abs(i-j)/(abs(i-j)+1);
end
end
sigma=zeros(n,1);
mu=zeros(n,1);
sigma(1)=2;
mu(1)=3;
for i=1:n-1
sigma(i+1)=sigma(i)+2*rand;
mu(i+1)=mu(i)+1;
end
D=diag(sigma);
C2=D*Corr*D;
C=0.5*(C2+C2');
%%
% End of code from the exercise description
% Exercise 1
mux1 = [];
sigmax1 = [];
for r = (3:0.25:9)
    Ain = [];
    bin = [];
    Aeq = [mu'; ones(1, n)];
    beq = [r 1]';
    f = [];
    lb = zeros(8, 1);
    X = quadprog(C, f, Ain, bin, Aeq, beq, lb, []);
    mux1 = [mux1 mu'*X];
    sigmax1 = [sigmax1 sqrt(X'*C*X)];
end
% Exercise 2
mux2 = [];
sigmax2 = [];
for r = (3:0.25:9)
    Ain = ones(1, n);
    bin = 1;
    Aeq = mu';
    beq = r;
    f = [];
    lb = zeros(8, 1);
    X = quadprog(C, f, Ain, bin, Aeq, beq, lb, []);
    mux2 = [mux2 mu'*X];
    sigmax2 = [sigmax2 sqrt(X'*C*X)];
    disp(sum(X))
end
figure(1);
hold on
plot(sigmax1, mux1);
plot(sigmax2, mux2);
title('Exercise 1 & 2');
legend('1', '2')
% Exercise 3
mux3 = [];
sigmax3 = [];
for r = (3:0.25:9)
    Ain = -1*mu';
    bin = -r;
    Aeq = ones(1, n);
    beq = 1;
    f = [];
    lb = zeros(8, 1);
    X = quadprog(C, f, Ain, bin, Aeq, beq, lb, []);
    mux3 = [mux3 mu'*X];
    sigmax3 = [sigmax3 sqrt(X'*C*X)];
end
figure(2);
plot(sigmax1, mux1);
plot(sigmax3, mux3);
legend('1', '3')
title('Exercise 1 & 3');
% Exercise 4
mux4 = [];
sigmax4 = [];
for r = (3:0.25:9)
    Ain = [];
    bin = [];
    Aeq = [mu'; ones(1, n)];
    beq = [r 1]';
    f = [];
    lb = [];
    X = quadprog(C, f, Ain, bin, Aeq, beq, lb, []);
    mux4 = [mux4 mu'*X];
    sigmax4 = [sigmax4 sqrt(X'*C*X)];
end
figure(3);
hold on
plot(sigmax1, mux1);
plot(sigmax4, mux4);
legend('1', '4')
title('Exercise 4');