Matlab Code for the Gauss Legendre Quadrature
564 次查看(过去 30 天)
显示 更早的评论
I need help contsructing the code for the gauss legendre quadrature using n =2 , 4, and 6. i was able to get the value for n =2 but after that im kind of lost.
%% parameters
a = -1; % lower bound
b = 1; % upper bound
n = [-1/sqrt(3) 1/sqrt(3)]; %location values for n=2
n1 = [-0.86113 -0.33998 0.33998 0.86113]; %location values for n=4
n2 = [-0.93246 -0.66120 -0.23861 0.23861 0.66120 0.993246]; % location values for n=6
w = 1; %weight for n=2
w1 = [0.34785 0.625214]; %weights for n=4
w2 = [0.17132 0.36076 0.6791]; % weights for n=6
f = @(x) (1)/(x.^2 * (sqrt(x.^2 +1)));%function
%% Gauss_Legendre Quadrature
% n=2
for i = 1:length(n)
h = n(i);
y(i) = w*f(h);
end
gaussleg1 = sum(y); % sum of wi*fi for n=2
%n=4
%n=6
0 个评论
采纳的回答
Tommy
2020-4-16
编辑:Tommy
2020-4-16
Have your vectors of weights be the same size as your vectors of locations, ordered so that corresponding weights and locations are in the same position:
%% parameters
a = -1; % lower bound
b = 1; % upper bound
n = [-1/sqrt(3) 1/sqrt(3)]; %location values for n=2
n1 = [-0.86113 -0.33998 0.33998 0.86113]; %location values for n=4
n2 = [-0.93246 -0.66120 -0.23861 0.23861 0.66120 0.993246]; % location values for n=6
w = [1 1]; %weight for n=2
w1 = [0.34785 0.625214 0.625214 0.34785]; %weights for n=4
w2 = [0.17132 0.36076 0.6791 0.6791 0.36076 0.17132]; % weights for n=6
f = @(x) (1)./(x.^2 .* (sqrt(x.^2 +1)));%function
%% Gauss_Legendre Quadrature
% n=2
gaussleg = sum(w.*f(n));
%n=4
gaussleg1 = sum(w1.*f(n1));
%n=6
gaussleg2 = sum(w2.*f(n2));
Higher precision in your weights and locations will give more accurate results.
(edit) I also changed f slightly so that you could pass in the entire vector of locations at once.
4 个评论
John D'Errico
2022-3-28
Almost certainly you have created a variable with the name sum. DON'T DO THAT!!!!!
Look at the list of variables in your workspace. Think about what happens now when MATLAB tries to use the function sum. It gets confused, and then you get that equally confusing error message.
更多回答(2 个)
Manikanta Balaji
2022-6-9
Y=input('Enter function directly f(x) =','s');
f=inline(Y);
a=input('Enter Initial interval point ''a'' =');
b=input('Enter final interval point ''b'' =');
n=input('Enter Gauss point ''n'' =');
if (n>=2 && n<=9)
%% data
%abscissa values
X=[-0.57735 -0.7746 -0.86114 -0.90618 -0.932470 -0.949108 -0.96029 -0.96816
0.57735 0.0000 -0.33998 -0.53847 -0.661209 -0.741531 -0.79666 -0.83603
0.00000 0.7746 0.33998 0.00000 -0.238619 -0.405845 -0.52553 -0.61337
0.00000 0.0000 0.86114 0.53847 0.238619 0.000000 -0.18343 -0.32425
0.00000 0.0000 0.00000 0.90618 0.661209 0.405845 0.18343 0.00000
0.00000 0.0000 0.00000 0.00000 0.932470 0.741531 0.52553 0.32425
0.00000 0.0000 0.00000 0.00000 0.000000 0.949108 0.79666 0.61337
0.00000 0.0000 0.00000 0.00000 0.000000 0.000000 0.96029 0.83603
0.00000 0.0000 0.00000 0.00000 0.000000 0.000000 0.00000 0.96816];
%weight values
w=[1.00000 0.5555 0.34785 0.23693 0.171324 0.129484 0.10123 0.08127
1.00000 0.8888 0.65214 0.47863 0.360761 0.279705 0.22238 0.18064
0.00000 0.5555 0.65214 0.56889 0.467914 0.381830 0.31370 0.26061
0.00000 0.0000 0.34785 0.47863 0.467914 0.417959 0.36268 0.31234
0.00000 0.0000 0.00000 0.23693 0.360761 0.381830 0.36268 0.33024
0.00000 0.0000 0.00000 0.00000 0.171324 0.279705 0.31370 0.31234
0.00000 0.0000 0.00000 0.00000 0.000000 0.129484 0.22238 0.26061
0.00000 0.0000 0.00000 0.00000 0.000000 0.000000 0.10123 0.18064
0.00000 0.0000 0.00000 0.00000 0.000000 0.000000 0.00000 0.08127];
%% Gauss Quadrature
G=0;
F=@(t) f(((b-a)*t+(b+a))/2);
for i=1:n
G=G+w(i,n-1)*F(X(i,n-1))*(b-a)/2;
end
fprintf('Integration by %d Gauss Quadrature = %d',n,G);
else
fprintf('Formula not available for %d point',n)
end
0 个评论
Mushimiyimana
2024-1-15
Y=input('Enter function directly f(x) =','s');
f=inline(Y);
a=input('Enter Initial interval point ''a'' =');
b=input('Enter final interval point ''b'' =');
n=input('Enter Gauss point ''n'' =');
if (n>=2 && n<=9)
%% data
%abscissa values
X=[-0.57735 -0.7746 -0.86114 -0.90618 -0.932470 -0.949108 -0.96029 -0.96816
0.57735 0.0000 -0.33998 -0.53847 -0.661209 -0.741531 -0.79666 -0.83603
0.00000 0.7746 0.33998 0.00000 -0.238619 -0.405845 -0.52553 -0.61337
0.00000 0.0000 0.86114 0.53847 0.238619 0.000000 -0.18343 -0.32425
0.00000 0.0000 0.00000 0.90618 0.661209 0.405845 0.18343 0.00000
0.00000 0.0000 0.00000 0.00000 0.932470 0.741531 0.52553 0.32425
0.00000 0.0000 0.00000 0.00000 0.000000 0.949108 0.79666 0.61337
0.00000 0.0000 0.00000 0.00000 0.000000 0.000000 0.96029 0.83603
0.00000 0.0000 0.00000 0.00000 0.000000 0.000000 0.00000 0.96816];
%weight values
w=[1.00000 0.5555 0.34785 0.23693 0.171324 0.129484 0.10123 0.08127
1.00000 0.8888 0.65214 0.47863 0.360761 0.279705 0.22238 0.18064
0.00000 0.5555 0.65214 0.56889 0.467914 0.381830 0.31370 0.26061
0.00000 0.0000 0.34785 0.47863 0.467914 0.417959 0.36268 0.31234
0.00000 0.0000 0.00000 0.23693 0.360761 0.381830 0.36268 0.33024
0.00000 0.0000 0.00000 0.00000 0.171324 0.279705 0.31370 0.31234
0.00000 0.0000 0.00000 0.00000 0.000000 0.129484 0.22238 0.26061
0.00000 0.0000 0.00000 0.00000 0.000000 0.000000 0.10123 0.18064
0.00000 0.0000 0.00000 0.00000 0.000000 0.000000 0.00000 0.08127];
%% Gauss Quadrature
G=0;
F=@(t) f(((b-a)*t+(b+a))/2);
for i=1:n
G=G+w(i,n-1)*F(X(i,n-1))*(b-a)/2;
end
fprintf('Integration by %d Gauss Quadrature = %d',n,G);
else
fprintf('Formula not available for %d point',n)
end
0 个评论
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Numerical Integration and Differential Equations 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!