NAG Toolbox: nag_ode_withdraw_ivp_rk_interp (d02px)

Purpose

Syntax

Description

References

Parameters

Compulsory Input Parameters

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

$n$, the number of ordinary differential equations in the system to be solved by the integration function.

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

2:     $\mathrm{twant}$ – double scalar

$t$, the value of the independent variable where a solution is desired.

3:     $\mathrm{reqest}$ – string (length ≥ 1)

Determines whether the solution and/or its first derivative are to be computed.

$\mathbf{reqest}=\text{'S'}$

Compute the approximate solution only.

$\mathbf{reqest}=\text{'D'}$

Compute the approximate first derivative of the solution only.

$\mathbf{reqest}=\text{'B'}$

Compute both the approximate solution and its first derivative.

Constraint: $\mathbf{reqest}=\text{'S'}$, $\text{'D'}$ or $\text{'B'}$.

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

The number of components of the solution to be computed. The first nwant components are evaluated.

Constraint: $1\le \mathbf{nwant}\le \mathit{n}$, where $\mathit{n}$ is specified by neq in the prior call to nag_ode_ivp_rk_setup (d02pv).

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

f must evaluate the functions $f_i$ (that is the first derivatives $y_i^{\prime }$) for given values of the arguments $t,y_i$. It must be the same procedure as supplied to nag_ode_ivp_rk_onestep (d02pd).

[yp] = f(t, y)

Input Parameters

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

$t$, the current value of the independent variable.

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

The current values of the dependent variables, $y_{\mathit{i}}$, for $\mathit{i}=1,2,\dots ,\mathit{n}$.

Output Parameters

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

The values of $f_{\mathit{i}}$, for $\mathit{i}=1,2,\dots ,\mathit{n}$.

6:     $\mathrm{work}\left(:\right)$ – double array

The dimension of the array work must be at least $\mathbf{lenwrk}$ (see nag_ode_ivp_rk_setup (d02pv))

This must be the same array as supplied to nag_ode_ivp_rk_onestep (d02pd) and must remain unchanged between calls.

7:     $\mathrm{wrkint}\left(\mathbf{lenint}\right)$ – double array

Must be the same array as supplied in previous calls, if any, and must remain unchanged between calls to nag_ode_ivp_rk_interp (d02px).

Optional Input Parameters

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

Default: the dimension of the array wrkint.

The dimension of the array wrkint.

Constraints:

• $\mathbf{lenint}\ge 1$ if $\mathbf{method}=1$ in the prior call to nag_ode_ivp_rk_setup (d02pv);
• $\mathbf{lenint}\ge \mathit{n}+5\times \mathbf{nwant}$ if $\mathbf{method}=2$ and $\mathit{n}$ is specified by neq in the prior call to nag_ode_ivp_rk_setup (d02pv).

Output Parameters

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

The dimension of the array ywant will be $\mathbf{nwant}$ if $\mathbf{reqest}=\text{'S'}$ or $\text{'B'}$ and $1$ otherwise

An approximation to the first nwant components of the solution at twant if $\mathbf{reqest}=\text{'S'}$ or $\text{'B'}$. Otherwise ywant is not defined.

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

The dimension of the array ypwant will be $\mathbf{nwant}$ if $\mathbf{reqest}=\text{'D'}$ or $\text{'B'}$ and $1$ otherwise

An approximation to the first nwant components of the first derivative at twant if $\mathbf{reqest}=\text{'D'}$ or $\text{'B'}$. Otherwise ypwant is not defined.

3:     $\mathrm{work}\left(:\right)$ – double array

The dimension of the array work will be $\mathbf{lenwrk}$ (see nag_ode_ivp_rk_setup (d02pv))

Contains information about the integration for use on subsequent calls to nag_ode_ivp_rk_onestep (d02pd) or other associated functions.

4:     $\mathrm{wrkint}\left(\mathbf{lenint}\right)$ – double array

The contents are modified.

5:     $\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$

On entry, an invalid input value for nwant or lenint was detected or an invalid call to nag_ode_ivp_rk_interp (d02px) was made, for example without a previous call to the integration function nag_ode_ivp_rk_onestep (d02pd), or after an error return from nag_ode_ivp_rk_onestep (d02pd), or if nag_ode_ivp_rk_onestep (d02pd) was being used with $\mathbf{method}=3$. You cannot continue integrating the problem.

   $\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: d02px_example

function d02px_example


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

% Set initial conditions and input
method = int64(2);
tstart = 0;
tend   = 2*pi;
yinit  = [0;1];
hstart = 0;
thresh = [1e-08; 1e-08];
n      = int64(2);
task   = 'Complex Task';
errass = false;
lenwrk = int64(32*n);
reqest = 'Both';
nwant  = int64(1);
tol    =  1.0E-2;
npts   = 32;
work   = zeros(n+5*nwant, 1);
twant  = zeros(npts+1, 1);
ywant  = zeros(npts, nwant);
ypwant = zeros(npts, nwant);
lenint = n+5*nwant;
wrkint = zeros(lenint,1);

% Set control for printing solution
tinc   = (tend-tstart)/npts;

% Run through the calculation twice with two tolerance values
for i = 1:2

  tol = tol*0.1;

  % Call setup function
  [work, ifail] = d02pv(tstart, yinit, tend, tol, thresh, method, task, ...
                        errass, lenwrk);

  % Set up first point at which solution is required
  twant(1) = tstart;
  twant(2) = tstart + tinc;
  ywant(1,1) = yinit(1);
  ypwant(1,1) = yinit(2);
  t = tstart;
  j=1;

  % Integrate by steps until tend is reached or error is encountered.
  while t < tend

  % Integrate one step from t, updating t.
    [t, y, yp, work, ifail] = d02pd(@f, n, work);

    % Interpolate at required additional points up to t
    while twant(j+1) < t
      j=j+1;
      % Interpolate and print solution at t = twant.
      [ywant(j,:), ypwant(j,:), work, wrkint, ifail] = d02px(n, twant(j), ...
                                             reqest, nwant, @f, work, wrkint);

      % Set next required solution point
      twant(j+1) = twant(j) + tinc;
    end
  end

  % Get integration statistics.
  [fevals, stepcost, waste, stepsok, hnext, fail] = d02py;
  fprintf('\nCalculation with TOL = %8.1e\n\n', tol);
  fprintf('     Number of evaluations of f = %d\n', fevals);

end

% Plot results
fig1 = figure;
title('Simple Sine Solution, tol=1e-4');
hold on;
xlabel('t');
ylabel('Solution');
plot(twant(1:npts), ywant(:, 1), '-xr');
text(3, 0.25, 'Solution', 'Color', 'r');
plot(twant(1:npts), ypwant(:, 1), '-xg');
text(1.45, 0.25, 'Derivative', 'Color', 'g');
hold off



function [yp] = f(t, y)
% Evaluate derivative vector.
yp = zeros(2, 1);
yp(1) =  y(2);
yp(2) = -y(1);

d02px example results


Calculation with TOL =  1.0e-03

     Number of evaluations of f = 68

Calculation with TOL =  1.0e-04

     Number of evaluations of f = 105