Solve boundary value problems for ordinary differential equations
sol = bvp4c(odefun,bcfun,solinit)
sol = bvp4c(odefun,bcfun,solinit,options)
solinit = bvpinit(x, yinit, params)
A function handle that evaluates the differential equations f(x,y). It can have the form
dydx = odefun(x,y) dydx = odefun(x,y,parameters)
For a scalar
A function handle that computes
the residual in the boundary conditions. For two-point boundary value
conditions of the form bc(y(a),y(b)),
res = bcfun(ya,yb) res = bcfun(ya,yb,parameters)
See Multipoint Boundary Value Problems for a description of
A structure containing the
initial guess for a solution. You create
Ordered nodes of the initial mesh. Boundary conditions
are imposed at a =
Initial guess for the solution such that
The structure can have any
name, but the fields must be named
Optional integration argument.
A structure you create using the
sol = bvp4c(odefun,bcfun,solinit) integrates
a system of ordinary differential equations of the form
y′ = f(x,y)
on the interval [a,b] subject to two-point boundary value conditions bc(y(a),y(b)) = 0.
bcfun are function
handles. See Create Function Handle for more information.
Parameterizing Functions explains how to provide additional
parameters to the function
odefun, as well as the
boundary condition function
bcfun, if necessary.
bvp4c can also solve multipoint boundary
value problems. See Multipoint Boundary Value Problems. You can use the function
bvpinit to specify the boundary points,
which are stored in the input argument
See the reference page for
bvpinit for more information.
bvp4c solver can also find unknown
parameters p for problems of the form
y′ = f(x,y, p)
0 = bc(y(a),y(b),p)
where p corresponds to
bvp4c an initial guess for any unknown
returns the final values of these unknown parameters in
bvp4c produces a solution that is continuous
[a,b] and has a continuous first derivative
there. Use the function
evaluate the solution at specific points
sxint = deval(sol,xint)
sol returned by
the following fields:
Mesh selected by
Approximation to y(x) at
the mesh points of
Approximation to y′(x) at
the mesh points of
Values returned by
|Computational cost statistics (also displayed when the |
sol can have any name, and
sol = bvp4c(odefun,bcfun,solinit,options) solves
as above with default integration properties replaced by the values
options, a structure created with the
bvpset for details.
solinit = bvpinit(x, yinit,
params) forms the initial guess
params of guesses for the unknown parameters.
bvp4c solves a class of singular boundary
value problems, including problems with unknown parameters
of the form
y′ = S · y/x + F(x,y,p)
0 = bc(y(0),y(b),p)
The interval is required to be [0, b] with b >
0. Often such problems arise when computing a smooth solution of ODEs
that result from partial differential equations (PDEs) due to cylindrical
or spherical symmetry. For singular problems, you specify the (constant)
S as the value of the
The boundary conditions must be consistent with the necessary condition S·y(0) = 0 and the initial
guess should satisfy this condition.
bvp4c can solve multipoint boundary value
problems where a = a0 < a1 < a2 < ...< an = b are
boundary points in the interval [a,b].
The points a1,a2,...,an−1 represent
interfaces that divide [a,b]
bvp4c enumerates the regions from
left to right (from a to b),
with indices starting from 1. In region k, [ak−1,ak],
the derivative as
yp = odefun(x,y,k)
In the boundary conditions function
yleft(:,k) is the solution at the left boundary
yright(:,k) is the solution at the right
boundary of region k. In particular,
yleft(:,1) = y(a)
yright(:,end) = y(b)
When you create an initial guess with
solinit = bvpinit(xinit,yinit),
use double entries in
xinit for each interface
point. See the reference page for
yinit is a function,
= yinit(x,k) to get an initial guess for the solution at
k. In the solution structure
sol.x has double entries
for each interface point. The corresponding columns of
the left and right solution at the interface, respectively.
To see an example that solves a three-point boundary value problem,
threebvp at the MATLAB® command prompt.
bvp5c function is
used exactly like
bvp4c, with the exception of
the meaning of error tolerances between the two solvers. If S(x)
approximates the solution y(x),
the residual |S′(x) – f(x,S(x))|.
This controls indirectly the true error |y(x)
the true error directly.
Boundary value problems can have multiple solutions. One purpose of the initial guess is to indicate which solution, among several, that you want.
The second-order differential equation
has exactly two solutions that satisfy the boundary conditions
bvp4c to solve the problem, you need to rewrite the differential equation as a system of two first-order ODEs,
where and . This system has the required form
The function, f, and the boundary conditions, bc, are coded in MATLAB as the functions
function dydx = twoode(x,y) dydx = [ y(2); -abs(y(1)) ];
function res = twobc(ya,yb) res = [ ya(1); yb(1) + 2 ];
Form a guess structure consisting of an initial mesh of five equally spaced points in [0,4] and a guess of the constant values
solinit = bvpinit(linspace(0,4,5),[1 0]);
Solve the problem using
sol = bvp4c(@twoode,@twobc,solinit);
Evaluate the numerical solution at 100 equally spaced points and plot y(x).
x = linspace(0,4); y = deval(sol,x); plot(x,y(1,:))
To obtain the other solution of this problem, use the initial guess
solinit = bvpinit(linspace(0,4,5),[-1 0]); sol = bvp4c(@twoode,@twobc,solinit); x = linspace(0,4); y = deval(sol,x); plot(x,y(1,:))
This boundary value problem involves an unknown parameter. The task is to compute the fourth (q = 5) eigenvalue λ of Mathieu's equation
y” + (λ – 2q cos2x)y = 0.
Because the unknown parameter λ is present, this second-order differential equation is subject to three boundary conditions:
y′(0) = 0
y′(π) = 0
y(0) = 1
It is convenient to use local functions to place all the functions
bvp4c in a single file.
function mat4bvp lambda = 15; solinit = bvpinit(linspace(0,pi,10),@mat4init,lambda); sol = bvp4c(@mat4ode,@mat4bc,solinit); fprintf('The fourth eigenvalue is approximately %7.3f.\n',... sol.parameters) xint = linspace(0,pi); Sxint = deval(sol,xint); plot(xint,Sxint(1,:)) axis([0 pi -1 1.1]) title('Eigenfunction of Mathieu''s equation.') xlabel('x') ylabel('solution y') % ------------------------------------------------------------ function dydx = mat4ode(x,y,lambda) q = 5; dydx = [ y(2) -(lambda - 2*q*cos(2*x))*y(1) ]; % ------------------------------------------------------------ function res = mat4bc(ya,yb,lambda) res = [ ya(2) yb(2) ya(1)-1 ]; % ------------------------------------------------------------ function yinit = mat4init(x) yinit = [ cos(4*x) -4*sin(4*x) ];
The differential equation (converted to a first-order system)
and the boundary conditions are coded as local functions
respectively. Because unknown parameters are present, these functions
must accept three input arguments, even though some of the arguments
are not used.
The guess structure
solinit is formed with
bvpinit. An initial guess for the solution
is supplied in the form of a function
We chose y = cos 4x because
it satisfies the boundary conditions and has the correct qualitative
behavior (the correct number of sign changes). In the call to
the third argument (
lambda = 15) provides an initial
guess for the unknown parameter λ.
After the problem is solved with
sol.parameters returns the value λ
= 17.097, and the plot shows the eigenfunction
associated with this eigenvalue.
bvp4c is a finite difference code that
implements the three-stage Lobatto IIIa formula. This is a collocation
formula and the collocation polynomial provides a C1-continuous
solution that is fourth-order accurate uniformly in
Mesh selection and error control are based on the residual of the
 Shampine, L.F., M.W. Reichelt, and J. Kierzenka, “Solving Boundary Value Problems for Ordinary Differential Equations in MATLAB with bvp4c”