hide long namesshow long names
hide short namesshow short names
Integer type:  int32  int64  nag_int  show int32  show int32  show int64  show int64  show nag_int  show nag_int

PDF version (NAG web site, 64-bit version, 64-bit version)
Chapter Contents
Chapter Introduction
NAG Toolbox

NAG Toolbox: nag_linsys_real_toeplitz_update (f04mf)

 Contents

    1  Purpose
    2  Syntax
    7  Accuracy
    9  Example

Purpose

nag_linsys_real_toeplitz_update (f04mf) updates the solution of the equations Tx=b, where T is a real symmetric positive definite Toeplitz matrix.

Syntax

[x, p, work, ifail] = f04mf(t, b, x, work, 'n', n)
[x, p, work, ifail] = nag_linsys_real_toeplitz_update(t, b, x, work, 'n', n)

Description

nag_linsys_real_toeplitz_update (f04mf) solves the equations
Tnxn=bn,  
where Tn is the n by n symmetric positive definite Toeplitz matrix
Tn= τ0 τ1 τ2 τn-1 τ1 τ0 τ1 τn-2 τ2 τ1 τ0 τn-3 . . . . τn-1 τn-2 τn-3 τ0  
and bn is the n-element vector bn=β1β2βnT, given the solution of the equations
Tn-1xn-1=bn-1.  
This function will normally be used to successively solve the equations
Tkxk=bk,   k= 1,2,,n.  
If it is desired to solve the equations for a single value of n, then function nag_linsys_real_toeplitz_solve (f04ff) may be called. This function uses the method of Levinson (see Levinson (1947) and Golub and Van Loan (1996)).

References

Bunch J R (1985) Stability of methods for solving Toeplitz systems of equations SIAM J. Sci. Statist. Comput. 6 349–364
Bunch J R (1987) The weak and strong stability of algorithms in numerical linear algebra Linear Algebra Appl. 88/89 49–66
Cybenko G (1980) The numerical stability of the Levinson–Durbin algorithm for Toeplitz systems of equations SIAM J. Sci. Statist. Comput. 1 303–319
Golub G H and Van Loan C F (1996) Matrix Computations (3rd Edition) Johns Hopkins University Press, Baltimore
Levinson N (1947) The Weiner RMS error criterion in filter design and prediction J. Math. Phys. 25 261–278

Parameters

Compulsory Input Parameters

1:     t: – double array
The dimension of the array t must be at least max1,n
ti+1 must contain the value τi, for i=0,1,,n-1.
Constraint: t1>0.0. Note that if this is not true, then the Toeplitz matrix cannot be positive definite.
2:     b: – double array
The dimension of the array b must be at least max1,n
The right-hand side vector bn.
3:     x: – double array
The dimension of the array x must be at least max1,n
With n>1 the (n-1) elements of the solution vector xn-1 as returned by a previous call to nag_linsys_real_toeplitz_update (f04mf). The element xn need not be specified.
4:     work: – double array
The dimension of the array work must be at least max1,2×n-1
With n>2 the elements of work should be as returned from a previous call to nag_linsys_real_toeplitz_update (f04mf) with (n-1) as the argument n.

Optional Input Parameters

1:     n int64int32nag_int scalar
Default: the dimension of the arrays t, b, x.
The order of the Toeplitz matrix T.
Constraint: n0. When n=0, then an immediate return is effected.

Output Parameters

1:     x: – double array
The dimension of the array x will be max1,n
The solution vector xn.
2:     p – double scalar
The reflection coefficient pn-1. (See Further Comments.)
3:     work: – double array
The dimension of the array work will be max1,2×n-1
The first (n-1) elements of work contain the solution to the Yule–Walker equations
Tn-1yn-1=-tn-1, 
where tn-1=τ1τ2τn-1t.
4:     ifail int64int32nag_int scalar
ifail=0 unless the function detects an error (see Error Indicators and Warnings).

Error Indicators and Warnings

Errors or warnings detected by the function:

Cases prefixed with W are classified as warnings and do not generate an error of type NAG:error_n. See nag_issue_warnings.

   ifail=-1
On entry,n<0,
ort00.0.
W  ifail=1
The Toeplitz matrix Tn is not positive definite to working accuracy. If, on exit, p is close to unity, then Tn was probably close to being singular.
   ifail=-99
An unexpected error has been triggered by this routine. Please contact NAG.
   ifail=-399
Your licence key may have expired or may not have been installed correctly.
   ifail=-999
Dynamic memory allocation failed.

Accuracy

The computed solution of the equations certainly satisfies
r=Tnxn-bn,  
where r1 is approximately bounded by
r1cεCTn,  
c being a modest function of n, ε being the machine precision and CT being the condition number of T with respect to inversion. This bound is almost certainly pessimistic, but it seems unlikely that the method of Levinson is backward stable, so caution should be exercised when Tn is ill-conditioned. The following bound on Tn-1 holds:
max1i=1 n-11-pi2 , 1i=1 n-11-pi Tn-11i=1 n-1 1+pi 1-pi .  
(See Golub and Van Loan (1996).) The norm of Tn-1 may also be estimated using function nag_linsys_real_gen_norm_rcomm (f04yd). For further information on stability issues see Bunch (1985), Bunch (1987), Cybenko (1980) and Golub and Van Loan (1996).

Further Comments

The number of floating-point operations used by this function is approximately 8n.
If yi is the solution of the equations
Tiyi=-τ1τ2τiT,  
then the reflection coefficient pi is defined as the ith element of yi.

Example

This example finds the solution of the equations Tkxk=bk, k=1,2,3,4, where
T4= 4 3 2 1 3 4 3 2 2 3 4 3 1 2 3 4   and  b4= 1 1 1 1 .  
function f04mf_example


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

t = [4; 3; 2; 1];
b = [1; 1; 1; 1];
x = [0];
work = zeros(9,1);
fprintf('  order  refl. coeff   solution\n');
for k=1:4
  [x, p, work, ifail] = f04mf( ...
                               t(1:k), b(1:k), x, work);
  if k > 1
    fprintf('%6d%11.4f    ', k, p);
  else
    fprintf('%6d%11.4s    ', k, ' ');
  end
  fprintf('%8.4f',transpose(x));
  fprintf('\n');
  if k < 4
    x = [x; 0]; % Extend x by one element
  end
end


f04mf example results

  order  refl. coeff   solution
     1                 0.2500
     2    -0.7500      0.1429  0.1429
     3     0.1429      0.1667  0.0000  0.1667
     4     0.1667      0.2000  0.0000 -0.0000  0.2000

PDF version (NAG web site, 64-bit version, 64-bit version)
Chapter Contents
Chapter Introduction
NAG Toolbox

© The Numerical Algorithms Group Ltd, Oxford, UK. 2009–2015