===== Cvičení 5 =====
<code matlab>
%figure;
hold on;

%primka
a = rand(1);
b = rand(1);
xLine = (-5:5);
yLine = a*xLine + b;

plot(xLine, yLine, 'y--')

%body v normalnim rozdeleni
epsilon = 0.05;

v = randn(2, 100);

sx = v(1, 1:length(v));
x = sx(:);

y = v(2, 1:length(v));
y = a * x + b+ y(:);

plot(x, y, 'r+')

%robustni prokladani primarni ulohy
heigth = length(x);

y1 = y+epsilon; 
y2 = -y+epsilon;
y3 = 0;

A = [ x ones(heigth,1) -1*eye(heigth); -1*x -1*ones(heigth,1) -1*eye(heigth); zeros(heigth,1) zeros(heigth,1) -1*eye(heigth) ];

b = [y1;y2;zeros(heigth,1)];

f = [0; 0; ones(heigth,1)];

v = linprog(f, A, b );

PrimaryOptimum = f'*v


plot(x,y,'r+');

lineX = [-5:0.1:5];
lineX = lineX(:);
lineY = v(1)*lineX + v(2);

plot(lineX,lineY, 'b');
plot(lineX,lineY+epsilon, 'g--');
plot(lineX,lineY-epsilon, 'g--');

hold off;
%robustni prokladani dualni ulohy
fD = -b';
AD = A';
bD = f';

AeqD = (A'); %tady skromnej trik...nema to bejt omezeny na nerovnost, ale na rovnost....
AeqD = AeqD(1:2,:); %protoze jsem linej psat znova ty nerovnice (aD, bD), tak je tam necham, jak jsou.
beqD = f'; %Protoze je rovnost tvrdsi omezeni, tak se jejim pridanim nic neporusi a docilim
beqD = beqD(1:2); %stejnych vlasnosti, jako kdyby aD, bD ony dve rovnice nebyly a byly jen v AeqD a beqD

v2 = linprog(fD, AD, bD,AeqD,beqD,[],zeros(length(b), 1));
SecondaryOptimum = -fD*v2

%vykresleni divergenci a odpovidajicich dualnich promennych
figure
hold on
xAxis = [1:1:length(v2)];
plot(xAxis, v2, 'r+')
hold off;
</code>

==== Bonus 1 ====
 --- //[[sachamar@fel.cvut.cz|Marek Sacha]] 2009/10/29 17:16//
<code matlab>
clear
clc
load data1.mat
% Construct inequalities matrix
%       | x   1  -1 |
% A =   |-x  -1  -1 |
% 
A = [ x ones(50,1) -1*ones(50,1); -1*x -1*ones(50,1) -1*ones(50,1) ];
% Construct right side
%     |  y |
% b = | -y |
b = [ y; -y ];
% Construct function to be minimized
% f = 0x + 0y + lambda
f = [ 0; 0; 1];
% Run linprog
o = linprog(f,A,b);
% o = [ a b lambda ];


% Modify linprog parameters for secondary problem    
% Run linprog
f2 = b';
A2 = -1*A';
b2 = f';

[ awh a2w ] = size(A2);
[ f2h f2w ] = size(f2);
lb = zeros(f2h,f2w);
% Dont forget to add equality conditions (for first two rows of A2,b2)
% - a, b are from 'R' not R+
o2 = linprog(f2,A2,b2,A2(1:2,[1:a2w]),b2(1:2),lb);


% Verify results
result_prim = (-f)'*o
result_sec = b'*o2


hold on;
% Print separate points
plot(x,y,'y.');
% Approximation line
plot(x,x*o(1)+o(2),'r');
% Top border
plot(x,x*o(1)+o(2)+o(3),'g');
% Bottom border
plot(x,x*o(1)+o(2)-o(3),'g');
hold off;

% New Plot
figure;

% Visualize dual
plot([x ; x],[y ; -y].*o2,'ro');
</code>

==== Bonus 2 ====
 --- //[[sachamar@fel.cvut.cz|Marek Sacha]] 2009/10/29 17:16//
<code matlab>
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Author: Marek Sacha <sachamar@fel.cvut.cz> %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

clear all
clc

numberOfPoints = [];
primal = [];
dual = [];

options = optimset('LargeScale','off','Simplex','on','MaxIter',1e5);
 
for j=1:4
    
    t = j*100;    
    realloc = [ numberOfPoints t ];
    numberOfPoints = realloc;
    
    % initialize random matrix
    r = rand(t,2);
    x = r(:,1);
    y = r(:,2);

    % Initialize e


    e = 0.05;

    [ xh xw ] = size(x);
    % Construct inequalities matrix
    %         a    b   l1  l2 ... ln 
    %        ++++++++++++++++++++++++
    %       | x1   1  -1   0  ...  0 |
    %       | x2   1   0  -1  ...  0 |
    %       ...
    %       | xn   1   0   0  ... -1 |
    % A =   |-x1  -1  -1   0  ...  0 |
    %       |-x2  -1   0  -1  ...  0 |
    %       ...
    %       | 0    0   0   0  ... -1 |
    %        ++++++++++++++++++++++++
    % 
    A = [ x ones(xh,1) -1*eye(xh); -1*x -1*ones(xh,1) -1*eye(xh); zeros(xh,1) zeros(xh,1) -1*eye(xh) ];
    % Construct right side
    %     | e + y |
    % b = | e - y |
    %     | 0     |
    b = [ (e*ones(xh,1))+y; (e*ones(xh,1))-y; zeros(xh,1) ];
    % Construct function to be minimized
    % f = 0a + 0b + lambda1 + lambda2 + lambda3 + ... + lambdan
    f = [ 0 ; 0; ones(xh,1)];

    % Run linprog
    tic;
    o = linprog(f,A,b,[],[],[],[],[],options);
    time = toc;
    
    realloc = [ primal time ];
    primal = realloc;

    % Modify linprog parameters for secondary problem    
    f2 = b';
    A2 = -1*A';
    b2 = f';

    [ awh a2w ] = size(A2);
    [ f2h f2w ] = size(f2);
    lb = zeros(f2h,f2w);
    % Dont forget to add equality conditions (for first two rows of A2,b2)
    % - a, b are from 'R' not R+
    
    tic;
    o2 = linprog(f2,A2,b2,A2(1:2,[1:a2w]),b2(1:2),lb,[],[],options);
    time = toc;
    
    realloc = [ dual time ];
    dual = realloc;

    % Verify results
    result_prim = (-f)'*o;
    result_sec = b'*o2;

end

hold on;
plot(numberOfPoints,primal,'g-');
plot(numberOfPoints,dual,'r-');
hold off;

</code>

~~DISCUSSION~~