Solving shooting method with ode45

7 views (last 30 days)
beginner
beginner on 21 Mar 2016
Answered: Torsten on 22 Mar 2016
Hi there, I want to solve 2nd order ODE system with the Shooting Method. With my code I can't find any solution...
I have H(inf) = dH/dX(inf) = 0 and dH/dX(0) = 0. What I want to find out is H(0). I used fsolve and ode45.
Is my code correct at all? It doesnt give any error message but it doesn't find any solution...
function dhdt = wihelmy(X,h,Bo)
%h(1) = H
%h(2) = H'
% X is the input.
% tspan 1:3
dhdt_1 = h(2); % dhdt_1 = H'
dhdt_2 = Bo*(h(1)*(1+h(2)))^(3/2); %dhdt_2 = H''
dhdt = [dhdt_1; dhdt_2];
function HdH = test1(h0, theta, Bo)
%Bo = 1; % Bo = rho*g*l^2/gamma
%h0 = ?
% h_inf = 0;
% dh0 = 0;
% dh_inf = 0;
%h0 = 10e-5;
%thetadeg = 90; %theta in degrees
%theta = thetadeg*pi/180;
init = [10e-5 -cot(theta)];
function_handle = @(X,h) wihelmy(X,h,Bo);
[Xout,Hout] = ode45(function_handle, [0:0.1:3], init);
H_x_inf = Hout(end,1);
dH_x_inf = Hout(end,2);
HdH = [H_x_inf; dH_x_inf];
Main function is:
Bo=1;
thetadeg = 85; %theta in degrees
theta = thetadeg*pi/180;
Hopt=fsolve(@(h0) test1(h0, theta, Bo), 10e-5);
Thank you sooo much, if you can help me!

Answers (1)

Torsten
Torsten on 22 Mar 2016
There are so many errors in your code that you should start simple.
Your problem is not a standard boundary value problem.
If your equation reads
y''=Bo*(y*(1+y'))^(3/2)
and you know y'(0)=0, I'd first examine some solutions for different values of y(0) using ODE45 alone.
Then, after you've got a feeling on how your function behaves, I'd switch to bvp4c to impose boundary conditions at x=0 and at x=L for a large value of L.
After you've been successful with these two foregoing steps, you can dare to apply your own solution method, namely shooting.
Best wishes
Torsten.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!