square matrix containing Chebyshev polynomials becomes singular when its size becomes greater than a certain value
이전 댓글 표시
Hi, I've been trying to approximate functions with Chebyshev polynomials. In this case, I'm looking at 
, and I've approximated 
as 
where 
are coefficients and 
is the k-th Chebyshev polynomial. I've got N+1 points 
inside the interval [0,4] such that 
. My goal is to find what the coefficients are, and I did this by creating a matrix equation of the form 
and solving for x, where x is a vector of length N+1 containing the coefficients 
. So I applied the previous Chebyshev approximation on all N+1 points and got a matrix equation where A = 
, x = 
and B= 
, and I solved for x. Now, in order to simplify things, I chose my points to be evenly spaced within [0,4]. When trying out my code, I would get a certain vector for x and compare it via chebcoeffs with the actual coefficients. My estimated coefficients from x got increasingly closer to those obtained from chebcoeffs as I increased N starting from 1, as one would expect. However, when I hit N = 56 and onwards, my matrix A began to not be full rank, and hence be singular. I'm not sure why this happens, even though my definition of all matrices and vectors appears correct from a theoretical standpoint. Here's the full code:
, and I've approximated as where are coefficients and is the k-th Chebyshev polynomial. I've got N+1 points inside the interval [0,4] such that . My goal is to find what the coefficients are, and I did this by creating a matrix equation of the form and solving for x, where x is a vector of length N+1 containing the coefficients . So I applied the previous Chebyshev approximation on all N+1 points and got a matrix equation where A = , x = and B= , and I solved for x. Now, in order to simplify things, I chose my points to be evenly spaced within [0,4]. When trying out my code, I would get a certain vector for x and compare it via chebcoeffs with the actual coefficients. My estimated coefficients from x got increasingly closer to those obtained from chebcoeffs as I increased N starting from 1, as one would expect. However, when I hit N = 56 and onwards, my matrix A began to not be full rank, and hence be singular. I'm not sure why this happens, even though my definition of all matrices and vectors appears correct from a theoretical standpoint. Here's the full code:t0 = 0;
tf = 4;
N = 56;
t = chebfun('t',[t0 tf]);
A = zeros(N+1,N+1);
d = zeros(N+1,1);
for c=0:N
d(c+1,1)=tf*c./(N);
end
h = chebpoly(0:N,[t0 tf]);
for r=0:N
A(r+1,:)= h(d(r+1,1));
end
A
B = zeros(N+1,1)
B(:,1) = (d).^3-(d)+(d).^2;
B
X= A\B
b=chebcoeffs((t).^3-(t)+(t).^2) %for comparison with actual chebyshev coefficients
rank(A)
댓글 수: 6
Torsten
2023년 11월 27일
Why do you use the Chebyshev polynomials for approximation as if they were the usual polynomial basis 1,x,x^2,...,x^n ? Chebyshev polynomials are orthogonal with respect to a special scalar product, and this property is used to determine the coefficients a_n in your above expansion.
Torsten
2023년 11월 28일
If you use the MATLAB function "chebyshevT" to create the polynomials and convert them to numerical functions using "matlabFunction", your code also works for n > 56.
Walter Roberson
2023년 11월 28일
Chang-Yu
2023년 11월 28일
Chang-Yu
2023년 11월 28일
채택된 답변
t0 = 0;
tf = 4;
N = 75;
syms x
h = chebyshevT(0:N, x);
A = zeros(N+1,N+1);
d = zeros(N+1,1);
for c=1:N+1
d(c,1) = 0.5*(t0+tf)+0.5*(tf-t0)*cos((2*c-1)/(2*(N+1))*pi);
end
for c=1:N+1
A(c,:)= vpa(subs(h,x,cos((2*c-1)/(2*(N+1))*pi)));
end
A
B = zeros(N+1,1);
B(:,1) = (d).^3-(d)+(d).^2;
X= A\B
hold on
plot(d,(d).^3-(d)+(d).^2)
fplot(subs(h,x,(x-0.5*(t0+tf))/(0.5*(tf-t0)))*X,[t0 tf])
hold off
grid on
댓글 수: 2
That is exactky what I have suggested in my answer: use the Chebyschev nodes
t0 = 0;
tf = 4;
N = 75;
syms x
h = chebyshevT(0:N, x);
A = zeros(N+1,N+1);
for c=1:N+1
A(c,:)= vpa(subs(h,x,cos((2*c-1)/(2*(N+1))*pi)));
end
and the matrix is absolutely well conditioned for any (large) N
cond(A)
In comparison to equidistant points:
t0 = 0;
tf = 4;
N = 75;
syms x
h = chebyshevT(0:N, x);
A = zeros(N+1,N+1);
d = zeros(N+1,1);
for c=0:N
d(c+1,1) = 0.5*(t0+tf)+0.5*(tf-t0)*(-1+2*c/N);
end
for c=0:N
A(c+1,:)= vpa(subs(h,x,-1+2*c/N));
end
A
B = zeros(N+1,1);
B(:,1) = (d).^3-(d)+(d).^2;
X= A\B
hold on
plot(d,(d).^3-(d)+(d).^2)
fplot(subs(h,x,(x-0.5*(t0+tf))/(0.5*(tf-t0)))*X,[t0 tf])
hold off
grid on
추가 답변 (1개)
Bruno Luong
2023년 11월 27일
0 개 추천
If you select discrete points d the Chebychev nodes see definition wiki, and not uniform, it will become well posed.
댓글 수: 2
Chang-Yu
2023년 11월 27일
Bruno Luong
2023년 11월 27일
Try it, if it works then I give you a mathematic explanation.
카테고리
도움말 센터 및 File Exchange에서 Polynomials에 대해 자세히 알아보기
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
