I want to extract w and x from this code, but unfortunately i cant and it gives me only 1 answer.

2 次查看（过去 30 天）
Alireza
Alireza 2022-4-11
回答： Riccardo Scorretti 2022-4-11
function [sol]=numintla(fun,n)
[xx,w]=gaussla(6); % Call the gaussla.m file
sx=size(xx,1);
sol = 0;
for i=1:sx,
x=xx(i);
fx=eval(fun);
sol=sol+w(i)*fx;
end
end
function [xx,ww]=gaussla(n)
format long;
b(1)=1;
for i=1:n, % Calculation of the binomial coefficients
b(i+1)=(factorial(n))/(factorial(i)*(factorial(n+1-i)));
end
for i=1:n+1, % The polynomial coefficients
poly(i)=((-1)^(n+1-i))*b(i)*(factorial(n)/(factorial(n+1-i)));
end
xx=roots(poly);% The polynomial roots
for i=1:n, % Coefficients of the first derivative of the polynomial
polycd(i)=poly(i)*(n+1-i);
end
for i=1:n, % Evaluation
x=xx(i);
solde=0;
for k=1:n,
solde=solde+polycd(k)*(x^(n-k));
end
ww(i,1)=(factorial(n)^2)/(xx(i)*(solde^2));
end
end

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

回答（2 个）

Walter Roberson
Walter Roberson 2022-4-11
function [sol, w, x] = numintla(fun,n)
I would recommend that you reconsider your requirement. Your x variable is whatever was last written into x by your for i loop. It would make more sense if you were to return xx instead of x.
Riccardo Scorretti
Riccardo Scorretti 2022-4-11
Hi Alireza,
I'm afraid that the fact your code returns only one value is the last of your problems. More specifically, you provides two functions: one returns two values and the other returns two. So, it depends on what you want to do: if your purpose is to get the quadrature formula, you can invoke guass1a in this way (the number 3 is just an example: you can use any positive integer):
[x, w] = gaussla(3)
In my opinion (but I may be wrong) the true problems are more fundamental. You seems to be willing to compute (gaussla) the nodes and the weights of a Gauss-Legendre quadrature formula, which is usually defined in the range [-1 1]. The function numintla integrates numerically the function fun over a not precised domain (which I suppose is [-1 1]).
I'm afraid the function gaussla returns the wrong result. For instance, you can check the first coefficients here: https://en.wikipedia.org/wiki/Gaussian_quadrature. For n = 1 and n = 2 you get respectively:
x =
1
w =
1
and
x =
1.000000000000000 + 1.000000000000000i
1.000000000000000 - 1.000000000000000i
w =
-0.500000000000000 + 0.500000000000001i
-0.500000000000000 - 0.500000000000001i
One observes that coefficients are complex (this is wrong). The error is likely ot be in the computation of the coefficients of Legendre polynomials (see https://en.wikipedia.org/wiki/Legendre_polynomials).
If I'm not wrong, I suggest you to have a look this code: https://fr.mathworks.com/matlabcentral/fileexchange/4540-legendre-gauss-quadrature-weights-and-nodes

类别

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

标签

产品

Community Treasure Hunt

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

Start Hunting!

Translated by