Azzera filtri
Azzera filtri

issue with fmincon function.

4 visualizzazioni (ultimi 30 giorni)
Nasrin
Nasrin il 25 Lug 2024
Commentato: Nasrin il 25 Lug 2024
I want to maximize rate using fmincon function. the problem is I have an error as "Arrays have incompatible sizes for this operation." and "Caused by: Failure in initial objective function evaluation. FMINCON cannot continue.". I provide a part of my code and the objective function here. I guess the problem is in the reshape part and it might be wrong. I would really appreciate it if someone could help me to write the objective function correctly.
when . I should mention that is a vector and is a matrix for .
here is my code.
%initial guess for P2 and W2
P2_0 = ones(Ns,1) * (tilde_P2 / Ns) / 2; % initial guess for P2
W2_0 = ones(Ms,Ns); %initial guess for W2
%combine initial guess in a single vector
x0 = [P2_0 ; W2_0(:)];
% Define the linear constraint (C2)
A = zeros(Np, Ns + Ms * Ns);
for l = 1:Np
A(l, 1:Ns) = xi'; % Place xi' in the first Ns columns of each row
end
% Define the bounds(C1)
lb = [zeros(Ns, 1); -Inf(Ms * Ns, 1)]; % Lower bound (P2 >= 0, no lower bound for W2)
ub = [repmat(tilde_P2 / Ns, Ns, 1); Inf(Ms * Ns, 1)]; % Upper bound (P2 <= tilde_P2_max / Ns, no upper bound for W2)
%%%objective function
obj_fun = @(x) -((Nup - Nt) / Nup * sum(log(1 + x(1:Ns) ./ (sigma2 * reshape(x(Ns+1:end), Ms, Ns).^2))));
% Set options for fmincon
options = optimoptions('fmincon', 'Display', 'iter', 'Algorithm', 'interior-point');
% Run the optimization
[x_opt, fval] = fmincon(obj_fun, x0, A, hat_zeta1, [], [], lb, ub, [], options);
% Extract optimized values for P2 and W2
P2_opt = x_opt(1:Ns);
W2_opt = reshape(x_opt(Ns+1:end), Ms, Ns);
  1 Commento
Nasrin
Nasrin il 25 Lug 2024
@Torstensure. it is long so I write the objective function.
Here is the excutable code.
%parameters:
Mp = 32; %number of antennas at PBS
Ms = 20; %number of antennas at SBS
Np = 4; %number of primary users
Ns = 4; %number of secondary users
dp = 1000; %distance parameter in meters
ds = 500; %distance parameter in meters
tilde_P2_dB = 20; %transmit power in dBW
tilde_P1_dB = 35; %example for tilde_P1 in dBW
p1_max_dB = tilde_P1_dB / Np; % example for p1_max in dBW
eta = 0.5; %parameter eta
rho = 0.9; %parameter rho
Nt = 12; %number of training symbols
Nup = 200; %number of uplink symbols
Ndp = 200; %number of downlink symbols
sigma2_dBm = -104; %noise power in dBm
a_dB = -137; %path loss at reference distance in dB
path_loss_exp = 3.5; %path loss exponent
d0 = 1000; %reference distance in meter
%convert dBW and dBm to linear scale
tilde_P2 = 10^(tilde_P2_dB / 10);
tilde_P1 = 10 .^ (tilde_P1_dB / 10 );
sigma2 = 10^(sigma2_dBm / 10 -3 ); %dBm to mW then to W
a = 10 ^ (a_dB / 10);
%power parameters:
alpha = tilde_P2 / Ns; %parameter alpha
p2_max = tilde_P2 /Ns;
p1_max = tilde_P1 / Np; % example for p1_max
%%%%%%%%%%%%%%%%%%%%Channel modeling%%%%%%%%%%%%%%%%%%%%%%
% Geometry of the system
PBS_location = [0, 0]; %PBS location
SBS_location = [(-dp+ds)/2, (dp-ds)/2]; %SBS location
%generate distances
distance_pu_pbs = dp * rand(Np, 1); %distances of PUs to PBS
distance_su_pbs = dp * rand(Ns, 1); %distances of SUs to PBS
distance_su_sbs = ds * rand(Ns,1); %distances of SUs to SBS
distance_pu_sbs = ds * rand(Ns, 1); %distances of PUs to SBS
%calculate variance for the chanel model
variance_pu_pbs = a * (distance_pu_pbs/ d0).^(-path_loss_exp);
variance_su_pbs = a * (distance_su_pbs/ d0).^(-path_loss_exp);
variance_su_sbs = a * (distance_su_sbs/ d0).^(-path_loss_exp);
variance_pu_sbs = a * (distance_pu_sbs/ d0).^(-path_loss_exp);
%generate chanel matrices
H_1 = (randn(Mp, Np) +1i *randn(Mp, Np)) .*sqrt(variance_pu_pbs.') / sqrt(2); %matrix 32*4
H_12 = (randn(Mp, Ns) +1i *randn(Mp, Ns)) .*sqrt(variance_su_pbs.') / sqrt(2); %matrix 32*4
H_2 = (randn(Ms, Ns) +1i *randn(Ms, Ns)) .*sqrt(variance_su_sbs.') / sqrt(2); %matrix 20*4
H_21 = (randn(Ms, Np) +1i *randn(Ms, Np)) .*sqrt(variance_pu_sbs.') / sqrt(2); %matrix 20*4
%%%%%%%%%%% mu_l %%%%%%%%%%%%%%%%%%%
ZFBF_mu_l = H_1' * H_1; % matrix 4*4
[ZFBF_mu_l_inv, ZFBF_mu_l_pinv] = compute_inverse_pseudoinverse(ZFBF_mu_l);
mu_l_matrix = H_1 * ZFBF_mu_l_inv; % matrix 32*4
mu_l = zeros(Np, 1); %matrix 4*1
for l = 1:Np
column_l = mu_l_matrix(:, l);
mu_l(l) = (norm(column_l))^(-1);
end
D1 = diag(mu_l); %matrix 4*4
W1 = mu_l_matrix * D1; % matrix 32*4
%%%%%%%%%%%%%%%%%%%%%% r1_l_max %%%%%%%%%%%%%%%%%
r1_max = zeros(Np,1); %matrix Np*1
for l = 1:Np
r1_l_max = log(1 + ((mu_l(l)^2 * p1_max) / sigma2));
r1_max(l) = r1_l_max;
end
r1 = rho * r1_max;
%%%%%%%%%%%%%%%% waterfilling for tilde_p1_l %%%%%%%%%%%%%%%
%computing phi
phi_l = mu_l .^2 / sigma2;
%%computing psi
% Objective function for psi
obj_fun = @(psi) sum(max(0, psi - 1 ./ phi_l)) - tilde_P1;
% Initial guess for psi
psi_init = tilde_P1 / Np;
% Options
options = optimoptions('fsolve', 'Display', 'off');
% Solve psi
psi = fsolve(obj_fun, psi_init, options);
%tilde_p1_l
tilde_p1_l = max(0, psi - 1./phi_l);
%%%%%%%%%%%%%%%%%%%%% xi %%%%%%%%%%%%%%%%%%%%%%%%%%%%%
xi = zeros(Ns,1); % matrix Ns*1
for k = 1:Ns
xi_k = 0;
for l = 1:Np
w1_l = W1(:, l);
h12_k = H_12(:, k);
% Compute the inner product and absolute value squared
abs_val_squared = abs(w1_l' * h12_k)^2;
xi_k = xi_k + tilde_p1_l(l) * abs_val_squared;
end
xi(k) = xi_k ./ sum(tilde_p1_l);% matrix 4*1
end
%%%%%%%%%%%%%% hat_pi %%%%%%%%%%%%%%%%%%%
hat_pi= zeros(Np,1); %matrix Np*1
for l = 1:Np
hat_pi_l = Nt * log(1 + (mu_l(l)^2 * p1_max) / (alpha * sum(xi) + sigma2));
hat_pi(l) = hat_pi_l;
end
%%%%%%%%%%%%%% hat_zeta1 %%%%%%%%%%%%%%%%%%%
hat_zeta1 = zeros(Np,1);% matrix Np*1
for l = 1:Np
num = mu_l(l)^2 * p1_max;
denom = 2 ^((Nup * r1(l) - hat_pi(l))/(Nup - Nt));
hat_zeta1_l = (num / (denom -1)) - sigma2;
hat_zeta1(l) = hat_zeta1_l;
end
%%%%%%%%%%%%%%% Optimization problem %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%initial guess for P2 and W2
P2_0 = ones(Ns,1) * (tilde_P2 / Ns) / 2; % initial guess for P2
W2_0 = ones(Ms,Ns); %initial guess for W2
%combine initial guess in a single vector
x0 = [P2_0 ; W2_0(:)];
% Define the linear constraint (C2)
A = zeros(Np, Ns + Ms * Ns);
for l = 1:Np
A(l, 1:Ns) = xi'; % Place xi' in the first Ns columns of each row
end
% Define the bounds(C1)
lb = [zeros(Ns, 1); -Inf(Ms * Ns, 1)]; % Lower bound (P2 >= 0, no lower bound for W2)
ub = [repmat(tilde_P2 / Ns, Ns, 1); Inf(Ms * Ns, 1)]; % Upper bound (P2 <= tilde_P2_max / Ns, no upper bound for W2)
%%%objective function
obj_fun = @(x) -((Nup - Nt) / Nup * sum(log(1 + x(1:Ns) ./ (sigma2 * reshape(x(Ns+1:end), Ms, Ns).^2))));
% Set options for fmincon
options = optimoptions('fmincon', 'Display', 'iter', 'Algorithm', 'interior-point');
% Run the optimization
[x_opt, fval] = fmincon(obj_fun, x0, A, hat_zeta1, [], [], lb, ub, [], options);
Arrays have incompatible sizes for this operation.

Error in solution>@(x)-((Nup-Nt)/Nup*sum(log(1+x(1:Ns)./(sigma2*reshape(x(Ns+1:end),Ms,Ns).^2)))) (line 148)
obj_fun = @(x) -((Nup - Nt) / Nup * sum(log(1 + x(1:Ns) ./ (sigma2 * reshape(x(Ns+1:end), Ms, Ns).^2))));

Error in objfunEvaluator (line 5)
fval = feval(Objfun, x, self.FunArgs.AdditionalParameters{:});

Error in OptimFunctions/objective (line 271)
[fval_, fgrad_, hess_] = self.ObjectiveFunAndGrad(self,self.FunFcn{3},...

Error in OptimFunctions/objective_first_eval (line 614)
[fval,self] = self.objective(X0);

Error in fmincon (line 500)
[initVals.f,initVals.g,HESSIAN,funObj] = funObj.objective_first_eval(X);

Caused by:
Failure in initial objective function evaluation. FMINCON cannot continue.
% Extract optimized values for P2 and W2
P2_opt = x_opt(1:Ns);
W2_opt = reshape(x_opt(Ns+1:end), Ms, Ns);
% Display the result
disp('Optimized values of P2:');
disp(P2_opt);
disp('Optimized values of W2:');
disp(W2_opt);
disp('Optimized objective value:');
disp(-fval);
% Compute the inverse and pseudoinverse of a matrix
function [matrix_inv, matrix_pinv] = compute_inverse_pseudoinverse(matrix)
% Initialize the inverse as empty, to be used if matrix is not invertible
matrix_inv = [];
% Check if the matrix is square
[rows, cols] = size(matrix);
if rows == cols
% Check if the matrix is invertible by determining its determinant
if det(matrix) ~= 0
% Compute the inverse of the matrix
matrix_inv = inv(matrix);
else
warning('Matrix is singular and cannot be inverted. Only pseudoinverse will be returned.');
end
else
warning('Matrix is not square. Only pseudoinverse will be returned.');
end
% Compute the pseudoinverse of the matrix
matrix_pinv = pinv(matrix);
end

Accedi per commentare.

Risposta accettata

Torsten
Torsten il 25 Lug 2024
Modificato: Torsten il 25 Lug 2024
I guess the objective function should read
obj_fun = @(x) -(Nup - Nt) / Nup * sum(log2(1 + x(1:Ns) ./ (sigma2 * (vecnorm(reshape(x(Ns+1:end), Ms, Ns)).^2).')));
If you still encounter errors, please give us executable code. Above, many variables are undefined.
  1 Commento
Nasrin
Nasrin il 25 Lug 2024
@Torsten It works now. Thank you. you helped me alot.

Accedi per commentare.

Più risposte (0)

Community Treasure Hunt

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

Start Hunting!

Translated by