**You are now following this question**

- You will see updates in your followed content feed.
- You may receive emails, depending on your communication preferences.

# How to solve "Conversion to logical from optim.problemdef.OptimizationEquality is not possible." ?

17 views (last 30 days)

Show older comments

Hello!

Hi, I have meet the problem in MATLAB using problem based optimization, my optimization variable is q_int , i have used to generate a vector "qu", but an error displays

"Conversion to logical from optim.problemdef.OptimizationEquality is not possible.

Error in Generate_position (line 4)

if q_int == 0"

a = floor(L/UAV_Speed*dt);

q_int = optimvar("q_int","LowerBound",0,"UpperBound",a);

function [qu] = Generate_position(Flight_direction,q_int,UAV_Speed,dt,L)

direction = Flight_direction;

if q_int == 0 % here the error

direction = 1;

end

if q_int == L

direction = -1;

end

for k=2:T

qu(k) = qu(k-1) + direction*UAV_Speed*dt;

if qu(k)==0 || qu(k)==L

direction = -direction;

end

end

Appreciate for your help!

##### 8 Comments

Jan
on 29 Jan 2023

Torsten
on 29 Jan 2023

You didn't initialize qu(1) in your function. This will be another reason it will error.

Matt J
on 29 Jan 2023

Edited: Matt J
on 29 Jan 2023

The error is produced by the part of the code that invokes Generate_position, which you haven't shown us. So, we can only guess as to how it is used.

If Generate_position is the objective function, it doesn't look like a reasonable one, since qu is a discontinuous function of q_int.

Maria
on 29 Jan 2023

q_int is the optimized variable should have a value between 0 and a , but there is a condition that it should be put in qu_int because it related to all the vector qu

here is the function that used the vector qu to calculate a formula, and i gave above the code to generate qu

qu = Generate_position(T,Flight_direction,q_int,UAV_Speed,dt,L);

arg1 = argument1(N,T,BP,H,qu,qc);

function [arg1] = argument1(N,T,BP,H,qu,qc)

M = size(N,1);

arg1 = optimexpr(M,T);

for i=1:size(N,1)

for k=1:T

if N(i,k)==1 && BP(i,k)~=0

arg1(i,k) = BP(i,k).*(H^2+(norm((qu(k))-qc(i,k))^2));

end

end

end

end

Torsten
on 29 Jan 2023

i initialize qu(1) = q_int*UAV_Speed*dt, but the same error is displayed,

I know. That's why I wrote: This will be another reason it will error.

Maria
on 29 Jan 2023

@Torsten the problem is when i change q_int with such value , it works. That's why i think that the problem in the initialization of my optimization variable q_int , but i see in documents that we wrote variables as shown, or there is something missing?

q_int = optimvar("q_int","LowerBound",0,"UpperBound",a);

### Accepted Answer

Matt J
on 29 Jan 2023

Edited: Matt J
on 29 Jan 2023

This will solve the "Conversion to logical" error, but not all of the other issues raised by Torsten and me. In particular, Generate_position will still fail if q_int has any value other than 0 or L.

q_int = optimvar("q_int","LowerBound",0,"UpperBound",a);

fun=fcn2optimexpr(@(INPUT)computeThing(INPUT, T,Flight_direction,UAV_Speed,dt,L, N,T,BP,H,qc),...

q_int );

function arg1=computeThing(q_int, T,Flight_direction,UAV_Speed,dt,L, N,T,BP,H,qc)

qu = Generate_position(T,Flight_direction,q_int,UAV_Speed,dt,L);

arg1 = argument1(N,T,BP,H,qu,qc);

end

##### 122 Comments

Torsten
on 29 Jan 2023

Edited: Torsten
on 30 Jan 2023

There is an error in the code I supplied. You will have to use

a = floor(L/(UAV_Speed*dt));

q_int = optimvar('q_int','Type','Integer','LowerBound',0,'UpperBound',a);

fun = fcn2optimexpr(@(INPUT)computeThing(INPUT,T,Flight_direction,UAV_Speed,dt,L,N,T,BP,H,qc),...

q_int);

function arg1 = computeThing(q_int,T,Flight_direction,UAV_Speed,dt,L,N,T,BP,H,qc)

qu = Generate_position(T,Flight_direction,q_int,UAV_Speed,dt,L,T);

arg1 = argument1(N,T,BP,H,qu,qc);

end

function [qu] = Generate_position(Flight_direction,q_int,UAV_Speed,dt,L,T)

direction = Flight_direction;

if q_int == 0

direction = 1;

end

if q_int == floor(L/(UAV_Speed*dt));

direction = -1;

end

qu(1) = q_int*UAV_Speed*dt;

for k=2:T

qu(k) = qu(k-1) + direction*UAV_Speed*dt;

if qu(k)==0 || qu(k)==L

direction = -direction;

end

end

end

function [arg1] = argument1(N,T,BP,H,qu,qc)

M = size(N,1);

arg1 = optimexpr(M,T);

for i=1:size(N,1)

for k=1:T

if N(i,k)==1 && BP(i,k)~=0

arg1(i,k) = BP(i,k).*(H^2+(norm((qu(k))-qc(i,k))^2));

end

end

end

end

See whether Matt's code suggestion works as required.

Maria
on 29 Jan 2023

Error in fcn2optimexpr

Error in MyoptimizationProblem1 (line 71)

fun=fcn2optimexpr(@(INPUT)computeThing(INPUT, T,Flight_direction,UAV_Speed,dt,L, N,T,BP,H,qc),...

Caused by:

Function evaluation failed while attempting to determine output size. The function might contain an error, or might not be well-defined at the automatically-chosen point. To specify output size without function evaluation, use 'OutputSize'.

Matt J
on 29 Jan 2023

We can run that line right here in the forum. As you can see, no errors:

a=5;

q_int = optimvar("q_int","LowerBound",0,"UpperBound",a),

q_int =

OptimizationVariable with properties:
Name: 'q_int'
Type: 'continuous'
IndexNames: {{} {}}
LowerBound: 0
UpperBound: 5
See variables with show.
See bounds with showbounds.

Maria
on 30 Jan 2023

@Torsten i found a solution for the error , but i didn't know why q_int is 0, i ran the code many times and it doesn't change

Solving problem using intlinprog.

LP: Optimal objective value is -75.000000.

Optimal solution found.

Intlinprog stopped at the root node because the objective value is within a gap tolerance of the optimal value,

options.AbsoluteGapTolerance = 0 (the default value). The intcon variables are integer within tolerance,

options.IntegerTolerance = 1e-05 (the default value).

sol =

struct with fields:

S: [96×90 double]

q_int: 0

fval =

75

exitflag =

OptimalSolution

Torsten
on 30 Jan 2023

Why is S a double array ? I thought it should be a binary integer array.

And do you have any reason to believe that q_int should be different from 0 ?

Matt J
on 30 Jan 2023

so why it doesn't accept a condition "if q_int ==0".

There's no reason to think it wont't. Below is a simple example where it does:

q_int = optimvar("q_int","LowerBound",0,"UpperBound",1);

prob=optimproblem;

prob.Objective=fcn2optimexpr(@entropyFcn, q_int);

sol0.q_int=0.4;

sol=solve(prob,sol0)

Solving problem using fmincon.
Local minimum found that satisfies the constraints.
Optimization completed because the objective function is non-decreasing in
feasible directions, to within the value of the optimality tolerance,
and constraints are satisfied to within the value of the constraint tolerance.

sol = struct with fields:

q_int: 0.3679

function out=entropyFcn(q_int)

if q_int==0

out=0;

else

out= q_int*log(q_int);

end

end

Maria
on 30 Jan 2023

Edited: Maria
on 30 Jan 2023

Torsten
on 30 Jan 2023

Edited: Torsten
on 30 Jan 2023

But if S is a solution matrix, you somewhere define it as you define q_int with the "Type", "Integer" , lower and upper bounds options, don't you ?

Concerning q_int: I don't know why it doesn't change its value. Did you try to initialize it to an integer value different from 0 as Matt does it in his example code ?

Matt J
on 30 Jan 2023

Edited: Matt J
on 30 Jan 2023

Torsten
on 30 Jan 2023

Edited: Torsten
on 30 Jan 2023

The value of 'x0' is invalid. Initial point must contain values for variable 'S'.

I think - also with problem-based optimization - you cannot supply an initial value for only a part of your solution variables. So supplying a value only for q_int and not for S is not possible, I guess.

So I think you will additionally have to pass a matrix of 0's and 1's for x0.S (and maybe additional solution variables I don't know of).

Torsten
on 30 Jan 2023

Edited: Torsten
on 30 Jan 2023

In function "Generate_position":

qu = zeros(1,T);

direction = Flight_direction;

if qu(1) == 0

direction = 1;

end

If you set qu = zeros(1,T), especially qu(1) is 0.

Since you define the constraint

constr1 = qu(1) == q_int*UAV_Speed*dt ;

and qu(1) = 0, q_int must be 0.

Why did you change the code to something wrong ?

qu(1) must be set to q_int*UAV_Speed*dt, not to 0, with a variable integer value for q_int in the limits from 0 to floor(L/(UAV_Speed*dt)).

Maria
on 30 Jan 2023

this constraint works only when i initialize the vector qu that's why i set qu = zeros(1,T).

constr1 = qu(1) == q_int*UAV_Speed*dt ;

also when i changed "generation position" the code works well, just for q_int = 0.

i didn't know what's the role of the solver here?

Maria
on 30 Jan 2023

Torsten
on 30 Jan 2023

Edited: Torsten
on 30 Jan 2023

as i know the solver will find the best solution according to such variables, here my optimized variable is q_int and i put it between 0 and a , but it take the result for just one point "0"

Yes, your problem definition is wrong.

The result from your setup will automatically be q_int = 0 because solving your constraint

qu(1) = q_int*UAV_Speed*dt

for q_int gives

q_int = qu(1)/(UAV_Speed*dt) = 0/((UAV_Speed*dt) = 0

So the solver doesn't have the possibility to vary q_int. You force the solver to set q_int = 0.

Torsten
on 30 Jan 2023

Edited: Torsten
on 30 Jan 2023

Did you read about the two approaches ?

If not, read this before proceeding:

The problem-based formulation is transformed internally to a solver-based formulation.

Thus if your problem can be solved as a problem-based formulation, it can also be solved as a solver-based formulation.

Maria
on 30 Jan 2023

Torsten
on 31 Jan 2023

The error message is due to the problem-based approach. I cannot help you here.

But wasn't the problem already solved by Matt's code ?

Maria
on 31 Jan 2023

@Torsten @Matt J i want to use piecewise instead of if statement, but an error was displayed in this case if qu(k)==0 || qu(k)==L

function [qu] = Generate_position(Flight_direction,q_int,UAV_Speed,dt,L,T,a)

direction = Flight_direction;

syms direction(q_int)

direction(q_int) = piecewise(q_int == 0,1, q_int ==a,-1);

% if (q_int == 0)

% direction = 1;

% end

% if q_int == floor(L/(UAV_Speed*dt))

% direction = -1;

%end

qu(1) = q_int*UAV_Speed*dt;

for k=2:T

qu(k) = qu(k-1) + direction*UAV_Speed*dt;

%mask = qu(k)==0 || qu(k)==L;

syms direction(qu)

direction(qu(k)) = piecewise((qu(k) == 0)|| (qu(k)==L),-direction, direction); %here error

% if qu(k)==0 || qu(k)==L

% direction = -direction;

%end

end

end

Matt J
on 31 Jan 2023

But wasn't the problem already solved by Matt's code ?...no it wasn't solved therefore i changed q_int by qu(1) and it takes only 0 and it is not what i want

I have reverted to your original definition of Generate_position, fixed some errors, and applied my suggestion to use fcn2optimexpr(). As you can see below, it produces output for qu successfully.

However, as I pointed out from the very beginning, it is also discontinuous, as you can see from when we evaluate with q_int=L-1e-12 instead of q_int=L. This small change in q_int produces a very different output. You need to make "direction" a continuous and differentiable function of q_int.

lambda = 1/10;

T = 3600; % Mission Time

dt = 1; % controller time step

L = 1000;

Vmin = 50;

Vmax = 100;

UAV_Speed = 20;

%qu_start = 0;

Flight_direction = 1;

H = 80;

G0 = -50;

BW = 20e6;

Pmax = 35;

N0 = -130;

C = 1e3 ;

f = 2e9;

alpha = 6e6;

Beta = 10e-28;

E_prop = 388944;

E_max = 400e3;

Q = 6;

K = 10e8;

a = floor(L/UAV_Speed*dt);

q_int = optimvar("q_int","LowerBound",0,"UpperBound",a);

fcn=@(INPUT)Generate_position(Flight_direction,INPUT,UAV_Speed,dt,L,T);

generatePosition=fcn2optimexpr(fcn,q_int);

s.q_int=0;

qu=evaluate(generatePosition,s); qu(1:5)

ans = 1×5

0 20 40 60 80

s.q_int=L;

qu=evaluate(generatePosition,s); qu(1:5)

ans = 1×5

0 -20 -40 -60 -80

s.q_int=L-1e-12;

qu=evaluate(generatePosition,s); qu(1:5)

ans = 1×5

0 20 40 60 80

function [qu] = Generate_position(Flight_direction,q_int,UAV_Speed,dt,L,T)

direction = Flight_direction;

if q_int == 0 % here the error

direction = 1;

end

if q_int == L

direction = -1;

end

qu = zeros(1,T); %<----Mattt J inserted

for k=2:T

qu(k) = qu(k-1) + direction*UAV_Speed*dt;

if qu(k)==0 || qu(k)==L

direction = -direction;

end

end

end

Torsten
on 31 Jan 2023

Edited: Torsten
on 31 Jan 2023

However, as I pointed out from the very beginning, it is also discontinuous, as you can see from when we evaluate with q_int=L-1e-12 instead of q_int=L. This small change in q_int produces a very different output. You need to make "direction" a continuous and differentiable function of q_int.

"q_int" is an integer optimization variable with value in between 0 and floor(L/(UAV_Speed*dt)) and "Flight_Direction" is another integer optimization variable with value +1 or -1, as far as I understood the optimization problem to be solved.

Matt J
on 31 Jan 2023

Edited: Matt J
on 31 Jan 2023

"q_int" is an integer optimization variable with value in between 0 and

@Torsten Well, we need Maria to confirm that. So far, that condition does not show up anywhere in her remarks or in her code. Noteably, she has deliberately set S to be integer-valued, but not q_int.

Maria
on 31 Jan 2023

Edited: Maria
on 31 Jan 2023

@Matt J @Torsten until now i want to fix Flight_Direction at the begining in 1 or -1 and optimize just q_int , after this i will take the better solution from both.

for more clarification, Flight_Direction is the decision i will take at the start of q_int, it should be fixed one time, but its value will be changed through the condition of 0 and L

Matt J
on 31 Jan 2023

but its value will be changed through the condition of 0 and L

This is critical to understand, and the English here is not clear enough. Should q_int be able to take any value between 0 and L? Or, should it only be able take one of two values, 0 or L?

Maria
on 31 Jan 2023

Maria
on 31 Jan 2023

Edited: Maria
on 31 Jan 2023

here the code works well according to the initialization of q_int , but i'm asking if each time should i change q_int value by myself or it will be done by the solver?

fcn=@(INPUT)Generate_position(Flight_direction,INPUT,UAV_Speed,dt,L,T);

generatePosition=fcn2optimexpr(fcn,q_int);

s.q_int=5;

qu=evaluate(generatePosition,s); qu(1:T)

100 120 140 160 180 200 220 240 260 280 300 320 340 360 380 400 420 440 ....

Torsten
on 31 Jan 2023

A solution variable is a variable the solver should optimize on its own. So interventions from your side are not necessary (and contraproductive).

You can see if something has changed concerning MATLAB's interpretation of your problem if "ga" is chosen as solver instead of "intlinprog". Your problem cannot be solved by "intlinprog".

Torsten
on 31 Jan 2023

Edited: Torsten
on 31 Jan 2023

If you use "solve" within the problem-based optimization approach, the suitable MATLAB solver is chosen automatically by the software (depending on the "difficulty" of your problem).

As you can see from the output so far:

Solving problem using intlinprog.

LP: Optimal objective value is -75.000000.

Optimal solution found.

Intlinprog stopped at the root node because the objective value is within a gap tolerance of the optimal value,

options.AbsoluteGapTolerance = 0 (the default value). The intcon variables are integer within tolerance,

options.IntegerTolerance = 1e-05 (the default value).

sol =

struct with fields:

S: [96×90 double]

q_int: 0

fval =

75

exitflag =

OptimalSolution

MATLAB chose "intlinprog" as solver because your problem formulation did not meet the formulation of your "real" problem.

If

Solving problem using ga.

appears here, you've made a step forward.

Matt J
on 31 Jan 2023

Attached is my rewrite of your code. It runs and produces non-trivial output. Unfortunately, you have a problem. Your problem is a MINLP with over 1 million integer constrained variables. The only algorithms MATLAB has to deal with MINLP's do not handle that many unknown integers without exceeding memory limits. So, the demo below was only possible when I reduce T to something much smaller, in this case 36. You will probably need to consider reformulation of the problem.

[sol,fval,exitflag] =optimizeIt()

Solving problem using ga.
Optimization terminated: average change in the penalty fitness value less than options.FunctionTolerance
and constraint violation is less than options.ConstraintTolerance.

sol = struct with fields:

S: [4×36 double]
q_int: 472

fval = 3

exitflag =

SolverConvergedSuccessfully

Maria
on 31 Jan 2023

@Torsten here it solved but q_int doesn't considered

Solving problem using ga.

Optimization terminated: average change in the penalty fitness value less than options.FunctionTolerance

and constraint violation is less than options.ConstraintTolerance.

sol =

struct with fields:

S: [16×90 double]

fval =

10

Unrecognized field name "q_int".

Error in MyoptimizationProblem1 (line 108)

sol.q_int

Torsten
on 31 Jan 2023

Maria
on 31 Jan 2023

@Matt J " Loop over all possible q " yes this is what i need exactly , each time S will be solved for fixed q_int and extract the result, and at final it display the best q_int that gives the max of S.

i already did that without using optimization, i have tested each value that q_int could take and calculated S, but now i need to use an optimizer to solve this.

Matt J
on 31 Jan 2023

i already did that without using optimization

What you describe can't be done without optimization. Even for fixed q_int, the solution for S is a non-trivial integer linear program.

In any case, I attached a solution that jointly optimizes in S and q_int, so you should check that out. However, I don't think the joint optimization can be scaled to large T.

Maria
on 31 Jan 2023

@Matt J Do you mean each row S(i,:) of S should contain only a single 1 and the rest should be zero?Yes exactly

i have N is the initial matrix that contains also a single 1 by row, so if this "1" satisfied the two constraints so S will take "1" otherwise "0"

That's why in my code i add the condition "if N(i,j)==1"

Matt J
on 31 Jan 2023

Edited: Matt J
on 31 Jan 2023

@Maria What you say implies that an optimal solution exists such that S(i,j)=0 if N(i,j)=0. However, that does not mean it is the unique optimal solution. There can be other optimal solutions which violate this, so long as (C6) and (C7) are satisfied. You may be accustomed to seeing simpler solutions for S because linear programming tends to find sparse solutions. However, now that we are jointly optimizing over both S and q_int, a linear programming algorithm can't be used.

Incidentally, your objective function implementation is unnecessarily complicated. Since N is a binary matrix, then the following is equivalent to the single expression f=S(:).'*N(:). You should probably switch to this, as it is much more efficient. In fact, the other double loops in your functions can also be easily vectorized.

function [f] = Objective_function(S)

f=zeros(M,T);

for i = 1:M

for j = 1:T

if N(i,j) == 1

f1(i,j) = S(i,j).*N(i,j);

end

end

end

f = sum(f1,"all");

end

Torsten
on 31 Jan 2023

Torsten
on 31 Jan 2023

in my objective function i can put just sum of S

Not according to your objective function definition.

Matt J
on 31 Jan 2023

Edited: Matt J
on 31 Jan 2023

I think you mean f=sum(S.*N,'all'). That's true, but f=S(:).'*N(:) might be a bit faster.

Incidentally also, if you know that the only S(i,j) that matter are the ones where N(i,j)=1, then you really have only M unknown S values, whereas you have currently formulated the problem with M*T unknown values. This greatly increases the burden on the solver. You probably should remove from the problem the S(i,j) you already know can be set to 0. That would probably make it so that ga() can be scaled up to larger T. However, I still think a loop over q_int may be more reliable.

Maria
on 31 Jan 2023

@Matt J i put the condition N(i,j)=1, in order to consider just ones in N and based on it S will be built.

M is unknown , it is generated randomly so its value can be changed but T is fix.

how can i use a loop over q_int

I really wonder how to solve this problem, I have passed much times to solve it but I always failed.

Matt J
on 31 Jan 2023

Edited: Matt J
on 31 Jan 2023

i put the condition N(i,j)=1, in order to consider just ones in N and based on it S will be built.

As soon as you wrote,

S = optimvar("S",[M,T],"Type","integer","LowerBound",0,"UpperBound",1);

you told the solver that there will be M*T unknown S values that it needs to solve for. The solver knows nothing about how the sparsity of N is supposed to reduce the number of unknowns. To communicate that to the solver, you could try this,

LB=zeros(M,T);

S = optimvar("S",[M,T],"Type","integer","LowerBound",LB,"UpperBound",N);

without changing anything else.

M is unknown , it is generated randomly so its value can be changed but T is fix.

No, M is fixed and known once the optimization starts, which is all that matters.

how can i use a loop over q_int

In the mathematics that you posted, constraint (C6) can be re-arranged in the form

S(i,j)<=UglyExpression(qu)

When q_int is fixed, you can precompute the right hand side. This gives a simple upper bound on S(i,j). You can then launch an optimization problem with

S = optimvar("S",[M,T],"Type","integer","LowerBound",LB,...

"UpperBound",min(N,UglyExpression(qu(q_int))));

The only other constraint on S is the single linear inequaltiy (C7). The whole thing reduces to a pretty tractable looking integer linear program.

Maria
on 1 Feb 2023

Edited: Maria
on 1 Feb 2023

@Torsten i'm sorry, i didn't see your comment,

so if this "1" satisfied the two constraints

??Where do we find this in your list of constraints ?

in constraint (C6) and (C7),i multiplied by S(i,j) in order to take 1 if constraints are verified otherwise "0"

Matt J
on 1 Feb 2023

Edited: Matt J
on 1 Feb 2023

you mean that i replace the constraint by this expression, sorry i didn't get it well.

Let's take a simpler example. Suppose s is the unknown variable in a constraint that looks like,

Q - s + F(qu) >= K*s

where K is positive. Can you see that this constraint is the same as,

s <= (Q+F(qu))/(K+1)

Maria
on 1 Feb 2023

Edited: Maria
on 1 Feb 2023

i tried to apply your suggestion, but i'm asking how can i define just one input in function that contains many other functions with many inputs, for example

function [myNlCON] = Myconstraint1(qu) % here i wrote just (qu) as input

qu = Generate_position(Flight_direction,q_int,UAV_Speed,dt,L,T);

arg1 = argument1(M,N,T,H,qu,qm);

arg2 = argument2(arg1);

Z = sumArgument(arg2);

Transmission_Delay_cons = nonlcons(Y,N,T,alpha,BP,Z);

Delay_exec = DelayAndProcessing(Y,N,D,T);

myNlCON = Transmission_Delay_cons + Delay_exec ;

end

and i called it in my code as, here it demands arguments,should i wrote all inputs or there is a way to just write one specific

constr1 = S >= Q/(Myconstraint1(qu));

Torsten
on 1 Feb 2023

Could you check

class(Myconstraint1(qu))

size(Myconstraint1(qu))

class(Q)

size(Q)

when you run your code ?

Matt J
on 1 Feb 2023

Edited: Matt J
on 1 Feb 2023

@Maria i tried to apply your suggestion, but i'm asking how can i define just one input in function that contains many other functions with many inputs

Nobody said it had to have just one input. However, if you study my ealier optimizeIt.m file, you will see I avoid passing the fixed parameters as function inputs by applying the technique described here. In particular, if Myconstraint1 is nested inside some outer function, it can see (and change!) all the variables in the workspace of that function, simplifying the code to,

function outerFunction()

Flight_direction,UAV_Speed=... %The nested functions below can use these

M,N,T,H,qm,dt,L=....

Delay_exec = DelayAndProcessing(Y,N,D,T);

function [myNlCON] = Myconstraint1(q_int) %nested in outerFunction()

qu = Generate_position(q_int);

arg1 = argument1(qu);

arg2 = argument2(arg1);

Z = sumArgument(arg2);

Transmission_Delay_cons = nonlcons(Z);

myNlCON = Transmission_Delay_cons + Delay_exec ;

end

function qu = Generate_position(q_int) %also nested

...

end

function arg1=argument1(qu) %also nested

...

end

function arg2=argument2(arg1) %also nested

...

end

function Z=sumArgument(arg2) %also nested

...

end

function Transmission_Delay_cons = nonlcons(Z) %also nested

...

end

end

Torsten
on 1 Feb 2023

Maria
on 1 Feb 2023

Edited: Maria
on 1 Feb 2023

@Matt J i put all functions in outerFunction() but they still also in ohter scripts, so an error is appeared "Unrecognized function or variable " i already see your file but you are already put all functions in one script. But for me each function is in a script.

the nested functions should not exist also in other script right? only in outerFunction()?

Matt J
on 1 Feb 2023

Edited: Matt J
on 1 Feb 2023

the nested functions should not exist also in other script right? only in outerFunction()?

By definition, if a function is nested inside outerFunction(), it will be in the same mfile as outerFunction.

Using nested functions was a recommendation only. You can pass all variables around as function inputs as you were doing. However, as you can see, that can get very cumbersome, and the code hard to read.

As you can see also, not all of the functions in optimizeIt.m were nested, only the ones which were specific to computing the objective function and constraints. All of the set-up functions like DelayAndProcessing were not nested, and could have been put in separate files. However, I don't know why you would want to have things like argument1() and argument2() in separate files. They seem very specific to the optimization that you are doing.

Maria
on 1 Feb 2023

Edited: Maria
on 1 Feb 2023

@Matt J i did it now, I returned to the code you provided and i have changed size of M , the problem is solved with large T, i want to verify if the q_int displayed is exactly the first point of my vector qu or not, also i want to show my constraints in order to verify the result.

i wrote show (qu) nothing was displayed.

I also want to know if there is a number of iterations was done before extracting the final result,

Matt J
on 1 Feb 2023

@Maria It sounds like you have run an optimization and now have solution structure, sol with

sol.q_int=...

sol.S=...

Is that correct? If so, you would get qu by doing

qu = Generate_position(sol.q_int)

However, I expect you will not find that qu(1)=sol.q_int, but rather qu(1)=sol.q_int*UAV_Speed*dt

Matt J
on 1 Feb 2023

Edited: Matt J
on 1 Feb 2023

@Torsten Well, it depends. It won't be found if you query qu from the base workspace, but it will be found if you query it inside OuterFunction. Or, OuterFunction can return a handle to Generate_position() to the base workspace, where it could be invoked freely:

function [sol,h]=outerFunction()

Flight_direction,UAV_Speed=... %The nested functions below can use these

M,N,T,H,qm,dt,L=....

Delay_exec = DelayAndProcessing(Y,N,D,T);

h.Generate_position =@Generate_position;

function [myNlCON] = Myconstraint1(q_int) %nested in outerFunction()

qu = Generate_position(q_int);

...

end

function qu = Generate_position(q_int) %also nested

...

end

....

end

Or, you could move Generate_position back to its own separate file and use a nested wrapper

function [sol,h]=outerFunction()

Flight_direction,UAV_Speed=... %The nested functions below can use these

M,N,T,H,qm,dt,L=....

Delay_exec = DelayAndProcessing(Y,N,D,T);

function [myNlCON] = Myconstraint1(q_int) %nested in outerFunction()

qu = get_qu(q_int);

...

end

function qu = get_qu(q_int) %nested wrapper

qu=Generate_position(Flight_direction,q_int,UAV_Speed,dt,L,T);

end

....

end

Maria
on 1 Feb 2023

Edited: Maria
on 1 Feb 2023

@Torsten @Matt J qu = Generate_position(sol.q_int) works well and it displayed the vector qu where q(1) = q_int*UAV_speed*dt.

I know I asked a lot of questions, and thanks to you I managed to solve my problem, but just if you could tell me if I can add a command line to see the number of iterations performed.

Maria
on 1 Feb 2023

@Matt J Does the number of generations refer to the number of iterations?

report =

struct with fields:

problemtype: 'integerconstraints'

rngstate: [1×1 struct]

generations: 68

funccount: 6566

message: 'Optimization terminated: average change in the penalty fitness value less than options.FunctionTolerance↵and constraint violation is less than options.ConstraintTolerance.'

maxconstraint: 0

hybridflag: []

solver: 'ga'

Maria
on 9 Feb 2023

I would like to ask you if I want to change the algorithm to solve my problem to make a comparison, which algorithm could I switch on?

Thank you,

Best regards

Matt J
on 9 Feb 2023

Matt J
on 9 Feb 2023

Matt J
on 9 Feb 2023

Matt J
on 9 Feb 2023

Edited: Matt J
on 9 Feb 2023

Your code, in outline, will look like below. Not all the LPs that you loop over will have a solution. That is why it is important to use the exitflag, as indicated, to check whether a solution was found.

for q_int=0:L

i=q_int+1;

S=optimvar(...)

probLP=....

[sol{i}, fval{i}, exitflag] = solve(probLP)

if exitflag<=0 %solver failed

fval{i}=-inf;

end

end

[best_fval,imax]=max([fval{:}]);

solution=sol{imax};

Torsten
on 9 Feb 2023

q_int goes from 0 to floor(L/(UAV_Speed*dt)) because qu(1) is defined as qu(1) = q_int*UAV_Speed*dt.

Matt J
on 9 Feb 2023

Torsten
on 9 Feb 2023

Edited: Torsten
on 9 Feb 2023

qu is supposed to reach a maximum of L and then goes back directed towards 0 by changing the variable "direction" to "-direction".

q_int as an integer variable has to satisfy 0 <= q_int <= floor(L/(UAV_Speed*dt)) for that qu(1) = q_int*UAV_Speed*dt does not exceed L.

Maria
on 9 Feb 2023

@Matt J @Torsten it works now but i don't know how to hide this, and all statements are finished by(;)

Solving problem using intlinprog.

LP: Optimal objective value is -53.000000.

Optimal solution found.

Intlinprog stopped at the root node because the objective value is within a gap tolerance of the optimal value,

options.AbsoluteGapTolerance = 0 (the default value). The intcon variables are integer within tolerance,

options.IntegerTolerance = 1e-05 (the default value).

ans = % how to hide this?

1×54 cell array

Columns 1 through 7

{0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double}

Columns 8 through 14

{0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double}

Columns 15 through 21

{0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double}

Columns 22 through 28

{0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double}

Columns 29 through 35

{0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double}

Columns 36 through 42

{0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double}

Columns 43 through 49

{0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double} {0×0 double}

Columns 50 through 54

{0×0 double} {0×0 double} {0×0 double} {0×0 double} {1×1 struct}

Matt J
on 9 Feb 2023

Does the result agree with (or improve upon) the result that ga was giving you?

Torsten
on 9 Feb 2023

Edited: Torsten
on 9 Feb 2023

And is it reasonable that for only one value of q_int, a solution exists ? At least that's how I interprete the {0x0 double}.

I would have thought there should always be a solution, but that you have to find the best one after you tried each possible value for q_int: [best_fval,imax]=max([fval{:}]);

Matt J
on 9 Feb 2023

qu is supposed to reach a maximum of L and then goes back directed towards 0 by changing the variable "direction" to "-direction".

Here is the code for qu that I submitted as my answer, and which Maria accepted. If qu is taking steps of UAV_Speed*dt, then a direction change could only occur if L were landed upon exactly. For that to occur, L would need to be an exact integer multiple of UAV_Speed*dt. If L is always an exact multiple of UAV_Speed*dt, then there would be no need for the use of the floor() command in computing the upper limit of q_int.

function [qu] = Generate_position(q_int)

qu = zeros(1,T);

qu(1) = q_int;

direction = Flight_direction;

for k=2:T

qu(k) = qu(k-1) + direction;

if qu(k)==0 || qu(k)==L

direction = -direction;

end

end

qu=qu*UAV_Speed*dt;

end

Torsten
on 9 Feb 2023

Yes, L is an exact multiple of UAV_Speed*dt. I only took the floor function to avoid floating point errors.

I don't know what code exactly you mean that Maria accepted, but I supplied this one:

function [qu] = Generate_position(Flight_direction,q_int,UAV_Speed,dt,L,T)

direction = Flight_direction;

if q_int == 0 % here the error

direction = 1;

end

if q_int == floor(L/(UAV_Speed*dt))

direction = -1;

end

qu = zeros(1,T); %<----Mattt J inserted

qu(1) = q_int*UAV_Speed*dt;

for k=2:T

qu(k) = qu(k-1) + direction*UAV_Speed*dt;

if qu(k)==0 || qu(k)==L

direction = -direction;

end

end

end

Maria
on 9 Feb 2023

Torsten
on 9 Feb 2023

I don't have a problem with qu or q_int now.

Well, it will make a difference if you loop q_int =0:L or q_int = 0:floor(L/(UAV_Speed*dt)) or if you use Matt's or my calculation of qu.

Matt J
on 9 Feb 2023

Maria
on 9 Feb 2023

Maria
on 21 Feb 2023