NAG Toolbox: nag_ode_ivp_rkts_interp (d02ps)

Purpose

Syntax

Description

References

Parameters

Compulsory Input Parameters

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

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

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

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

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

3:     $\mathrm{ideriv}$ – int64int32nag_int scalar

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

$\mathbf{ideriv}=0$

compute approximate solution.

$\mathbf{ideriv}=1$

compute approximate first derivative.

$\mathbf{ideriv}=2$

compute approximate solution and first derivative.

Constraint: $\mathbf{ideriv}=0$, $1$ or $2$.

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 \mathbf{n}$.

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_rkts_onestep (d02pf).

[yp, user] = f(t, n, y, user)

Input Parameters

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

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

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

$\mathit{n}$, the number of ordinary differential equations in the system to be solved.

3:     $\mathrm{y}\left(\mathbf{n}\right)$ – double array

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

4:     $\mathrm{user}$ – Any MATLAB object

f is called from nag_ode_ivp_rkts_interp (d02ps) with the object supplied to nag_ode_ivp_rkts_interp (d02ps).

Output Parameters

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

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

2:     $\mathrm{user}$ – Any MATLAB object

6:     $\mathrm{wcomm}\left(\mathbf{lwcomm}\right)$ – double array

This array stores information that can be utilized on subsequent calls to nag_ode_ivp_rkts_interp (d02ps).

7:     $\mathrm{iwsav}\left(130\right)$ – int64int32nag_int array

8:     $\mathrm{rwsav}\left(32\times \mathbf{n}+350\right)$ – double array

These must be the same arrays supplied in a previous call nag_ode_ivp_rkts_onestep (d02pf). They must remain unchanged between calls.

Optional Input Parameters

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

Default: the dimension of the array wcomm.

Length of wcomm.

If in a previous call to nag_ode_ivp_rkts_setup (d02pq):

• $\mathbf{method}=1$ or $-1$ then lwcomm must be at least $1$.
• $\mathbf{method}=2$ or $-2$ then lwcomm must be at least $\mathbf{n}+\max \left(\mathbf{n},5\times \mathbf{nwant}\right)$.
• $\mathbf{method}=3$ or $-3$ then wcomm and lwcomm are not referenced.

2:     $\mathrm{user}$ – Any MATLAB object

user is not used by nag_ode_ivp_rkts_interp (d02ps), but is passed to f. Note that for large objects it may be more efficient to use a global variable which is accessible from the m-files than to use user.

Output Parameters

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

An approximation to the first nwant components of the solution at twant if $\mathbf{ideriv}=0$ or $2$. Otherwise ywant is not defined.

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

An approximation to the first nwant components of the first derivative at twant if $\mathbf{ideriv}=1$ or $2$. Otherwise ypwant is not defined.

3:     $\mathrm{wcomm}\left(\mathbf{lwcomm}\right)$ – double array

Communication array, used to store information between calls to nag_ode_ivp_rkts_interp (d02ps).

4:     $\mathrm{user}$ – Any MATLAB object

5:     $\mathrm{iwsav}\left(130\right)$ – int64int32nag_int array

6:     $\mathrm{rwsav}\left(32\times \mathbf{n}+350\right)$ – double array

Information about the integration for use on subsequent calls to nag_ode_ivp_rkts_onestep (d02pf), nag_ode_ivp_rkts_interp (d02ps) or other associated functions.

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

Constraint: $1\le \mathbf{nwant}\le \mathbf{n}$.

Constraint: for $\mathbf{method}=-1$ or $1$, $\mathbf{lwcomm}\ge 1$.

Constraint: for $\mathbf{method}=-2$ or $2$, $\mathbf{lwcomm}\ge \mathbf{n}+\max \left(\mathbf{n},5\times \mathbf{nwant}\right)$.

Constraint: $\mathbf{ideriv}=0$, $1$ or $2$.

$\mathbf{method}=-3$ or $3$ in setup, but interpolation is not available for this method. Either use $\mathbf{method}=-2$ or $2$ in setup or use reset function to force the integrator to step to particular points.

On entry, a previous call to the setup function has not been made or the communication arrays have become corrupted, or a catastrophic error has already been detected elsewhere.
You cannot continue integrating the problem.

On entry, $\mathbf{n}=\_$, but the value passed to the setup function was .

You cannot call this function after the integrator has returned an error.

You cannot call this function before you have called the step integrator.

You cannot call this function when you have specified, in the setup function, that the range integrator will be used.

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

function d02ps_example


fprintf('d02ps 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);
nwant  = int64(1);
tol    =  1.0E-2;
npts   = 32;
wcomm  = zeros(n+5*nwant, 1);
twant  = zeros(npts+1, 1);
ywant  = zeros(npts, nwant);
ypwant = zeros(npts, nwant);

% 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
  [iwsav, rwsav, ifail] = d02pq(tstart, tend, yinit, tol, thresh, method);

  % 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, user, iwsav, rwsav, ifail] = d02pf(@f, n, iwsav, rwsav);

    % Interpolate at required additional points up to t
    while twant(j+1) < t
      j=j+1;
      % Interpolate and print solution at t = twant.
      ideriv = int64(2);
      [ywant(j, :), ypwant(j, :), wcomm, user, iwsav, rwsav, ifail] = ...
        d02ps(n, twant(j), ideriv, nwant, @f, wcomm, iwsav, rwsav);

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

  % Get integration statistics.
  [fevals, stepcost, waste, stepsok, hnext, iwsav, ifail] = d02pt(iwsav,rwsav);
  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, user] = f(t, n, y, user)
  yp = [y(2); -y(1)];

d02ps example results


Calculation with TOL =  1.0e-03

     Number of evaluations of f = 152

Calculation with TOL =  1.0e-04

     Number of evaluations of f = 231