NAG Toolbox: nag_ode_bvp_shoot_bval (d02ha)

Purpose

Syntax

Description

References

Parameters

Compulsory Input Parameters

1:     $\mathrm{u}\left(\mathbf{n},2\right)$ – double array

$\mathbf{u}\left(\mathit{i},1\right)$ must be set to the known or estimated value of $y_{\mathit{i}}$ at $a$ and $\mathbf{u}\left(\mathit{i},2\right)$ must be set to the known or estimated value of $y_{\mathit{i}}$ at $b$, for $\mathit{i}=1,2,\dots ,\mathit{n}$.

2:     $\mathrm{v}\left(\mathbf{n},2\right)$ – double array

$\mathbf{v}\left(\mathit{i},\mathit{j}\right)$ must be set to $0.0$ if $\mathbf{u}\left(\mathit{i},\mathit{j}\right)$ is a known value and to $1.0$ if $\mathbf{u}\left(\mathit{i},\mathit{j}\right)$ is an estimated value, for $\mathit{i}=1,2,\dots ,\mathit{n}$ and $\mathit{j}=1,2$.

Constraint: precisely $\mathit{n}$ of the $\mathbf{v}\left(i,j\right)$ must be set to $0.0$, i.e., precisely $\mathit{n}$ of the $\mathbf{u}\left(i,j\right)$ must be known values, and these must not be all at $a$ or all at $b$.

3:     $\mathrm{a}$ – double scalar

$a$, the initial point of the interval of integration.

4:     $\mathrm{b}$ – double scalar

$b$, the final point of the interval of integration.

5:     $\mathrm{tol}$ – double scalar

Must be set to a small quantity suitable for:

(a)	testing the local error in $y_i$ during integration,
(b)	testing for the convergence of $y_i$ at $b$,
(c)	calculating the perturbation in estimated boundary values for $y_i$, which are used to obtain the approximate derivatives of the residuals for use in the Newton iteration.

You are advised to check your results by varying tol.

Constraint: $\mathbf{tol}>0.0$.

6:     $\mathrm{fcn}$ – function handle or string containing name of m-file

fcn must evaluate the functions $f_{\mathit{i}}$ (i.e., the derivatives $y_{\mathit{i}}^{\prime }$), for $\mathit{i}=1,2,\dots ,\mathit{n}$, at a general point $x$.

[f] = fcn(x, y)

Input Parameters

1:     $\mathrm{x}$ – double scalar

$x$, the value of the argument.

2:     $\mathrm{y}\left(:\right)$ – double array

$y_{\mathit{i}}$, for $\mathit{i}=1,2,\dots ,\mathit{n}$, the value of the argument.

Output Parameters

1:     $\mathrm{f}\left(:\right)$ – double array

The values of $f_{\mathit{i}}\left(x\right)$, for $\mathit{i}=1,2,\dots ,\mathit{n}$.

7:     $\mathrm{m1}$ – int64int32nag_int scalar

A value which controls output.

$\mathbf{m1}=1$

The final solution is not evaluated.

$\mathbf{m1}>1$

The final values of $y_{\mathit{i}}$ at interval $\left(b-a\right)/\left(\mathbf{m1}-1\right)$ are calculated and stored in the array soln by columns, starting with values $y_{\mathit{i}}$ at $a$ stored in $\mathbf{soln}\left(\mathit{i},1\right)$, for $\mathit{i}=1,2,\dots ,\mathit{n}$.

Constraint: $\mathbf{m1}\ge 1$.

Optional Input Parameters

1:     $\mathrm{n}$ – int64int32nag_int scalar

Default: the first dimension of the arrays u, v. (An error is raised if these dimensions are not equal.)

$\mathit{n}$, the number of equations.

Constraint: $\mathbf{n}\ge 1$.

Output Parameters

1:     $\mathrm{u}\left(\mathbf{n},2\right)$ – double array

The known values unaltered, and corrected values of the estimates, unless an error has occurred. If an error has occurred, u contains the known values and the latest values of the estimates.

2:     $\mathrm{soln}\left(\mathbf{n},\mathbf{m1}\right)$ – double array

The solution when $\mathbf{m1}>1$.

3:     $\mathrm{w}\left(\mathbf{n},\mathit{sdw}\right)$ – double array

$\mathit{sdw}=3\mathbf{n}+17+\max \left(11,\mathbf{n}\right)$.

If $\mathbf{ifail}=\mathbf{2}$, $\mathbf{3}$, $\mathbf{4}$ or $\mathbf{5}$, $\mathbf{w}\left(\mathit{i},1\right)$, for $\mathit{i}=1,2,\dots ,\mathit{n}$, contains the solution at the point where the integration fails and the point of failure is returned in $\mathbf{w}\left(1,2\right)$.

4:     $\mathrm{ifail}$ – int64int32nag_int scalar

$\mathbf{ifail}=\mathbf{0}$ unless the function detects an error (see Error Indicators and Warnings).

Error Indicators and Warnings

   $\mathbf{ifail}=1$

One or more of the arguments v, n, m1, sdw, or tol is incorrectly set.

   $\mathbf{ifail}=2$

The step length for the integration is too short whilst calculating the residual (see Further Comments).

   $\mathbf{ifail}=3$

No initial step length could be chosen for the integration whilst calculating the residual.

   $\mathbf{ifail}=4$

As for $\mathbf{ifail}=\mathbf{2}$ but the error occurred when calculating the Jacobian of the derivatives of the residuals with respect to the parameters.

   $\mathbf{ifail}=5$

As for $\mathbf{ifail}=\mathbf{3}$ but the error occurred when calculating the derivatives of the residuals with respect to the parameters.

   $\mathbf{ifail}=6$

The calculated Jacobian has an insignificant column.

   $\mathbf{ifail}=7$

The linear algebra function (nag_lapack_dgesvd (f08kb)) used has failed. This error exit should not occur and can be avoided by changing the estimated initial values.

   $\mathbf{ifail}=8$

The Newton iteration has failed to converge.

   $\mathbf{ifail}=9$

   $\mathbf{ifail}=10$

   $\mathbf{ifail}=11$

   $\mathbf{ifail}=12$

   $\mathbf{ifail}=13$

Indicates that a serious error has occurred in an internal call. Check all array subscripts and function argument lists in calls to nag_ode_bvp_shoot_bval (d02ha). Seek expert help.

   $\mathbf{ifail}=-99$

An unexpected error has been triggered by this routine. Please contact NAG.

   $\mathbf{ifail}=-399$

Your licence key may have expired or may not have been installed correctly.

   $\mathbf{ifail}=-999$

Dynamic memory allocation failed.

Accuracy

Further Comments

Example

Open in the MATLAB editor: d02ha_example

function d02ha_example


fprintf('d02ha example results\n\n');

% Initialize variables and arrays.
u = [0, 0; 0.5, 0.46; 1.15, -1.2];
v = [0, 0; 0, 1; 1, 1];
a = 0;
b = 5;
m1 = 6;
n = 3;
ykeep = zeros(m1,3);
xkeep = zeros(1,m1);

for id = 3:4
    tol = 5*10^(-id);

    [uOut, soln, w, ifail] = d02ha(u, v, a, b, tol, @fcn, int64(m1));
    if (ifail ~= 0)
        % Problems in integration, or incorrect parameters.  Print
        % message and exit.
        error('Warning: d02ha returned with ifail = %d ',ifail);
    else
        % Output results.
        fprintf('Results with tol = %1.3e\n\n',tol);
        fprintf('X-value and final solution\n\n');
        for i = 1:m1
            disp((sprintf('%d   %8.4f   %8.4f   %8.4f', i-1,...
                soln(1,i),soln(2,i),soln(3,i))));
            % Store results for plotting.
            ykeep(i,1) = soln(1,i);
            ykeep(i,2) = soln(2,i);
            ykeep(i,3) = soln(3,i);
            xkeep(i) = i-1;
        end
    end
    disp(' ');
end
% Plot results.
fig1 = figure;
display_plot(xkeep, ykeep)


function f = fcn(x,y)
% Evaluate the derivatives.
f = zeros(3,1);
f(1) = tan(y(3));
f(2) = -0.032*tan(y(3))/y(2) - 0.02*y(2)/cos(y(3));
f(3) = -0.032/y(2)^2;

function  display_plot(xkeep, ykeep)
% Plot curves.
plot(xkeep, ykeep(:,1), '-+',...
     xkeep, ykeep(:,2), '--x',...
     xkeep, ykeep(:,3), ':*');
% Add title.
title({'Solution of Two-point Boundary-value Problem',...
    ['using Runge-Kutta-Merson and Newton Correction in a ', ...
    'Shooting Method']});
% Label the axes.
xlabel('x');
ylabel('Solution');
% Add legend.
legend('height','angle','velocity','location','Best');

d02ha example results

Results with tol = 5.000e-03

X-value and final solution

0     0.0000     0.5000     1.1685
1     1.9197     0.3341     0.9752
2     2.9304     0.2066     0.4917
3     2.9767     0.1956    -0.4228
4     2.0258     0.3093    -0.9754
5    -0.0005     0.4598    -1.2020
 
Results with tol = 5.000e-04

X-value and final solution

0     0.0000     0.5000     1.1681
1     1.9177     0.3343     0.9749
2     2.9280     0.2070     0.4929
3     2.9769     0.1955    -0.4194
4     2.0210     0.3095    -0.9751
5    -0.0000     0.4597    -1.2014