Integer type:  int32  int64  nag_int  show int32  show int32  show int64  show int64  show nag_int  show nag_int

Chapter Contents
Chapter Introduction
NAG Toolbox

# NAG Toolbox: nag_ode_bvp_ps_lin_cgl_vals (d02ub)

## Purpose

nag_ode_bvp_ps_lin_cgl_vals (d02ub) evaluates a function, or one of its lower order derivatives, from its Chebyshev series representation at Chebyshev Gauss–Lobatto points on [a,b]$\left[a,b\right]$. The coefficients of the Chebyshev series representation required are usually derived from those returned by nag_ode_bvp_ps_lin_coeffs (d02ua) or nag_ode_bvp_ps_lin_solve (d02ue).

## Syntax

[f, ifail] = d02ub(n, a, b, q, c)
[f, ifail] = nag_ode_bvp_ps_lin_cgl_vals(n, a, b, q, c)

## Description

nag_ode_bvp_ps_lin_cgl_vals (d02ub) evaluates the Chebyshev series
 S (x) = (1/2) c1 T0 (x) + c2 T1 (x) + c3T2 (x) + ⋯ + cn + 1 Tn (x) , $S (x-) = 12 c1 T0 (x-) + c2 T1 (x-) + c3T2 (x-) +⋯+ cn+1 Tn (x-) ,$
or its derivative (up to fourth order) at the Chebyshev Gauss–Lobatto points on [a,b]$\left[a,b\right]$. Here Tj(x)${T}_{j}\left(\stackrel{-}{x}\right)$ denotes the Chebyshev polynomial of the first kind of degree j$j$ with argument x$\stackrel{-}{x}$ defined on [1,1]$\left[-1,1\right]$. In terms of your original variable, x$x$ say, the input values at which the function values are to be provided are
 xr = − (1/2) (b − a) cos(π(r − 1) / n) + (1/2) (b + a) ,   r = 1,2, … ,n + 1 , ​ $xr = - 12 ( b - a ) cos( π(r-1) /n ) + 1 2 ( b + a ) , r=1,2,…,n+1 , ​$
where b$b$ and a$a$ are respectively the upper and lower ends of the range of x$x$ over which the function is required.
The calculation is implemented by a forward one-dimensional discrete Fast Fourier Transform (DFT).

## References

Canuto C (1988) Spectral Methods in Fluid Dynamics 502 Springer
Canuto C, Hussaini M Y, Quarteroni A and Zang T A (2006) Spectral Methods: Fundamentals in Single Domains Springer
Trefethen L N (2000) Spectral Methods in MATLAB SIAM

## Parameters

### Compulsory Input Parameters

1:     n – int64int32nag_int scalar
n$n$, where the number of grid points is n + 1$n+1$. This is also the largest order of Chebyshev polynomial in the Chebyshev series to be computed.
Constraint: n > 0${\mathbf{n}}>0$ and n is even.
2:     a – double scalar
a$a$, the lower bound of domain [a,b]$\left[a,b\right]$.
Constraint: a < b${\mathbf{a}}<{\mathbf{b}}$.
3:     b – double scalar
b$b$, the upper bound of domain [a,b]$\left[a,b\right]$.
Constraint: b > a${\mathbf{b}}>{\mathbf{a}}$.
4:     q – int64int32nag_int scalar
The order, q$q$, of the derivative to evaluate.
Constraint: 0q4$0\le {\mathbf{q}}\le 4$.
5:     c(n + 1${\mathbf{n}}+1$) – double array
The Chebyshev coefficients, ci${c}_{\mathit{i}}$, for i = 1,2,,n + 1$\mathit{i}=1,2,\dots ,n+1$.

None.

None.

### Output Parameters

1:     f(n + 1${\mathbf{n}}+1$) – double array
The derivatives S(q) xi ${S}^{\left(q\right)}{x}_{\mathit{i}}$, for i = 1,2,,n + 1$\mathit{i}=1,2,\dots ,n+1$, of the Chebyshev series, S$S$.
2:     ifail – int64int32nag_int scalar
${\mathrm{ifail}}={\mathbf{0}}$ unless the function detects an error (see [Error Indicators and Warnings]).

## Error Indicators and Warnings

Errors or warnings detected by the function:
ifail = 1${\mathbf{ifail}}=1$
Constraint: n > 0${\mathbf{n}}>0$.
Constraint: n is even.
ifail = 2${\mathbf{ifail}}=2$
Constraint: a < b${\mathbf{a}}<{\mathbf{b}}$.
ifail = 3${\mathbf{ifail}}=3$
Constraint: 0q4$0\le {\mathbf{q}}\le 4$.
ifail = 999${\mathbf{ifail}}=-999$
Dynamic memory allocation failed.

## Accuracy

Evaluations of DFT to obtain function or derivative values should be an order n$n$ multiple of machine precision assuming full accuracy to machine precision in the given Chebyshev series representation.

The number of operations is of the order n log(n) $n\mathrm{log}\left(n\right)$ and the memory requirements are O(n)$\mathit{O}\left(n\right)$; thus the computation remains efficient and practical for very fine discretizations (very large values of n$n$).

## Example

```function nag_ode_bvp_ps_lin_cgl_vals_example
n = int64(16);
a = -pi/2;
b =  pi/2;

% Set up boundary condition on left side of domain
y = [a, a, b];
% Set up Dirichlet condition using exact solution at x=a.
bmat = zeros(3, 4);
bmat(1, 1) = 1;
bmat(2, 1:3) = [1, 2, 3];
bmat(3, 1:3) = [1, 2, 3];
bvec = [0, 2, -2];

% Set up problem definition
f = [1, 2, 3, 4];

% Set up solution grid
[x, ifail] = nag_ode_bvp_ps_lin_cgl_grid(n, a, b);

% Set up problem right hand sides for grid and transform
f0 = 2*sin(x) - 2*cos(x);
[c, ifail] = nag_ode_bvp_ps_lin_coeffs(n, f0);

% Solve in coefficient space
[bmat, f, uc, resid, ifail] = nag_ode_bvp_ps_lin_solve(n, a, b, c, bmat, y, bvec, f);

% Transform solution and derivative back to real space.
u = zeros(17, 4);
for q=0:3
[u(:, q+1),  ifail] = nag_ode_bvp_ps_lin_cgl_vals(n, a, b, int64(q), uc(:, q+1));
end

% Print Solution
fprintf('\nNumerical solution U and its first three derivatives\n');
fprintf('      x          U          Ux         Uxx       Uxxx\n');
for i=1:17
fprintf('%10.4f %10.4f %10.4f %10.4f %10.4f\n', x(i), u(i, :));
end
```
```

Numerical solution U and its first three derivatives
x          U          Ux         Uxx       Uxxx
-1.5708    -0.0000     1.0000     0.0000    -1.0000
-1.5406     0.0302     0.9995    -0.0302    -0.9995
-1.4512     0.1193     0.9929    -0.1193    -0.9929
-1.3061     0.2616     0.9652    -0.2616    -0.9652
-1.1107     0.4440     0.8960    -0.4440    -0.8960
-0.8727     0.6428     0.7661    -0.6428    -0.7661
-0.6011     0.8247     0.5656    -0.8247    -0.5656
-0.3064     0.9534     0.3017    -0.9534    -0.3017
-0.0000     1.0000     0.0000    -1.0000    -0.0000
0.3064     0.9534    -0.3017    -0.9534     0.3017
0.6011     0.8247    -0.5656    -0.8247     0.5656
0.8727     0.6428    -0.7661    -0.6428     0.7661
1.1107     0.4440    -0.8960    -0.4440     0.8960
1.3061     0.2616    -0.9652    -0.2616     0.9652
1.4512     0.1193    -0.9929    -0.1193     0.9929
1.5406     0.0302    -0.9995    -0.0302     0.9995
1.5708    -0.0000    -1.0000     0.0000     1.0000

```
```function d02ub_example
n = int64(16);
a = -pi/2;
b =  pi/2;

% Set up boundary condition on left side of domain
y = [a, a, b];
% Set up Dirichlet condition using exact solution at x=a.
bmat = zeros(3, 4);
bmat(1, 1) = 1;
bmat(2, 1:3) = [1, 2, 3];
bmat(3, 1:3) = [1, 2, 3];
bvec = [0, 2, -2];

% Set up problem definition
f = [1, 2, 3, 4];

% Set up solution grid
[x, ifail] = d02uc(n, a, b);

% Set up problem right hand sides for grid and transform
f0 = 2*sin(x) - 2*cos(x);
[c, ifail] = d02ua(n, f0);

% Solve in coefficient space
[bmat, f, uc, resid, ifail] = d02ue(n, a, b, c, bmat, y, bvec, f);

% Transform solution and derivative back to real space.
u = zeros(17, 4);
for q=0:3
[u(:, q+1),  ifail] = d02ub(n, a, b, int64(q), uc(:, q+1));
end

% Print Solution
fprintf('\nNumerical solution U and its first three derivatives\n');
fprintf('      x          U          Ux         Uxx       Uxxx\n');
for i=1:17
fprintf('%10.4f %10.4f %10.4f %10.4f %10.4f\n', x(i), u(i, :));
end
```
```

Numerical solution U and its first three derivatives
x          U          Ux         Uxx       Uxxx
-1.5708    -0.0000     1.0000     0.0000    -1.0000
-1.5406     0.0302     0.9995    -0.0302    -0.9995
-1.4512     0.1193     0.9929    -0.1193    -0.9929
-1.3061     0.2616     0.9652    -0.2616    -0.9652
-1.1107     0.4440     0.8960    -0.4440    -0.8960
-0.8727     0.6428     0.7661    -0.6428    -0.7661
-0.6011     0.8247     0.5656    -0.8247    -0.5656
-0.3064     0.9534     0.3017    -0.9534    -0.3017
-0.0000     1.0000     0.0000    -1.0000    -0.0000
0.3064     0.9534    -0.3017    -0.9534     0.3017
0.6011     0.8247    -0.5656    -0.8247     0.5656
0.8727     0.6428    -0.7661    -0.6428     0.7661
1.1107     0.4440    -0.8960    -0.4440     0.8960
1.3061     0.2616    -0.9652    -0.2616     0.9652
1.4512     0.1193    -0.9929    -0.1193     0.9929
1.5406     0.0302    -0.9995    -0.0302     0.9995
1.5708    -0.0000    -1.0000     0.0000     1.0000

```