polyval error with recurrence relation

1 次查看（过去 30 天）
T
T 2013-10-6
评论： T 2013-10-11
function c = relation( n )
if n==0
c = 1;
elseif n==1
syms x
c = [x 0];
else
syms x
c = ((2*n-1)*x*[relation(n-1),0] - (n-1)*[0,0,relation(n-2)])/n;
sum(c)
end
n=8
m=n+1;
xval=linspace(-1,1,5000)
I=trapsum(polyval(relation(n),xval).*polyval(relation(m),xval),-1,1);
Gives the following output:
??? Error using ==> polyval
Inputs to polyval must be floats, namely single or double.
Error in ==> polyval at 63
y = zeros(siz_x, superiorfloat(x,p));
What needs to be done to ensure that the output for the recurrence relation is a double?
I mean, in the function, we have syms x but this is to acquire the polynomial in analytic form. It's not necessary in this case but I need to plot it later on. The main issue I have is evaluating the integral using trapsum but I would need to do both; plot for various n, evaluate the integral.

请先登录，再回答此问题。

采纳的回答

Walter Roberson
Walter Roberson 2013-10-6
In your code for relation(), why do you have
sum(c)
which is just displaying the sum rather than returning the sum ?
After you build the relation symbolically, use
RelationM = relation(m);
RelationN = relation(n);
funM = matlabFunction(RelationM);
funN = matlabFunction(RelationN);
Y = funM(xval) .* funN(xval);
trapsum(Y, -1, 1)
However, I do not know what your trapsum routine does: it is not a routine in any Mathworks toolbox, and the obvious trapsum in the File Exchange uses a different calling sequence.
  14 个评论
显示 12更早的评论
T
T 2013-10-8
That's right. And for some reason using c1 and c2 will not generate the same polynomials for n where if you used
c = ((2*n-1)*x*[relation(n-1),0] - (n-1)*[0,0,relation(n-2)])/n;
The point was to do this using recursive programming.
T
T 2013-10-11
编辑：T 2013-10-11
Actually I don't think I need to do any of the above:
for n=0:inf
m=n+1;
int( ((2*n-1)*x*[relation(n-1),0] - (n-1)*[0,0,relation(n-2)])/n )* ((2*m-1)*x*[relation(m-1),0] - (m-1)*[0,0,relation(m-2)])/m),x,-1,1)
end
This will just give me zero. But is it possible to determine the round-off error? Or at least get the maximum value of n before the loop crashes?
I can set the recursion limit using:
set(0, 'RecursionLimit', 100000) but it keeps crashing still.

请先登录，再进行评论。

更多回答（1 个）

Matt J
Matt J 2013-10-6
编辑：Matt J 2013-10-6
Get rid of the "syms x". You are doing only numerical manipulations with the polynomial so there is no apparent need to have it in symbolic form, and that includes for plotting. Example
t=linspace(-1,1,1000);
p=[1 2 1]; %polynomial coefficients
plot(t,polyval(p,t))
Also, why use trapsum instead of the builtin trapz?
  3 个评论
显示 1更早的评论
Matt J
Matt J 2013-10-6
It's the function relation() that can't be found. Check your path.
T
T 2013-10-6
Okay now I receive the error:
??? Undefined function or variable 'x'.
Error in ==> rlegendre at 9
c = ((2*n-1)*x*[rlegendre(n-1),0] - (n-1)*[0,0,rlegendre(n-2)])/n;
This happens when I remove syms x.

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Common Operations 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by