Trapezoidal numerical integration without use of function?

365 ビュー (過去 30 日間)
Weronika Ring
Weronika Ring 2017 年 11 月 3 日
コメント済み: Mpho 2023 年 4 月 11 日
I was trying to do trapezoid integration without actually creating a function, but I seem to have done something wrong.
a=0;
b=1;
f=@(x)x.*sin(x);
n=100;
x=linspace(a,b,n+1);
h=(b-a)/n;
qt=sum((h*f((x(1:n)+x(2:n+1)))/2))
This returns 0.4354, when the actual value of the integral is 0.3012. Where am I going wrong with the code?
  2 件のコメント
murat onay
murat onay 2019 年 5 月 21 日
編集済み: murat onay 2019 年 5 月 21 日
trapezoid matlab code is here. This code return 0.301180
clc;clear;
a=0;
b=1;
n=100;
h=(b-a)/n;
sum=0;
f=@(x) x.*sin(x);
for i=1:1:n-1
sum= sum + f(a+i*h);
end
result = h/2*(f(a)+f(b)+2*sum);
fprintf('%f',result);
Mpho
Mpho 2023 年 4 月 11 日
Hi, how would I achieve this if I want to integrate a constant matrix?

サインインしてコメントする。

回答 (3 件)

John D'Errico
John D'Errico 2017 年 11 月 3 日
編集済み: John D'Errico 2017 年 11 月 7 日
Trapezoidal rule is easy enough. It depends on whether the step is constant or not. The entire point of my response is you need to get the weights correct. If not, then of course your code must fail.
For a constant step size, you need to remember that the weights look like this:
h*[1 2 2 2 2 2 ... 1]/2
So in any interval, we have a trapezoid. The area of a trapezoid is the average of the heights at each end, then multiply by the width.
h*(f(x(i)) + f(x(i+1)))/2
But each trapezoid shares its endpoint with the neighbors. So that gives the first and last points half the weight of the rest. It also gives us a simple way to do trapezoidal rule.
n = 20;
x = linspace(0,pi,n);
dx = x(2) - x(1);
f_x = sin(x);
trapint = (sum(f_x) - (f_x(1) + f_x(end))/2)*dx
trapint =
1.9954
Did we get it correct?
trapz(x,f_x)
ans =
1.9954
Yes. And both agree with the actual result.
syms u
int(sin(u),0,pi)
ans =
2
Could I have done this for unequal spacing? Still easy enough. I'll pick a random spacing here, using a few more points because a random spacing can have some wide places where nothing fell.
n = 50;
x = [0,sort(rand(1,n-2))*pi,pi];
f_x = sin(x);
dx = diff(x);
trapint = dot(dx,(f_x(1:end-1) + f_x(2:end))/2)
trapint =
1.9966
Note that on a function like sin(x) over that interval, trapezoidal rule will tend to underestimate the integral. As you can see, this is exactly what happened, and will always happen for that function, on that interval.

Hiren Rana
Hiren Rana 2021 年 11 月 11 日
a=0; b=1; n=100; h=(b-a)/n; sum=0; f=@(x) x.*sin(x); for i=1:1:n-1 sum= sum + f(a+i*h);
end result = h/2*(f(a)+f(b)+2*sum); fprintf('%f',result);

VBBV
VBBV 2021 年 12 月 29 日
編集済み: VBBV 2022 年 2 月 5 日
a=0;
b=1;
f=@(x)x.*sin(x);
n=100;
x=linspace(a,b,n+1);
h=(b-a)/n;
qt=sum(h*(f(x(1:n))+f(x(2:n+1)))/2)
qt = 0.3012
You did not use the function correctly. Try above

カテゴリ

Help Center および File ExchangeCalculus についてさらに検索

Community Treasure Hunt

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

Start Hunting!

Translated by