Differential equation concerning echolocation I am in the process of trying to solve a differential equation using MATLAB. 
Given information is the following: 
$$
c(z) = 4800 - 20.2090 + \frac{17.3368z}{1000}+272.9057e^\frac{-0.7528z}{1000}
$$
$$
\frac{d^2z}{dx^2} = -q_0\frac{c'(z)}{c(z)^3}
$$
$$
q_0=(\frac{c(z_0)}{cos(B_0)})^2
$$
$$
x = 0, z = z_0, \frac{dz}{dx} = tan(B_0) 
$$
B0 is the angle between the line z = 2000 and the ray in the start point.
I have solved the equation using MATLAB's ode45 command, but now I'm supposed to use the secant method to find which initial values of $B_0$ also give me an end value of $z = 2500$ at $x = 151900$.
I've been told that there's 4, since I'm supposed to use the $\Delta t$ between the arrivals to pinpoint the initial position, but no matter how I try I cannot figure out how to use the secant method in this particular instance. 
Currently my program looks something like this: 
angle = (-10:0.0035:14); 

z_0 = 2000;



error_tolerance = 1;    % Väljer hur långt ljudvågen
it = 1;

figure
for i = 1:length(angle)

B_0 = tand(angle(i)); 
q_0 = (c_sound(z_0)/cosd(angle(i)))^2;    
[X,Z] = ode45(@soundwave_function, [0 distance], [z_0 B_0]);
if abs(Z(end,1)-2500) <= error_tolerance

    saved_angles(it,1) = angle(i); 

    it = it+1;                           

end

end

angles_str = num2str(saved_angles);

global saved_angles

if numel(saved_angles)>1

for i = 1:numel(saved_angles)
    disp([num2str(saved_angles(i)), ' degrees'])
end

else
    disp([num2str(saved_angles), ' grader'])
end

This does find the angles where $ z(151900) \approx 2500 $ in the interval specified, -10 to 14 degrees (with 0 degrees counting as z = 2000).
I don't know how to approach the problem using the secant method, and since I'm supposed to use that solution for the second bit I don't really know how to continue. 
Is there anyone out there that knows how to apply the secant method in this specific case? Would be much appreciated! 
 A: You need to use a shooting method. Let's abstract it a bit.
Let's define $F(x_{ic})$ to be a function that takes your initial conditions as its argument and outputs the value of your solution at $x = 151900$.
In other words, $F$ is a solution trajectory of your ODE, and we might think of $F$ as a single run of ode45, executed with initial conditions $x_{ic}$.
If we want to compute $F(x_{ic}) - 2500 = 0$, we can find initial values that satisfy the problem. Therefore, we simply want to write a root finding algorithm.
Newton's method is a good fit, but it is hard to compute $\frac{dF}{dx_{ic}}$. Instead, we might use the secant method, also known as Broyden's method.
It is not so hard: simply implement a function for Broyden's method, and wrap your ODE solver in a function that takes initial conditions as its arguments. For instance:
function z151900 = myODE(IC)
      [x,z] = ode45(myode,[0 151900],IC);
      z151900 = z(end);
end;

function xic = SecantMethod(F,IC1,IC2)
      % your secant method iteration using two initial guesses
end;

IC1 = ic1;
IC2 = ic2;
SecantMethod(z151900,ic1,ic2);

I can post a more concrete example later, or you can search my answer history for "Broyden's method".
