Using the given code, how can I optimize a matrix to minimize a cost function?

8 visualizaciones (últimos 30 días)
Given S1, (1,K) vector, I want to optimize B (N,M) matrix to minimize the following cost function:
Subject to:
Where:
S2 is (1,K) vector and a function of matrix B.
S2 can be calculated after optimizing matrix B using the following system parameters and equations:
clc;
clear;
% Given system parameters:
N = 2;
K = 4;
M=2;
C_l = 4;
H = [0.1185 0.2811; 0.3550 0.8224; 0.3260 0.9644; 0.5333 0.6083]; % 4*2 matrix
A = [-2 1; -1 1]; % 2*2 matrix
C = [7 -3; 7 -3; -2 1; -2 1]; % 4*2 matrix
P = [25000000 0; 0 25000000]; % 4*4 matrix
S1 = [3.1683 3.1686 1.8716 1.8898]; % 1*4 vector
S2 = zeros(1,K); % intial value
B = zeros(N,M); % intial value
% How can we optimize the value of B matrix to achieve our goal?
%calculate S2 from B and the other given inputs
for j=1:1:N
d(j) = (B(j,:)*P*B(j,:)')/((2^(2*C_l))-(norm(A(:,j))^2));
end
D_d = diag(d);
for i=1:1:K
V_d(i)=C(i,:)*P*B'*H(i,:)'*inv(1+H(i,:)*(A'*D_d*A+B*P*B')*H(i,:)');
sigma_d(i)=norm((V_d(i)*H(i,:)*B-C(i,:))*(P^(1/2)))^2+(V_d(i)^2)*(1+H(i,:)*A'*D_d*A*H(i,:)');
S2(i)=0.5*log2((P(1,1))/sigma_d(:,i));
end
  6 comentarios
Rik
Rik el 18 de Nov. de 2020
I edited your code. I implemented your constraint inside the cost function.
I am using fminsearch, which can be sensitive to a local minimum depending on the initial guess. In this case that doesn't matter, because it can't find any valid solution. Are you sure it exists? Are you sure you did not make any mistake when implementing the mathematics in Matlab? Either way, the code below could be helpful.
clc;clear;
% Given system parameters:
N = 2;
K = 4;
M=2;
C_l = 4;
H = [0.1185 0.2811; 0.3550 0.8224; 0.3260 0.9644; 0.5333 0.6083]; % 4*2 matrix
A = [-2 1; -1 1]; % 2*2 matrix
C = [7 -3; 7 -3; -2 1; -2 1]; % 4*2 matrix
P = [25000000 0; 0 25000000]; % 4*4 matrix
S1 = [3.1683 3.1686 1.8716 1.8898]; % 1*4 vector
enforce_S_constraint=true;
% Create a struct so the parameters can be easily captured by the anonymous function
params=struct;
[params.N,params.K,params.M,params.C_l,params.H,params.A,params.C,...
params.P,params.S1,params.enforce_S_constraint]=deal(...
N,K,M,C_l,H,A,C,P,S1,enforce_S_constraint);
params_pass1=params;
params_pass1.enforce_S_constraint=false;
% First try without enforcing S2>=S1, then use the result as initial guess
B_initial_guess=zeros(N,M);
B_pass1=fminsearch(@(B) calculate_cost(B,params_pass1),B_initial_guess);
B_pass2=fminsearch(@(B) calculate_cost(B,params ),B_pass1);
% This is only relevant once you find a valid solution
[~,S2]=calculate_cost(B_pass2,params);
function [cost,S2]=calculate_cost(B,params)
%calculate S2 from B and the other given inputs
[N,K,C_l,H,A,C,P,S1,enforce_S_constraint]=deal(...
params.N,params.K,params.C_l,params.H,params.A,params.C,...
params.P,params.S1,params.enforce_S_constraint);
d=zeros(N,1);
for j=1:1:N
d(j) = (B(j,:)*P*B(j,:)')/((2^(2*C_l))-(norm(A(:,j))^2));
end
D_d = diag(d);
V_d=zeros(1,K);
sigma_d=zeros(1,K);
S2 = zeros(1,K); % intial value
for i=1:1:K
V_d(i)=C(i,:)*P*B'*H(i,:)'*inv(1+H(i,:)*(A'*D_d*A+B*P*B')*H(i,:)');
sigma_d(i)=norm((V_d(i)*H(i,:)*B-C(i,:))*(P^(1/2)))^2+(V_d(i)^2)*(1+H(i,:)*A'*D_d*A*H(i,:)');
S2(i)=0.5*log2((P(1,1))/sigma_d(:,i));
end
if enforce_S_constraint && any(S1>S2)
cost=inf;%set the cost to inf if the constraint is not satisfied
else
cost=sum( (S2-S1).^2 );
end
end
Adi Nor
Adi Nor el 20 de Nov. de 2020
@Rik,Thanks for your help. your code helps me

Iniciar sesión para comentar.

Respuesta aceptada

Matt J
Matt J el 18 de Nov. de 2020
Well, the way you would approach it is using fmincon as below. Using your initial guess of zeros(N,M), fmincon is unable to find a solution so, assuming a feasible solution does exist, you will need to choose a less arbitrary initial guess. I urge you to get familiar with the fmincon documentation.
B_Optimal=doOptimization
Converged to an infeasible point. fmincon stopped because the size of the current step is less than the value of the step size tolerance but constraints are not satisfied to within the value of the constraint tolerance.
B_Optimal = 2×2
46.2078 -19.7949 8.1964 -3.4964
function B_Optimal=doOptimization
% Given system parameters:
N = 2;
K = 4;
M=2;
C_l = 4;
H = [0.1185 0.2811; 0.3550 0.8224; 0.3260 0.9644; 0.5333 0.6083]; % 4*2 matrix
A = [-2 1; -1 1]; % 2*2 matrix
C = [7 -3; 7 -3; -2 1; -2 1]; % 4*2 matrix
P = [25000000 0; 0 25000000]; % 4*4 matrix
S1 = [3.1683 3.1686 1.8716 1.8898]; % 1*4 vector
Binitial = zeros(N,M); % intial value
% How can we optimize the value of B matrix to achieve our goal?
B_Optimal=fmincon( @(B)norm( S2(B)-S1).^2, Binitial,[],[],[],[],[],[],@nlcon);
function [c,ceq]=nlcon(B)
ceq=[];
c=S1-S2(B);
end
function out=S2(B)
%calculate S2 from B and the other given inputs
for j=N:-1:1
d(j) = (B(j,:)*P*B(j,:)')/((2^(2*C_l))-(norm(A(:,j))^2));
end
D_d = diag(d);
for i=K:-1:1
V_d(i)=C(i,:)*P*B'*H(i,:)'*inv(1+H(i,:)*(A'*D_d*A+B*P*B')*H(i,:)');
sigma_d(i)=norm((V_d(i)*H(i,:)*B-C(i,:))*(P^(1/2)))^2+(V_d(i)^2)*(1+H(i,:)*A'*D_d*A*H(i,:)');
out(i)=0.5*log2((P(1,1))/sigma_d(:,i));
end
end
end
  2 comentarios
Stephan
Stephan el 18 de Nov. de 2020
Editada: Stephan el 18 de Nov. de 2020
i agree - i tried with ga and also ga was not able to find a feasible solution that satisfied the constraint:
myOptProb
function myOptProb
% Given system parameters:
N = 2;
K = 4;
M = 2;
C_l = 4;
H = [0.1185 0.2811; 0.3550 0.8224; 0.3260 0.9644; 0.5333 0.6083]; % 4*2 matrix
A = [-2 1; -1 1]; % 2*2 matrix
C = [7 -3; 7 -3; -2 1; -2 1]; % 4*2 matrix
P = [25000000 0; 0 25000000]; % 4*4 matrix
S1 = [3.1683 3.1686 1.8716 1.8898]; % 1*4 vector
S2 = zeros(1,K); % intial value
B0 = zeros(N,M); % intial value
% How can we optimize the value of B matrix to achieve our goal?
opts = optimoptions('ga','MaxGenerations',20000,'PopulationSize',5000,'PlotFcn', @gaplotbestf);
[B, fval, exitflag, output] = ga(@fun,numel(B0),[],[],[],[],[],[],@(x)nonlcon(x,S1,S2),opts)
fun(B)
%calculate S2 from B and the other given inputs
function F = fun(B)
B = reshape(B,N,M);
for j=1:1:N
d(j) = (B(j,:)*P*B(j,:)')/((2^(2*C_l))-(norm(A(:,j))^2));
end
D_d = diag(d);
for i=1:1:K
V_d(i)=C(i,:)*P*B'*H(i,:)'*inv(1+H(i,:)*(A'*D_d*A+B*P*B')*H(i,:)');
sigma_d(i)=norm((V_d(i)*H(i,:)*B-C(i,:))*(P^(1/2)))^2+(V_d(i)^2)*(1+H(i,:)*A'*D_d*A*H(i,:)');
S2(i)=0.5*log2((P(1,1))/sigma_d(:,i));
end
F = sum((S2-S1).^2);
end
function [c,ceq] = nonlcon(~,S1,S2)
c = S1(:)-S2(:); % Compute nonlinear inequalities at x.
ceq = []; % Compute nonlinear equalities at x.
end
end
ends up with:
Optimization terminated: no feasible point found.
B =
0.5862 0.1652 1.4282 -1.0264
fval =
54.7745
exitflag =
-2
output =
struct with fields:
problemtype: 'nonlinearconstr'
rngstate: [1×1 struct]
generations: 1
funccount: 252512
message: 'Optimization terminated: no feasible point found.'
maxconstraint: 3.1686
Rik
Rik el 18 de Nov. de 2020
How certain is it that such a B actually exists? It is not my field, so algebraically determining is such a matrix is possible in the first place is beyond my skills. However, I don't see anything to indicate that it would be impossible to find this out.

Iniciar sesión para comentar.

Más respuestas (0)

Categorías

Más información sobre Optimization Toolbox en Help Center y File Exchange.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by