When I try to solve optimization problems, the MATLAB returned a error, and this situation only happens after solving the convex problem with binary variable.
while(true)
%优化akn
cvx_begin
cvx_solver mosek
variable akn(K,N) binary
expression Efly2(N,1)
expression Ecomp2(K,N)
expression euser2(K,N)
%Efly、Ecomp、euser的表达式
for i=1:N
Efly2(i,1) = tn*(P0*(1+3norm(vnr(i))/Utip^2) + Piv0/norm(vnr(i)) + d0rous0Apower(norm(vnr(i)),3)/2);
end
for k=1:K
for i=1:N
Ecomp2(k,i) = tnKufumax^3akn(k,i);
end
end
for k=1:K
for i=1:N
euser2(k,i) = (tnKkfmax^3 + wknr(k,i)Pknr(k,i))(1-akn(k,i));
end
end
for k=1:K
for i=1:N
gamakn(k,i)=Pknr(k,i)beta0/zs;
end
end
minimize sum(Efly2)+sum(Ecomp2(:))+sum(euser2(:))
subject to
for k=1:K
for i=1:N
0 <= akn(k,i) <=1;
end
end
sum=0;
for k=1:K
for i=1:N
sum = sum + tnfmax(1-akn(k,i))/C + tn*fumax.*akn(k,i)/C;
end
end
sum >= Ik;
for k=1:K
for i=1:N
wknr(k,i).akn(k,i)Blog2(1+gamakn(k,i)/H^2+power(norm(qnr(:,i)-uk(:,k)),2)) >= tnfumax/C;
end
end
cvx_end
aknr=akn;
%优化时间分配和无人机轨迹
cvx_begin
cvx_solver mosek
variable wkn(K,N)
variable vn(N,1)
variable qn(2,N+1)
variable ykn(K,N)
variable un(N,1)
expression Efly(N,1)
expression Ecomp(K,N)
expression euser(K,N)
for i=1:N
Efly(i,1) = tn*(P0*(1+3*norm(vn(i))/Utip^2) + Pi*un(i) + d0*rou*s0*A*pow_pos(norm(vn(i)),3)/2);
end
for k=1:K
for i=1:N
Ecomp(k,i) = tn*Ku*fumax^3*aknr(k,i);
end
end
for k=1:K
for i=1:N
euser(k,i) = (tn*Kk*fmax^3 + wkn(k,i)*Pknr(k,i))*(1-aknr(k,i));
end
end
for k=1:K
for i=1:N
gamakn(k,i)=Pknr(k,i)*beta0/zs;
end
end
minimize sum(Efly)+sum(euser(:))+sum(Ecomp(:))
subject to
for i=1:N
norm(vn(i)) <= Vmax;
end
a = sum(wkn);
for i=1:N
a(1,i) <= tn;
end
for i=1:N
norm(qn(:,i+1)-qn(:,i)) <= 20;
end
qn(:,1) == qI;
qn(:,N+1) == qF;
for k=1:K
wkn(k,N) == 0;
end
for k=1:K
for i=1:N
quad_over_lin(ykn(k,i),wkn(k,i)) <= B*log2(1+gamakn(k,i)./(H^2+power(norm(qnr(:,i)-uk(:,k)),2))) - B.*gamakn(k,i)*log2(exp(1))*(pow_pos(norm(qn(:,i)-uk(:,k)),2)-power(norm(qnr(:,i)-uk(:,k)),2))./((H^2+power(norm(qnr(:,i)-uk(:,k)),2))*(H^2+power(norm(qnr(:,i)-uk(:,k)),2)+gamakn(k,i)));
end
end
for k=1:K
for i=1:(N-1)
square(yknr(k,i)) + 2*yknr(k,i)*(ykn(k,i)-yknr(k,i)) >= tn*fumax/C;
end
end
for i=1:N
square_pos(quad_over_lin(1,un(i))) <= unr(i)^2+2*unr(i)*(un(i)-unr(i))+pow_pos(norm(vnr(i,:)),2)/v0^2+2.*vnr(i,:)'*(vn(i,:)-vnr(i,:))/v0^2;
end
cvx_end
cvx_begin
variable Pkn(K,N)
expression E(K,N)
for k=1:K
for i=1:N
E(k,i)=wknr(k,i)*Pkn(k,i)*(1-aknr(k,i))+tn*Kk*fmax^3;
end
end
minimize sum(E(:))
subject to
for k=1:K
for i=1:(N-1)
wknr(k,i)*aknr(k,i)*B*log2(1+gamakn(k,i)/(H^2+pow_pos(norm(qnr(:,i)-uk(:,k)),2))) >= tn*fumax;
end
end
for k=1:K
for i=1:N
Pkn(k,i) >= 0;%Pmin设置为0
end
end
for k=1:K
for i=1:N
Pkn(k,i) <= 0.05;%Pmax设置为50毫瓦
end
end
cvx_end
Pknr=Pkn;
r=r+1;
if(((sum(Efly)+sum(Ecomp(:))+sum(euser(:)))^(r+1)-(sum(Efly)+sum(Ecomp(:))+sum(euser(:))^r)) < o)
break;
end
end
Please show the simplest reproducible problem which illustrates this error. That means, complete with all input data, and showing all solver and CVX output.
If you are solving a sequence of CVX problems, have you carefully looked at the results of each CVX problem? Does each problem solve to optimality? Did you look at each of the terms in the statement which triggers the error?
Calling Mosek 9.1.9: 4001 variables, 3001 equality constraints
MOSEK Version 9.1.9 (Build date: 2019-11-21 11:34:40)
Copyright © MOSEK ApS, Denmark. WWW: mosek.com
Platform: Windows/64-X86
Problem
Name :
Objective sense : min
Type : LO (linear optimization problem)
Constraints : 3001
Cones : 0
Scalar variables : 4001
Matrix variables : 0
Integer variables : 1000
Optimizer started.
Mixed integer optimizer started.
Threads used: 16
Presolve started.
Presolve terminated. Time = 0.00
Presolved problem: 0 variables, 0 constraints, 0 non-zeros
Presolved problem: 0 general integer, 0 binary, 0 continuous
Clique table size: 0
BRANCHES RELAXS ACT_NDS DEPTH BEST_INT_OBJ BEST_RELAX_OBJ REL_GAP(%) TIME
0 0 1 0 9.7986500000e+03 NA NA 0.1
0 1 1 0 9.7986500000e+03 9.7986500000e+03 0.00e+00 0.1
An optimal solution satisfying the relative gap tolerance of 1.00e-02(%) has been located.
The relative gap is 0.00e+00(%).
An optimal solution satisfying the absolute gap tolerance of 0.00e+00 has been located.
The absolute gap is 0.00e+00.
Objective of best integer solution : 9.798650000000e+03
Best objective bound : 9.798650000000e+03
Construct solution objective : Not employed
User objective cut value : Not employed
Number of cuts generated : 0
Number of branches : 0
Number of relaxations solved : 1
Number of interior point iterations: 0
Number of simplex iterations : 0
Time spend presolving the root : 0.00
Time spend optimizing the root : 0.01
Mixed integer optimizer terminated. Time: 0.13
Optimizer terminated. Time: 0.19
Integer solution solution summary
Problem status : PRIMAL_FEASIBLE
Solution status : INTEGER_OPTIMAL
Primal. obj: 9.7986500000e+03 nrm: 3e+09 Viol. con: 0e+00 var: 0e+00 itg: 0e+00
Optimizer summary
Optimizer - time: 0.19
Interior-point - iterations : 0 time: 0.00
Basis identification - time: 0.00
Primal - iterations : 0 time: 0.00
Dual - iterations : 0 time: 0.00
Clean primal - iterations : 0 time: 0.00
Clean dual - iterations : 0 time: 0.00
Simplex - time: 0.00
Primal simplex - iterations : 0 time: 0.00
Dual simplex - iterations : 0 time: 0.00
Mixed integer - relaxations: 1 time: 0.13
Status: Solved
Optimal value (cvx_optval): +297924
没有为类 ‘cvx’ 的值定义函数 ‘subsindex’。
出错 minimize (line 14)
x = evalin( ‘caller’, sprintf( '%s ', varargin{:} ) );
出错 Untitled1 (line 126)
minimize sum(Efly)+sum(euser(:))+sum(Ecomp(:))
I am not finding the sequence of events to be clear. You hvae not provided input data to make this problem reproducible. How many CVX problems have been solved prior to the one in which the error occurs? Please insert the following just prior to the statement triggering the error. CVX will display what CVX considers these terms to be.
Efly
euser
Ecomp
You should follow the principle of simplifying when debugging. What you now is just a confusing non-reproducible mess.
Also, this may not trigger an error message, but could cause your program to not work as intended. After cvx_end, CVX expressions are not guaranteed to be populated with their optimal values. CVX variables are populated with their optimal values. So in order to be assured of having optimal values of expressions, you must reconstruct them after cvx _end, starting with using CVX variable values (i.e., reeexecute the appropriate lines of code within your CVX program, except without declaring expressions, and you will be using the numerically populated values of what were the CVX variables, but which are not just double precision MATLAB variables.
The complete code is as follows：
clear all;
B=1e6;H=100;Ik=1e8;C=1000;T=100;fmax=3e8;fumax=6e9;
qI=[-700;0];qF=[700;0];Vmax=40;Kk=1e-28;Ku=1e-28;
rou=1.225;A=0.603;Utip=200;d0=0.301;v0=6.2089;
s0=0.0499;P0=225.79;Pi=426.07;tn=0.5;N=200;
beta0=-60;zs=-100;K=5;o=1e-4;r=0;
Q01={0,1};
%地面设备的位置
uk=[-510,-400,0,480,750;
20,300,250,600,110];
%定义优化变量的初始值，初始轨迹为一半圆
vn=ones(N,1);
for i=1:N
vn(i)=40;
end
qn=zeros(2,N+1);
R=700;
theta=pi:-pi/(N):0;
qn(1,:)=R*cos(theta);
qn(2,:)=R*sin(theta);
Pkn=ones(K,N);
ykn=ones(K,N);
un=ones(N,1);
wkn=ones(K,N);
akn=ones(K,N);
%关于vnr、unr、qnr、yknr、wknr、Pknr、fknr、fuknr(考虑到有的参数在分母位置，分母不能为0，所以用ones)
vnr=ones(N,1);
aknr=ones(K,N);
unr=ones(N,1);
qnr=ones(2,N+1);
yknr=ones(K,N);
wknr=ones(K,N);
Pknr=ones(K,N);
while(true)
%优化akn
cvx_begin
cvx_solver mosek
variable akn(K,N) binary
expression Efly2(N,1)
expression Ecomp2(K,N)
expression euser2(K,N)
%Efly、Ecomp、euser的表达式
for i=1:N
Efly2(i,1) = tn*(P0*(1+3*norm(vnr(i))/Utip^2) + Pi*v0/norm(vnr(i)) + d0*rou*s0*A*power(norm(vnr(i)),3)/2);
end
for k=1:K
for i=1:N
Ecomp2(k,i) = tn*Ku*fumax^3*akn(k,i);
end
end
for k=1:K
for i=1:N
euser2(k,i) = (tn*Kk*fmax^3 + wknr(k,i)*Pknr(k,i))*(1-akn(k,i));
end
end
for k=1:K
for i=1:N
gamakn(k,i)=Pknr(k,i)*beta0/zs;
end
end
minimize sum(Efly2)+sum(Ecomp2(:))+sum(euser2(:))
subject to
for k=1:K
for i=1:N
0 <= akn(k,i) <=1;
end
end
sum=0;
for k=1:K
for i=1:N
sum = sum + tn*fmax*(1-akn(k,i))/C + tn*fumax.*akn(k,i)/C;
end
end
sum >= Ik;
for k=1:K
for i=1:N
wknr(k,i).*akn(k,i)*B*log2(1+gamakn(k,i)/H^2+power(norm(qnr(:,i)-uk(:,k)),2)) >= tn*fumax/C;
end
end
cvx_end
aknr=akn;
%优化时间分配和无人机轨迹
cvx_begin
cvx_solver mosek
variable wkn(K,N)
variable vn(N,1)
variable qn(2,N+1)
variable ykn(K,N)
variable un(N,1)
expression Efly(N,1)
expression Ecomp(K,N)
expression euser(K,N)
for i=1:N
Efly(i,1) = tn*(P0*(1+3*norm(vn(i))/Utip^2) + Pi*un(i) + d0*rou*s0*A*pow_pos(norm(vn(i)),3)/2);
end
for k=1:K
for i=1:N
Ecomp(k,i) = tn*Ku*fumax^3*aknr(k,i);
end
end
for k=1:K
for i=1:N
euser(k,i) = (tn*Kk*fmax^3 + wkn(k,i)*Pknr(k,i))*(1-aknr(k,i));
end
end
for k=1:K
for i=1:N
gamakn(k,i)=Pknr(k,i)*beta0/zs;
end
end
minimize sum(Efly)+sum(euser(:))+sum(Ecomp(:))
subject to
for i=1:N
norm(vn(i)) <= Vmax;
end
a = sum(wkn);
for i=1:N
a(1,i) <= tn;
end
for i=1:N
norm(qn(:,i+1)-qn(:,i)) <= 20;
end
qn(:,1) == qI;
qn(:,N+1) == qF;
for k=1:K
wkn(k,N) == 0;
end
for k=1:K
for i=1:N
quad_over_lin(ykn(k,i),wkn(k,i)) <= B*log2(1+gamakn(k,i)./(H^2+power(norm(qnr(:,i)-uk(:,k)),2))) - B.*gamakn(k,i)*log2(exp(1))*(pow_pos(norm(qn(:,i)-uk(:,k)),2)-power(norm(qnr(:,i)-uk(:,k)),2))./((H^2+power(norm(qnr(:,i)-uk(:,k)),2))*(H^2+power(norm(qnr(:,i)-uk(:,k)),2)+gamakn(k,i)));
end
end
for k=1:K
for i=1:(N-1)
square(yknr(k,i)) + 2*yknr(k,i)*(ykn(k,i)-yknr(k,i)) >= tn*fumax/C;
end
end
for i=1:N
square_pos(quad_over_lin(1,un(i))) <= unr(i)^2+2*unr(i)*(un(i)-unr(i))+pow_pos(norm(vnr(i,:)),2)/v0^2+2.*vnr(i,:)'*(vn(i,:)-vnr(i,:))/v0^2;
end
cvx_end
cvx_begin
variable Pkn(K,N)
expression E(K,N)
for k=1:K
for i=1:N
E(k,i)=wknr(k,i)*Pkn(k,i)*(1-aknr(k,i))+tn*Kk*fmax^3;
end
end
minimize sum(E(:))
subject to
for k=1:K
for i=1:(N-1)
wknr(k,i)*aknr(k,i)*B*log2(1+gamakn(k,i)/(H^2+pow_pos(norm(qnr(:,i)-uk(:,k)),2))) >= tn*fumax;
end
end
for k=1:K
for i=1:N
Pkn(k,i) >= 0;%Pmin设置为0
end
end
for k=1:K
for i=1:N
Pkn(k,i) <= 0.05;%Pmax设置为50毫瓦
end
end
cvx_end
Pknr=Pkn;
r=r+1;
if(((sum(Efly)+sum(Ecomp(:))+sum(euser(:)))^(r+1)-(sum(Efly)+sum(Ecomp(:))+sum(euser(:))^r)) < o)
break;
end
endI get this:at the same location:
Efly =
cvx convex expression (200x1 vector)
euser =
cvx zero expression (5x200 matrix)
Ecomp =
cvx positive constant expression (5x200 matrix)
Unable to use a value of type cvx as an index.
I don’t know what is going on.
I will also point out that the numerical scaling is horrible and should be fixed. Also note that constraining binary to be 0 <= binary_variable <= 1 is redundant, but I don’t thin it causes harm.
Thanks for your help！
I’m new to CVX coding, and I really appreciate your patience!!!
Even
sum(Efly)
produces this error
Unable to use a value of type cvx as an index.
I don’t know why. Perhaps it is a bug in CVX.
Anyway，thanks a lot !!!
I discovered that CVX does process
ones(1,N)*Efly
without error. That is equivalent to sum(Efly).
The last 2 terms in
minimize sum(Efly2)+sum(Ecomp2(:))+sum(euser2(:))
are constants; so they can be removed from the objective (and added in after CVX completes if the “correct” optimal objective value is needed).
So the objective can be replaced by
minimize ones(1,N)*Efly
which is processed without error.
But then the same error occurs later in the program. So my workaround worked in one place, but whatever the problem is still manifested itself later.