optimization of a function with array of parameters using matlab solver ga

Hello, i i am trying to optimize function S, so it will get values between 0 and 10 using its E_opt parameter

coefficients which shown bellow in the expression of S.

I defined optimvar for each sub cell of the E_opt, and used solve as shown bellow.

but then i get an error on optimvar saying E_opt(1) is not a valid matlab name.

how do we run such optimization so i will get optimized array of parameters E_opt?

Thanks.

clc
clear all
f=6;
lambda=300/6;
theta_in_rad=(0/180)*pi;
theta_out_rad=(70/180)*pi;
x=linspace(-15*lambda,15*lambda,1000);
y=linspace(-15*lambda,15*lambda,1000);
[xx,yy]=meshgrid(x,y);
k=2*pi/lambda;
E_i=1;
E_r=E_i*sqrt(cos(theta_in_rad)/cos(theta_out_rad))
eta=377;
D_x=1.064*lambda;
k_x=k*sin(theta_in_rad)+2*pi/D_x;
k_y=sqrt(k^2-k_x.^2);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
n=[2,3,4,5,6,7,8];
k_x_n=k*sin(theta_in_rad)+(2*pi.*n)/D_x;
k_y_n=sqrt(k^2-k_x_n.^2);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%



E_opt(1)=optimvar('E_opt(1)','Lower',0,'Upper',+1);
E_opt(2)=optimvar('E_opt(2)','Lower',0,'Upper',+1);
E_opt(3)=optimvar('E_opt(3)','Lower',0,'Upper',+1);
E_opt(4)=optimvar('E_opt(4)','Lower',0,'Upper',+1);
E_opt(5)=optimvar('E_opt(5)','Lower',0,'Upper',+1);
E_opt(6)=optimvar('E_opt(6)','Lower',0,'Upper',+1);
E_opt(7)=optimvar('E_opt(7)','Lower',0,'Upper',+1);


S=@(E_opt)abs(E_i.*exp(-j*k*xx*sin(theta_in_rad)+j*k*yy*cos(theta_in_rad))+E_r.*exp(-j*k_x*xx-j*k_y*yy)+E_opt(1).*exp(-j*k_x_n(1)*xx-j*k_y_n(1)*yy)+E_opt(2).*exp(-j*k_x_n(2)*xx-j*k_y_n(2)*yy)+E_opt(3).*exp(-j*k_x_n(3)*xx-j*k_y_n(3)*yy)+E_opt(4).*exp(-j*k_x_n(4)*xx-j*k_y_n(4)*yy)+E_opt(5).*exp(-j*k_x_n(5)*xx-j*k_y_n(5)*yy)+E_opt(6).*exp(-j*k_x_n(6)*xx-j*k_y_n(6)*yy)+E_opt(7).*exp(-j*k_x_n(7)*xx-j*k_y_n(7)*yy).*conj((j*k*cos(theta_in_rad))*E_i.*exp(-j*k*xx*sin(theta_in_rad)+j*k*yy*cos(theta_in_rad))+(-j*k_y)*E_r.*exp(-j*k_x*xx-j*k_y*yy)+(-j*k_y_n(1))*E_opt(1).*exp(-j*k_x_n(1)*xx-j*k_y_n(1)*yy)+(-j*k_y_n(2))*E_opt(2).*exp(-j*k_x_n(1)*xx-j*k_y_n(2)*yy)+(-j*k_y_n(3))*E_opt(3).*exp(-j*k_x_n(3)*xx-j*k_y_n(3)*yy)+(-j*k_y_n(4))*E_opt(4).*exp(-j*k_x_n(4)*xx-j*k_y_n(4)*yy)+(-j*k_y_n(5))*E_opt(5).*exp(-j*k_x_n(5)*xx-j*k_y_n(5)*yy)+(-j*k_y_n(6))*E_opt(6).*exp(-j*k_x_n(6)*xx-j*k_y_n(6)*yy)+(-j*k_y_n(7))*E_opt(7).*exp(-j*k_x_n(7)*xx-j*k_y_n(7)*yy)));

S_max=@(E_opt)max(abs(E_i.*exp(-j*k*xx*sin(theta_in_rad)+j*k*yy*cos(theta_in_rad))+E_r.*exp(-j*k_x*xx-j*k_y*yy)+E_opt(1).*exp(-j*k_x_n(1)*xx-j*k_y_n(1)*yy)+E_opt(2).*exp(-j*k_x_n(2)*xx-j*k_y_n(2)*yy)+E_opt(3).*exp(-j*k_x_n(3)*xx-j*k_y_n(3)*yy)+E_opt(4).*exp(-j*k_x_n(4)*xx-j*k_y_n(4)*yy)+E_opt(5).*exp(-j*k_x_n(5)*xx-j*k_y_n(5)*yy)+E_opt(6).*exp(-j*k_x_n(6)*xx-j*k_y_n(6)*yy)+E_opt(7).*exp(-j*k_x_n(7)*xx-j*k_y_n(7)*yy).*conj((j*k*cos(theta_in_rad))*E_i.*exp(-j*k*xx*sin(theta_in_rad)+j*k*yy*cos(theta_in_rad))+(-j*k_y)*E_r.*exp(-j*k_x*xx-j*k_y*yy)+(-j*k_y_n(1))*E_opt(1).*exp(-j*k_x_n(1)*xx-j*k_y_n(1)*yy)+(-j*k_y_n(2))*E_opt(2).*exp(-j*k_x_n(1)*xx-j*k_y_n(2)*yy)+(-j*k_y_n(3))*E_opt(3).*exp(-j*k_x_n(3)*xx-j*k_y_n(3)*yy)+(-j*k_y_n(4))*E_opt(4).*exp(-j*k_x_n(4)*xx-j*k_y_n(4)*yy)+(-j*k_y_n(5))*E_opt(5).*exp(-j*k_x_n(5)*xx-j*k_y_n(5)*yy)+(-j*k_y_n(6))*E_opt(6).*exp(-j*k_x_n(6)*xx-j*k_y_n(6)*yy)+(-j*k_y_n(7))*E_opt(7).*exp(-j*k_x_n(7)*xx-j*k_y_n(7)*yy))));
S_min=@(E_opt)min(abs(E_i.*exp(-j*k*xx*sin(theta_in_rad)+j*k*yy*cos(theta_in_rad))+E_r.*exp(-j*k_x*xx-j*k_y*yy)+E_opt(1).*exp(-j*k_x_n(1)*xx-j*k_y_n(1)*yy)+E_opt(2).*exp(-j*k_x_n(2)*xx-j*k_y_n(2)*yy)+E_opt(3).*exp(-j*k_x_n(3)*xx-j*k_y_n(3)*yy)+E_opt(4).*exp(-j*k_x_n(4)*xx-j*k_y_n(4)*yy)+E_opt(5).*exp(-j*k_x_n(5)*xx-j*k_y_n(5)*yy)+E_opt(6).*exp(-j*k_x_n(6)*xx-j*k_y_n(6)*yy)+E_opt(7).*exp(-j*k_x_n(7)*xx-j*k_y_n(7)*yy).*conj((j*k*cos(theta_in_rad))*E_i.*exp(-j*k*xx*sin(theta_in_rad)+j*k*yy*cos(theta_in_rad))+(-j*k_y)*E_r.*exp(-j*k_x*xx-j*k_y*yy)+(-j*k_y_n(1))*E_opt(1).*exp(-j*k_x_n(1)*xx-j*k_y_n(1)*yy)+(-j*k_y_n(2))*E_opt(2).*exp(-j*k_x_n(1)*xx-j*k_y_n(2)*yy)+(-j*k_y_n(3))*E_opt(3).*exp(-j*k_x_n(3)*xx-j*k_y_n(3)*yy)+(-j*k_y_n(4))*E_opt(4).*exp(-j*k_x_n(4)*xx-j*k_y_n(4)*yy)+(-j*k_y_n(5))*E_opt(5).*exp(-j*k_x_n(5)*xx-j*k_y_n(5)*yy)+(-j*k_y_n(6))*E_opt(6).*exp(-j*k_x_n(6)*xx-j*k_y_n(6)*yy)+(-j*k_y_n(7))*E_opt(7).*exp(-j*k_x_n(7)*xx-j*k_y_n(7)*yy))));

cam=fcn2optimexpr( S,E_opt);
camMax=fcn2optimexpr( S_max,E_opt);
camMin=fcn2optimexpr( S_min,E_opt);
prob = optimproblem("Objective",camMax);
prob.Constraints.camUpper=camMax<=10;
prob.Constraints.camLower=camMin>=0;
[sol,fval] = solve(prob,"Solver","ga");

[fima v]

Hello ,my E_opt array of variables needed to be complex numbers

is there a way to change the command bellow so each cell will be complex so the amplitude will be between 1 to 7 and the phase is between 0 to 220 degrees?

Thanks.

E_opt = optimvar('E_opt',7,'Lower',0,'Upper',1);

Answer 1

[Torsten]

E_opt = optimvar('E_opt',7,'Lower',0,'Upper',+1);

instead of

E_opt(1)=optimvar('E_opt(1)','Lower',0,'Upper',+1);
E_opt(2)=optimvar('E_opt(2)','Lower',0,'Upper',+1);
E_opt(3)=optimvar('E_opt(3)','Lower',0,'Upper',+1);
E_opt(4)=optimvar('E_opt(4)','Lower',0,'Upper',+1);
E_opt(5)=optimvar('E_opt(5)','Lower',0,'Upper',+1);
E_opt(6)=optimvar('E_opt(6)','Lower',0,'Upper',+1);
E_opt(7)=optimvar('E_opt(7)','Lower',0,'Upper',+1);

and

[sol,fval] = solve(prob);

instead of

[sol,fval] = solve(prob,"Solver","ga");

[fima v]

Hello Torsten, three question regarding your great answer.

1. your optim var fix puts the same boundary on every cell of the array ,what if i want the second cell be betweem 0 and 5 ,how can it been done?
2. In the solve comand you proposed there is no mentioning that its "GA" algorithm.how the solve command chooses to use GA?
3. in my solution i want my S function to ne optimised so it will be from 0 to 10. why do we use the S_max in as shown bellow.

prob = optimproblem("Objective",camMax);

Thanks.

[Torsten]

your optim var fix puts the same boundary on every cell of the array ,what if i want the second cell be betweem 0 and 5 ,how can it been done?

Specify "Lower" and "Upper" as a vector of length 7, not as a scalar.

In the solve comand you proposed there is no mentioning that its "GA" algorithm.how the solve command chooses to use GA?

Running your modified code with "Solver","ga" gave an error. You might try without the specification of a solver and see which one MATLAB chooses.

in my solution i want my S function to ne optimised so it will be from 0 to 10. why do we use the S_max in as shown bellow.

prob = optimproblem("Objective",camMax);

I don't understand what you mean. You chose the objective to be camMax.

[Walter Roberson]

solve() would choose to use ga if it detected any integer constraints. I am not sure if there are other situations it will choose to use ga

The max() in your functions has a discontinuous derivative, and so mathematically cannot be used with fmincon() or some of the other minimizers. I do not know, however, whether solve() will detect the max() and know to avoid using those solvers.

[fima v]

Hello walter ,

1/how did you know that Smax has discontinues derivative?

2.how would you formulate the last part shown bellow.

i just want to get S be 0 and 10.

i dont know why we have to use extra smin smax function and put its in constraint

why not put S (cam as constraint)?

Thanks.

[code]

prob = optimproblem("Objective",camMax);

prob.Constraints.camUpper=camMax<=10;

prob.Constraints.camLower=camMin>=0;

[sol,fval] = solve(prob,"Solver","ga");

[/code]

[Walter Roberson]

The derivative of max() is always discontinuous except in the case where it can be proven that the max() is always selecting the same entry.

Consider for example,

f(x) = max(5, x)

then df/dx is 1 when x > 5 where the max() is always going to be selecting x. When x < 5, then max(x,5) is always 5, and the derivative of that is 0. So if we look at the plot of the derivative we see constant 0 until x = 5, then we see constant 1. What is the derivative at x = 5 exactly ? In any closed compact region around x = 5, the derivative is 0 on one side and 1 on the other side: df at (5-eps) is 0 and df(5+eps) is 1 for all non-zero eps including infinitestimal eps. Any finite value you assign to df at x = 5, it there will be a jump in the derivative as you move to at least one side of x = 5. Therefore the derivative is discontinuous.

Discontinuous derivative is not a problem for ga() but is a problem for fsolve() and fminsearch() and fmincon()

[fima v]

cam function gave me an error

why its logical to take max of this function?

prob = optimproblem("Objective",camMax);

[Torsten]

Did you copy the code somewhere or why do you ask questions that you had to answer ?

[Walter Roberson]

You wish to "optimize" a function but have constaints on what the value of the function is.

If you are minimizing a function S and wish to constrain the minimized value of the function is at least 0, then there are two possibilities:

• The function is naturally non-negative and the >= 0 constraint has no effect. In this case, minimizing the function without the constraint will return a positive minima for the function; OR
• The function crosses 0, and the >= 0 constraint does have an effect. In this case you wish to know the location where it reaches 0. But you can get that by solving the function for 0.

This leads to the following strategy if you are doing minimization

1. minimize the function without the constraints
2. if the value at the minima is greater than 10 then there is no solution that is consistent with your constraints. Solving for the function equalling 10 will not work because if the minima were 10 or lower it would have been found by the minimizer
3. If the value of the unconstrained minima is between 0 and 10 (inclusive) then that fits within your constraints, so accept the location without further processing
4. If the value of the unconstrained minima is less than 0, then attempt to solve for the equations equalling 0 exactly. If the equations always return negative than this will fail. If the equations cross 0 then solving for 0 will find such a location, and that will satisfy the constraints that 0 <= S <= 10 because it will be a point at which S == 0 by design.

You can use a very similar strategy if you are maximizing instead of minimizing, except that in that case you would potentially be solving for S == 10 rather than S == 0

[fima v]

Hello ,my E_opt array of variables needed to be complex numbers

is there a way to change the command bellow so each cell will be complex so the amplitude will be between 1 to 7 and the phase is between 0 to 220 degrees?

Thanks.

E_opt = optimvar('E_opt',7,'Lower',0,'Upper',1);

[Torsten]

Read the warning in the MATLAB documentation:

https://de.mathworks.com/help/optim/ug/complex-numbers-in-optimization-toolbox-solvers.html

Warning

The problem-based approach does not support complex values in an objective function, nonlinear equalities, or nonlinear inequalities. If a function calculation has a complex value, even as an intermediate value, the final result might be incorrect.