main.m

clc
clear
format short%long


%% 

n = 3; % Ðàçìåðíîñòü "òî÷êè"
m = 4000; % Êîëè÷åñòâî òî÷åê

%cc
range = 10;
cc = range - rand(n,1) * 2 * range;

%Q 
% Çäåñü ïûòàåìñÿ "æåëåçíî" ñôîðìèðîâàòü ïîëíîðàíãîâóþ ìàòðèöó
isFullRank = false;
while isFullRank == false
    M = 1 - rand(n,n) * 2; %cîçäàåì ñëó÷àéíóþ ìàòðèöó M, ðàçìåðà n x n
    M = triu(M) + triu(M,1)'; %  áåðåì òîëüêî âåðõíþþ ïîëîâèíêó è îòîáðàæàåì åå îòíîñèòåëüíî ãë. äèàãîíàëè
    if rank(M) == n;
        isFullRank = true;
        Q = M;
        disp('Q - ïîëíîðàíãîâàÿ ìàòðèöà')
    end
end

%A
A = Q * Q;

%S
S = inv(Q);

%M
M = [];
for i = 1:1:m
    z = 1 - rand(n,1) * 2;
    nz = sqrt(z'*z);
    p = (Q \ z) * rand(1) / nz;   % (Q \ z) <==> (inv(Q) * z)
    Mj = [p + cc ; 1]';
    M = [M ; Mj];
end
disp(['Ðàçìåðíîñòü ìàòðèöû M: ' int2str(size(M))])

%TT
TT = Tt(M');

%L
L = ones(m,1);

%O
O = zeros(m,1);

%cF
cF = Coef(n+1);

%% Ðåøàòåëü
AA = [TT;-TT];
bb = [L;O];
[aa fval] = linprog(-cF,AA,bb);


%% Ïðîâåðêè

maxL = max(TT*aa);
minL = min(TT*aa);

disp(['max(TT*aa) = ' num2str(maxL)])
disp(['min(TT*aa) = ' num2str(minL)])

%%

%W
W = gA(n+1,aa);

%B
B = W(1:n,1:n);

%b
b = W(1:n,n+1:n+1);

%beta
beta = W(n+1,n+1);

%k
k = 1 / (b'*inv(B)*b + 1 - beta);

%Ar
Ar = k * B;

%Delta
Delta = A - Ar;

%c
c = -inv(B)*b;

%delta
delta = cc - c;

%d0
d0 = 1 / sqrt(det(A));

%d1 
d1 = 1 / sqrt(det(Ar));

%eps
eps = 100 * (d1 - d0) / d0;


disp(['d0 = ' num2str(d0)])
disp(['d1 = ' num2str(d1)])
disp(['eps = ' num2str(eps)])


%% Ãðàôèêè
j = 1:1:m;
fi = 0;

hold on
plot(M(j,1),M(j,2),'+r');
plot(R(A,fi) * cos(fi) + cc(1), R(A,fi) * sin(fi) + cc(2),'ob');
plot(rr(Ar,fi) * cos(fi) + cc(1), rr(Ar,fi) * sin(fi) + cc(2),'-k');