Help with fminsearch to maximise non-analytic function via specific parameters

I have a simple trajectory model. Let's pretend its non-analytic as it would be with drag. Key variables of interest are launch elevation, theta, and impact range z(:,1).

Is it possible to parameterise the objective function to allow fminsearch to maximise (negative minimise) the function value for range with theta as the independent variable?

A plot of it would look similar to below but with theta (0 to 90 deg) replacing range, and range replacing height ( 0 to max range at theta=45 deg).

I'm struggling to set up theta as the driver and range as the function output to drive such a calculation. Basically finding theta = 45 deg is the maximum range.

Something like [x,fval] = fminsearch(@traj, x0) ???

Any help setting this up would be very much appreciated.

clear,clc

global rho d m g
g = 9.81; rho = 1.2; m = 0.145; d = 0.0732; v = 100;
tspan = 0:0.001:20;

for theta = 35:10:55
   z0 = [0 0 v*cosd(theta) v*sind(theta)];
   options = odeset('Events',@trajevents);
   [~,z,~,ze,ie] = ode45(@traj,tspan,z0,options);
   hold on, plot(z(:,1),z(:,2))
   plot(ze(ie(1),1),ze(ie(1),2),'k+')
   xlabel('range'), ylabel('height')
   axis equal, xlim([0 1100]); ylim([0 400]); grid
end

function zdot = traj(~,z)
   global g
   zdot = [z(3); z(4); 0; -g];
end

function [check,stop,direction] = trajevents(~,z)
   check = z(2); stop = 1; direction = -1;
end

Answer 1

[Torsten]

Here is the optimization solution:

theta0 = 10;
sol = fminsearch(@fun,theta0)

sol = 45

function range = fun(theta)
  global rho d m g
  g = 9.81; rho = 1.2; m = 0.145; d = 0.0732; v = 100;
  tspan = 0:0.001:20;
  z0 = [0 0 v*cosd(theta) v*sind(theta)];
  options = odeset('Events',@trajevents);
  [~,z,~,ze,~] = ode45(@traj,tspan,z0,options);
  range = -z(end,1) ;  
end
function zdot = traj(~,z)
  global g
  zdot = [z(3); z(4); 0; -g];
end
function [check,stop,direction] = trajevents(~,z)
  check = z(2); stop = 1; direction = -1;
end

[Brantosaurus]

Thank you very much for your advice on this.

I had struggled with this for some time.

Much appreciated :)

Answer 2

[Torsten]

Of course, you could use fminsearch for optimization. But using it, you will only get a single point of the below graph.

Theta = 1:89;
global rho d m g
g = 9.81; rho = 1.2; m = 0.145; d = 0.0732; v = 100;
tspan = 0:0.001:20;
range = zeros(size(Theta));
for i=1:numel(Theta)
   theta = Theta(i);
   z0 = [0 0 v*cosd(theta) v*sind(theta)];
   options = odeset('Events',@trajevents);
   [~,z,~,ze,~] = ode45(@traj,tspan,z0,options);
   range(i) = z(end,1) ;  
end
plot(Theta,range)
grid on

function zdot = traj(~,z)
   global g
   zdot = [z(3); z(4); 0; -g];
end
function [check,stop,direction] = trajevents(~,z)
   check = z(2); stop = 1; direction = -1;
end