clear
echo on
format compact
%Solution to Problem 30
%
%In this problem we are asked to determine the
%displacement thickness for laminar boundary
%layer flow past a flat plate.  The velocity
%profile is described by the Blasius equation,
%and the dimensionless displacement thickness is
%given by the integral of (1-f') where f' is
%the dimensionless velocity.  The integration should
%go to infinity, but instead we will truncate it at
%some finite value etalimit, and increase this until
%the answer for delta does not change.
pause
%We take y(1) to be f, y(2) to be f', y(3) to be f''
%whose value at the wall is unknown, and y(4) to be
%the integral of (1-f').  Thus:
global etaout yout etalimit
etalimit=1;
deltalast=0;
delta=1;
fpp=1;
figure(1)
while abs(delta-deltalast)>1e-5,
 deltalast=delta;
 etalimit=etalimit*2;
 fpp=fzero('bcfit',fpp)
 plot(etaout,yout)
 axis([0,5,0,2]);
 xlabel('eta')
 legend('f','fp','fpp','delta')
 drawnow
 n=length(etaout);
 delta=yout(n,4);
 echo off
end
echo on
%Which does it all!  We can plot up the final result:
figure(2)
plot(etaout,yout(:,2))
axis([0,5,0,1]);
hold on
plot([delta,delta],[0,1],'k')
hold off
xlabel('eta')
ylabel('velocity')
text(delta,.4,'<-- displacement thickness')
%The displacement thickness is thus:
delta
%The problem converged at a distance from the plate of:
etalimit
echo off

 


___________________________




function diff=bcfit(fpp)
%This function integrates the Blasius equation
%subject to the initial conditions f(0)=f'(0)=0,
%and the input guess of f''(0).  It returns the
%difference between f'(0) and unity, which should
%be zero.  The routine is called by the routine
%fzero in the main script, and calls in turn the
%system of differential equations given by blasius.m

global etaout yout etalimit

y0=zeros(4,1);
y0(3)=fpp;
[etaout,yout]=ode45('blasius',[0,etalimit],y0);
n=length(etaout);
diff=yout(n,2)-1;



_____________________________



function ydot=blasius(eta,y)
%This function returns the system of equations
%corresponding to the third order differential
%blasius equation.  We also track the integral
%of (1-f') as a fourth variable, as this is
%used to calculate the displacement thickness.

ydot=zeros(size(y));
ydot(1)=y(2);
ydot(2)=y(3);
ydot(3)=-y(1)*y(3);
ydot(4)=1-y(2);