MATLAB: I want a code for solving a coupled 3rd order and 2nd order ode using shooting method and RK-4 numerical technique , please if anyone could help

mathematicsMATLABnumerical analysisrk-4shooting method

(1+2M*eta)f''' + 2M*f"+ f*f" -f'^2 – k1*f' + lambda*theta=0 ———- (1)
(1+2M*eta)theta" + 2M*theta' + Pr(f*theta'-f'*theta)=0 ——-(2)
'f' and 'theta' are functions of 'eta', eta is an independent variable
3 initial conditions are given: eta=0, f(0)=0, f'(0)=1,theta(0)=1
Say I reduce these equations (1) and (2) to five ode (shooting method)
f'=z ; f(0)=0 —–(3)
z'=p ; z(0)=1 ——(4)
p'= (-2M*f"-f*f"+f'^2+k1*f'-lamda*theta)/(1+2*M*eta) ;p(0)= (guess value) —–(5)
theta'= q ; theta(0)=1 —–(6)
q' = (-2M*theta'-Pr(f*theta'-f'*theta))/(1+2*M*eta) ; q(0)= (guess value) ——(7)
The boundary conditions that needs to be satisfied are: f'(eta=10)= 0 and theta(eta=10)=0 as eta=10
Given:
M= 1
k1= 0.1
lamda= 0.1
Pr= 0.7
taking step length: h= 0.01

Best Answer

Try
function main
global Pr k1 M lambda
Pr=0.7; k1=0.1; M=1; lambda=0.1;
rlow=0;
rhigh=10;
N=1000;
options=bvpset('stats','on','RelTol',1e-5);
solinit=bvpinit(linspace(rlow,rhigh,N),[0, -1, 0, 1 ,0]); %
sol=bvp4c(@projode,@mybcs,solinit,options);
function dy = projode(n,y)
global Pr k1 M lambda
dy = [y(2); y(3); (-2*M*y(3)-y(1)*y(3)+y(2)^2+k1*y(2)-lambda*y(4))/(1+2*M*n); y(5); (-2*M*y(5)-Pr*y(1)*y(5)+Pr*y(2)*y(4))/(1+2*M*n)];
function res = mybcs(ya,yb)
res = [ya(1); ya(2)-1.0; ya(4)-1.0 ; yb(2); yb(4)];
Best wishes
Torsten.