-
Notifications
You must be signed in to change notification settings - Fork 0
Expand file tree
/
Copy pathchebpts.m
More file actions
120 lines (107 loc) · 3.63 KB
/
chebpts.m
File metadata and controls
120 lines (107 loc) · 3.63 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
function [x, w, v] = chebpts(n, dom, type)
%CHEBPTS Chebyshev points.
% CHEBPTS(N) returns N Chebyshev points of the 2nd-kind in [-1,1].
%
% CHEBPTS(N, D), where D is vector of length 2 and N is a scalar integer,
% scales the nodes and weights for the interval [D(1),D(2)]. If length(D) > 2
% and N a vector of length(D)-1, then CHEBPTS(N, D) returns a column vector of
% the stacked N(k) Chebyshev points on the subintervals D(k:k+1). If length(N)
% is 1, then D is treated as [D(1),D(end)].
%
% [X, W] = CHEBPTS(N, D) returns also a row vector of the (scaled) weights for
% Clenshaw-Curtis quadrature (computed using [1]). (For nodes and weights of
% Gauss-Chebyshev quadrature, use [X, W] = JACPTS(N, -.5, -.5, D))
%
% [X, W, V] = CHEBPTS(N, D) returns, in addition to X and W, the barycentric
% weights V corresponding to the Chebyshev points X. The weights are scaled to
% have infinity norm 1.
%
% [X, W, V] = CHEBPTS(N, KIND) or CHEBPTS(N, D, KIND) returns Chebyshev points
% and weights of the 1st-kind if KIND = 1 and 2nd-kind if KIND = 2 (default).
%
% [1] Jarg Waldvogel, "Fast construction of the Fejer and Clenshaw-Curtis
% quadrature rules", BIT Numerical Mathematics, 46, (2006), pp 195-202.
%
% See also FOURPTS, LEGPTS, JACPTS, LAGPTS, and HERMPTS.
% Copyright 2014 by The University of Oxford and The Chebfun Developers.
% See http://www.chebfun.org/ for Chebfun information.
% Parse inputs:
if ( nargin == 2 )
if ( length(dom) == 1 )
type = dom;
dom = chebfunpref().domain;
else
type = 2;
end
end
if ( nargin == 1 )
dom = chebfunpref().domain;
type = 2;
end
% Verify that the number of points requested and the domains match:
if ( length(n) == 1 && length(dom) > 1 )
% Ignore interior breaks in this instance.
dom = dom([1,end]);
elseif ( length(n) ~= length(dom) - 1 )
error('CHEBFUN:chebpts:mismatchND', 'Vector N does not match domain D.');
end
% Create a dummy CHEBTECH of appropriate type to access static methods.
f = feval(['chebtech', num2str(type)]);
if ( length(n) == 1 ) % Single grid.
% Call the static CHEBTECH.CHEBPTS() method:
x = f.chebpts(n);
% Scale the domain:
x = scaleNodes(x, dom);
if ( nargout > 1 )
% Call the static CHEBTECH.CHEBPTS() method:
w = f.quadwts(n);
% Scale the domain:
w = scaleWeights(w, dom);
end
if ( nargout > 2 )
v = f.barywts(n);
end
else % Piecewise grid.
% Initialise cell arrays for temporary storage:
x = cell(length(n), 1);
w = cell(1, length(n));
v = cell(length(n), 1);
% Loop over the number of subintervals: (as above)
for k = 1:numel(n)
x{k} = f.chebpts(n(k));
x{k} = scaleNodes(x{k}, dom(k:k+1));
if ( nargout > 1 )
w{k} = f.quadwts(n(k));
w{k} = scaleWeights(w{k}, dom(k:k+1));
end
if ( nargout > 2 )
v{k} = f.barywts(n(k));
end
end
% Convert the cell to an array for output:
x = cell2mat(x);
w = cell2mat(w);
v = cell2mat(v);
end
end
function y = scaleNodes(x, dom)
%SCALENODES Scale the Chebyshev nodes X from [-1,1] to DOM.
% TODO: Deal with unbounded domains
if ( dom(1) == -1 && dom(2) == 1 )
% Nodes are already on [-1, 1];
y = x;
return
end
% Scale the nodes:
y = dom(2)*(x + 1)/2 + dom(1)*(1 - x)/2;
end
function w = scaleWeights(w, dom)
%SCALEWEIGHTS Scale the Chebyshev weights W from [-1,1] to DOM.
% TODO: Deal with unbounded domains
if ( dom(1) == -1 && dom(2) == 1 )
% Nodes are already on [-1, 1];
return
end
% Scale the weights:
w = (diff(dom)/2)*w;
end