The main part of the lecture was devoted to BVP whose boundary conditions are given by derivatives. Our running example is example 7.6 of the text book (pages 538-539):
2
d u(x)
-------- = u(x) with u'(1) = 1.17520 and u'(3) = 10.0179.
2
d x
For this problem we use MATLAB to apply the shooting method and
the finite differences method.
To use the shooting method in MATLAB, we need to translate the second-order differential equation into a system of two first-order equations. The right-hand side of this system is evaluated with the following .m file (save as ex7_6.m):
function ode_rhs = ex7_6(x,u) % % ode_rhs = ex7_6(x,u) % % The second-order differential equation u" = u is % written as a system of two first-order equations: % u'(1) = u(2) % u'(2) = u(1) % This function returns the right-hand side of the system. % ode_rhs = [u(2); u(1)];Since we need to shoot at least twice, we use a function to solve the IVP:
function endpt = shoot_ex7_6(startpt) % % endpt = shoot_ex7_6(startpt) % % Solves the second-order differential equation u" = u, % with u(1) = startpt and u'(1) = 1.17520. % On return is the value for u'(3). % span = [1 3]; uone = [startpt; 1.17520]; [x,u] = ode45(@ex7_6,span,uone); plot(x,u(:,1)); endpt = u(size(u,1),2);Then here is the output of the MATLAB session to find the appropriate value for u(1) with the shooting method:
>> end1 = shoot_ex7_6(1) end1 = 8.04819372641815 >> end2 = shoot_ex7_6(2) end2 = 11.67505449896711 >> shoot_ex7_6(1.54308) ans = 10.01786927882996We found the value 1.54308 by linear interpolation. Consult your class notes or the text book for further explanation.
Below is a MATLAB function (save as findif_ex7_6.m) which generates and solves the linear system for this BVP:
function [x,u] = findif_ex7_6(n)
%
% [x,u] = findif_ex7_6(n)
% applies finite differences with n subintervals
% to elaborate Example 7.6 on pages 538-539
%
% While executing, the matrix M of the system, and right-hand side
% vector R are displayed on screen.
%
a = 1; % left bound
uda = 1.17520; % value for derivative of u at 1
b = 3; % right bound
udb = 10.0179; % value for derivative of u at 3
h = (b-a)/n; % step size
for i = 1:n+1 % the nodes
x(i) = a + (i-1)*h;
end;
M = zeros(n+1,n+1); % matrix of the linear system
for i = 1:n+1
M(i,i) = -(2+h^2);
if (i < n+1)
M(i,i+1) = 1;
end;
if (i > 1)
M(i,i-1) = 1;
end;
end;
M(1,2) = 2;
M(n+1,n) = 2;
R = zeros(n+1,1); % right-hand side vector
R(1) = 2*h*uda;
R(n+1) = -2*h*udb;
M
R
u = M\R;
x = x';
The output of one session to reproduce the numbers in the book is
>> [x,u] = findif_ex7_6(4)
M =
-2.2500 2.0000 0 0 0
1.0000 -2.2500 1.0000 0 0
0 1.0000 -2.2500 1.0000 0
0 0 1.0000 -2.2500 1.0000
0 0 0 2.0000 -2.2500
R =
1.1752
0
0
0
-10.0179
x =
1.0000
1.5000
2.0000
2.5000
3.0000
u =
1.5522
2.3338
3.6989
5.9887
9.7757
Note that to get the same accuracy as with the shooting method,
we need about 200 intervals.