MATLAB: How to use BVP5c (or 4c) iteratively to solve BVP for different values of a known parameter

bvp4cbvp5c

I have managed to solve a system of 4 (first-order autonomous) ODE's subject to 7 boundary conditions and for fixed values of my known parameter g. However, now I want to vary g and plot the different curves. Originally I specified the value of g inside my ODE and ODE_bc functions (see below) – e.g. I may write g=1; just below where I specified my unknown parameter r. However, now it seems that I need to specify g outside my functions in order to iterate over different values of g (using a for loop) and using the code below I get the error:
Error using bvp5c
Too many input arguments.
Error in Iterated_ODE (line 4)
sol1 = bvp5c(@simp_ODE,@simp_bc,solinit1,options,g);
I've tried taking the 4th argument g out entirely for the simp_ODE and simp_bc and then for the 5th argument in BVP5c (above) but with no success.
How can I correct this problem?
(see current code below)
options = bvpset('stats','on');
g = 1;
solinit1 = bvpinit(linspace(0,1,20),@ODE_init,[0.75 -0.7 1.47]); %initial guess
sol1 = bvp5c(@simp_ODE,@simp_bc,solinit1,options,g);
sol = bvp5c(@ODE,@ODE_bc,sol1,options,g);
plot(sol.y(1,:),sol.y(3,:)); %plots y(3)=theta vs y(1)=x $
axis equal
hold on
for i=1:10
g = g+0.5; %Increment known/prescribed parameter
sol = bvp5c(@ODE,@ODE_bc,sol,options); %use previous solution as inital guess
plot(sol.y(1,:),sol.y(3,:)); %plots y(3)=theta vs y(1)=x $
axis equal
hold on
end
function G = ODE_init(s)
G = [s %guess x(s) = s%
0
0
0];
end
function dydx = simp_ODE(x,y,array, g) %first-order simplified ODE using cos(x) = 1+O(x), sin(x)=x+O(x^3)%
b = array(1);
n=array(2);
r = array(3);
dydx = [x
y(2)
y(3)
n*y(2)-g*(x+y(1)*y(2))/r]; % b,n,r are unknown parameters to solve for%
end
%tried to specify unknown parameters using array and g is a known parameter%
function res = simp_bc(ya,yb,array,g) %bc's to simplified ODE using a 'first order approximation%
b = array(1) ;
n = array(2);
r = array(3);
res = [ya(1) %x(0)=0%
ya(2) %y(0) = 0%
ya(3) %theta(0) = 0%
yb(4) %psi(1) = theta'(1) =0 in non-dim bc%
yb(3)+b-pi/3 %Wetting angle relation at boundary theta(1)+beta = pi/6%
1 – r*sin(b) %geometric condition%
yb(2) – r*(n +cos(b))]; %internal force y(1) -r*n-r*cos(b) == 0%
end
function dyds = ODE(s,y,array,g)
b = array(1);
n = array(2);
r = array(3);
% g = 1;
dyds = [cos(y(3));
sin(y(3));
y(4);
n*sin(y(3))-g*(y(1)+y(2)*sin(y(3)))/r];
end
function res = ODE_bc(ya,yb,array,g)
b = array(1) ;
n = array(2);
r = array(3);
res = [ya(1) %x(0)=0%
ya(2) %y(0) = 0%
ya(3) %theta(0) = 0%
yb(4) %psi(1) = theta'(1) =0 in non-dim bc%
yb(3)+b-pi/3 %Wetting angle relation at boundary theta(1)+beta = pi/6%
yb(1) – r*sin(b) %geometric condition%
yb(2) – r*(n +cos(b))]; %internal force bc y(1) -r*n-r*cos(b) == 0%
end

Best Answer

sol = bvp5c(@(x,y)ODE(x,y,ar,g),@(a,b)ODE_bc(a,b,ar,g),sol,options); %use previous solution as inital guess