How to find maximum of a function of 2 variables

I have a function z = cos(x^2 + y^2) and the assumption that both x and y belong to interval 1;5. How can I find the maximum value of the function z? I thought about using fminsearch for finding minimum of -z, but I’m new to Matlab and it doesn’t work.

[Alex Sha]

Obviously, the maximum value of cos() should be 1, in the range of [1,5], there are multi-solutions:

1:

x: 2.32564131982953
y: 4.44118606739161

2:

x: 4.43419681301284
y: 2.33893989844493

Answer 1

[Star Strider]

One approach:

z = @(x,y) cos(x.^2 + y.^2);

XY = fminsearch(@(b)norm(-z(b(1),b(2))), [3; 3])

producing:

XY =
          2.96385862016178
          2.91449853520917

[enter]

How did you get [3;3]?

[Star Strider]

It is necessary to choose inittial estimates for the parameters to optimise, since in most optimisation routines the result is sensitive (sometimes very sensitive) to the initial parameter estimates. I chose [3;3] because it is within the desired [1 5] values you want, and since fminsearch does not allow parameter constraints (as fmincon does, for example).

[enter]

Is it correct way to calculate max within a given interval? What if I had [1;10] interval?

[Star Strider]

Choose any initial parameter estimates for both parameters within your chosen interval. They do not have to be the same values or integers, so for example [5.5; 7.2] returns:

XY =
          5.62165225908844
          7.40208239906857

when I run it.

Feel free to experiment with the code!

[enter]

How do initial parameters influence the result?

[Star Strider]

The optimisation functions search for the closest result to the initial parameter estimates. That is true for fminsearch (that is derivative-free) and all gradient-descent algorithms. It is not true for global search functions such as genetic algorithms (the ga function).

[enter]

So there is a possibility that the result isn’t maximum?

[Star Strider]

There is the possibility that it is not the global maximum, yes. To find the global maximum, it would likely be necessary to use ga or one of the other global optimisers.

An example of usiing the ga function to do the optimisation:

z = @(x,y) cos(x.^2 + y.^2);
ftns = @(b)norm(-z(b(1),b(2)));
PopSz = 500;
Parms = 2;
opts = optimoptions('ga', 'PopulationSize',PopSz, 'InitialPopulationMatrix',randi(1E+4,PopSz,Parms)*1E-3, 'MaxGenerations',2E3, 'PlotFcn',@gaplotbestf, 'PlotInterval',1);
t0 = clock;
fprintf('\nStart Time: %4d-%02d-%02d %02d:%02d:%07.4f\n', t0)
[theta,fval,exitflag,output] = ga(ftns, Parms, [],[],[],[],[1 1],[10 10],[],[],opts)
t1 = clock;
fprintf('\nStop Time: %4d-%02d-%02d %02d:%02d:%07.4f\n', t1)

It performs much better, however it still comes up with sometimes inconsistent results (not always the same parameter values, although quite similar fitness values). Note that it allows constraints on the parameters, that here I limited to 1 and 10 for both.

[enter]

That's quite a sophisticated solution. is there a method like fminbnd, but for multivariable functions?

[Star Strider]

I referred to it earlier: fmincon. The would be more appropriate, since it is for multivariable functions, and it allows constraints on them.