Hello all. I am trying to solve the following algorithm however, MatLab did not return the correct result for it. Can anyone tell me what is the problem with my code?

Screenshot 2022-09-28 160718|423x499

Here, X is the N*3 dimensions array. f_r is the r-th column of the eigenvector of X*X^{T} which makes sure that X*X^{T} is a square matrix. V is the Vector in my code.

Here is the full code:

*close all

clear

clc

rng(‘default’)

q= rng;

tic

% Generating the sample

x = -0.5 + rand(4,1) ;

y = -0.5 + rand(4,1) ;

X4 = [];

X5 = [];

X8 = [];

X9 = [];

for i = 1:2

if ((x(i))^2 + (y(i))^2) <= 0.25

X4 = [X4,x(i)];

X5 = [X5,y(i)];

end

end

for j = 3:4

if ((x(j))^2 + (y(j))^2) <= 0.25

X8 = [X4,x(j)];

X9 = [X5,y(j)];

end

end

X6 = zeros(1,length(X4)) ;

X7 = 1 * ones(1,length(X8));

radius = 0.5 * ones(1,2);

theta = 2*rand(2,1);
X1 = diag(radius) * cos(theta*pi);

X2 = diag(radius) * sin(theta

*pi);*

X3 = 1rand(2,1);

X3 = 1

X1 = vertcat(X1,X4’,X8’);

X2 = vertcat(X2,X5’,X9’);

X3 = vertcat(X3,X6’,X7’);

X = horzcat(X1,X2,X3); % sample

%algorithm 1

S = X - mean(X,2);

S1 = S(:,1);

S2 = S(:,2);

S3 = S(:,3);

R_s = S * S’;

[Vector,Value] = eig(R_s);

Vector = real(Vector);

s = length(X);

resultX = zeros(length(X),length(X));

resultY = zeros(length(X),length(X));

resultZ = zeros(length(X),length(X));

for i = 1:s

parX = (transpose(Vector(:,i)) * X(:,1)) * (Vector(:,i));

resultX(:,i) = parX;

```
parY = (transpose(Vector(:,i)) * X(:,2)) * (Vector(:,i));
resultY(:,i) = parY;
parZ = (transpose(Vector(:,i)) * X(:,3)) * (Vector(:,i));
resultZ(:,i) = parZ;
```

end

result = zeros(length(X),length(X));

C = {};

```
for i = 1:length(Vector)
res = kron(Vector(:,i),transpose(Vector(:,i)));
for j = 1:length(Vector)
result = res*Vector(j,i);
C{j,i} = result;
end
end
```

D = {};

alpha = 1;

res_sum = zeros(length(Vector),length(Vector));

beta = 10000;

n = length(X1);

obj_min = 100;

z_min = [];

for k = 1:n

cvx_begin sdp

```
variable z(n)
W1 = sum(resultX *z,2);
W2 = sum(resultY *z,2);
W3 = sum(resultZ *z,2);
minimize (alpha*( (transpose(S1 - W1)*(S1 - W1)) + (transpose(S2 - W2)*(S2 - W2))+ (transpose(S3 - W3)*(S3 - W3)) ) + beta*z'*z);
subject to
z(k) == 1;
z(k:end,:) <= 1;
if k>1
z(1:k-1,:) <= 1;
end
z >= 0;
for i = 1:length(Vector)
int = zeros(length(Vector),length(Vector));
for j = 1:length(Vector)
int = int + C{i,j}*z(j);
end
D{i,1} = int;
0 <= D{i,1};
end
cvx_end
if cvx_optval < obj_min
obj_min = cvx_optval;
z_min = z;
end
```

end