ВУЗ: Не указан
Категория: Не указан
Дисциплина: Не указана
Добавлен: 06.04.2021
Просмотров: 156
Скачиваний: 1
# i, j — переменные циклов
# xp — середины частичных отрезков
# n — количество частичных отрезков
# x — массив узлов
# h — длины частичных отрезков
n:=1; x[1]:=a;
for i from 1 while (x[i]<=b) do
n:=n+1;
x[i+1]:=x[1]+sh*i; # заполнение массива узлов
end do;
if x[n]<b then x[n]:=b; end if;
for i from 1 to n-1 do
# вычисление длин частичных отрезков
h[i]:=abs(x[i+1]-x[i]);
# нахождение середин частичных отрезков
xp[i]:=h[i]/2+x[i];
end do;
slog:=0;
for j from 1 to n-1 do
slog:=slog+simplify(h[j]*f(xp[j])); # формула прямоугольников
end do;
slog:=evalf(slog);
return slog; # вывод результата
end proc;
Проверим результат работы процедуры
> prl(0,1,0.0001);
0.5222818376
> evalf(int((sin(x))-xˆ 3*ln(x),x=0..1));
0.5221976941
4.9.2.
Метод трапеций
> restart;
> f:=x-> evalf(sin(x)-xˆ 3*ln(x));
# вычисление определенного интеграла методом трапеций
> trp:=proc(a,b,sh)
11
# b, a — верхний и нижний пределы интегрирования соответственно
# sh — шаг разбиения отрезка [a,b]
local i,j,n,x,y,h,slog;
# i, j — переменные циклов
# n — количество частичных отрезков
# x — массив узлов
# y — массив значений функции в полученных узлах
# h — длины частичных отрезков
n:=1; x[1]:=a; y[1]:=f(x[1]);
for i from 1 while (x[i]<=b) do
n:=n+1;
x[i+1]:=x[1]+sh*i; # заполнение массива узлов
# заполнение массива значений функции в полученных узлах
y[i+1]:=f(x[i+1]);
end do;
if x[n]<b then x[n]:=b; y[n]:=f(x[n]); end if;
for i from 1 to n-1 do
h[i]:=abs(x[i+1]-x[i]); # вычисление длин частичных отрезков
end do;
slog:=0;
for j from 1 to n-1 do
slog:=slog+evalf(h[j]*(y[j]+y[j+1])); # формула трапеций
end do;
slog:=evalf(1/2*slog);
return slog; # вывод результата
end proc;
Проверим работу процедуры
> trp(0.001,1,0.001);
0.5230383120
> evalf(int((sin(x))-xˆ 3*ln(x),x=0..1));
0.5221976941
4.9.3.
Метод Симпсона
> restart;
> f:=x-> evalf(sin(x)+ln(x)*xˆ 3);
12
> Simpson:=proc(a,b,n)
# b, a — верхний и нижний пределы интегрирования соответственно
# n — количество частичных отрезков
local i,h,x,y,slog,slogn,slogch;
# i — переменная циклов
# h — шаг
# x — массив узлов
# y — массив значений функции в полученных узлах
# slogn — сумма всех значений функции в узлах с нечетным номером
# slogch — сумма всех значений функции в узлах с четным номером
(кроме x[0] и x[n])
if (n mod 2 = 1) then print(’n-нечетное_число’); exit; end if;
h:=abs(evalf(b)-evalf(a))/n; # определение шага, для заполнения
массива узлов
for i from 0 to n do
x[i]:=evalf(a)+i*h; # заполнение массива узлов
# заполнение массива значений функции в полученных узлах
y[i]:=f(x[i]);
end do;
slog:=0; slogn:=0; slogch:=0;
for i from 1 by 2 to (n-1) do
# вычисление суммы всех значений функции в узлах с нечетным
номером
slogn:=slogn+y[i];
end do;
for i from 2 by 2 to (n-2) do
# вычисление суммы всех значений функции в узлах с четным номе-
ром (кроме x[0] и x[n])
slogch:=slogch+y[i];
end do;
slog:=(y[0]+4*slogn+2*slogch+y[n])*h/3; # формула Симпсона
return slog; # вывод результата
end proc;
Проверим работу процедуры
> Simpson(0.0001,Pi,22);
23.78877518
> int(f(x),x=0.0001..Pi);
23.78870622
13
4.9.4.
Метод Гаусса
Таблица узлов и коэффициентов формулы Гаусса для n=1..6
n
x
i
A
i
1
x
1
= 0
A
1
= 2
2
x
2
=
−
x
1
= 0
.
5773502692
A
1
=
A
2
= 1
3
x
2
= 0
A
2
= 0
.
8888888889
x
3
=
−
x
1
= 0
.
7745966692
A
1
=
A
3
= 0
.
5555555556
4
x
3
=
−
x
2
= 0
.
3399810436
A
2
=
A
3
= 0
.
6521451549
x
4
=
−
x
1
= 0
.
8611363116
A
1
=
A
4
= 0
.
3478548451
5
x
3
= 0
A
3
= 0
.
5688888899
x
4
=
−
x
2
= 0
.
53384693101
A
2
=
A
4
= 0
.
4786286705
x
5
=
−
x
1
= 0
.
9061798459
A
1
=
A
5
= 0
.
2369268851
6
x
6
=
−
x
1
= 0
.
9324695142
A
1
=
A
6
= 0
.
1713244924
x
5
=
−
x
2
= 0
.
6612093865
A
2
=
A
5
= 0
.
3607615730
x
4
=
−
x
3
= 0
.
2386191861
A
3
=
A
4
= 0
.
4679139346
> restart;
> n:=6:
# табличные значения узлов формулы Гаусса для n=6
> x[1]:=-0.9324695142:
> x[2]:=-0.6612093865:
> x[3]:=-0.2386191861:
> x[4]:=-x[3]:
> x[5]:=-x[2]:
> x[6]:=-x[1]:
# табличные значения коэффициентов формулы Гаусса для n=6
> A[1]:=0.1713244924:
> A[2]:=0.3607615730:
> A[3]:=0.4679139346:
> A[4]:=A[3]:
> A[5]:=A[2]:
> A[6]:=A[1]:
# процедура для вычисления интеграла методом Гаусса с 6 узлами
> Gauss6:=proc(a,b)
# b, a — верхний и нижний пределы интегрирования
local i,sl,t;
# i — переменная циклов
14
# sl — переменная для накопления суммы слагаемых
# t — массив узлов
for i from 1 to n do
t[i]:=(b+a)/2+((b-a)/2)*x[i]; # заполнение массива узлов
end do;
sl:=0;
for i from 1 to n do
sl:=sl+evalf(A[i]*f(t[i])); # формула Гаусса
end do;
sl:=evalf((b-a)/2*sl); # вывод результата
end proc;
# процедура для вычисления определенного интеграла (метод Гаус-
са) с заданной точностью eps
> Gauss:=proc(a,b,eps)
# b, a — верхний и нижний пределы интегрирования
# eps — точность интегрирования
local int1,int2,int,i,j,c,k,r,l;
# int1, int2 — два последовательных приближения к точному реше-
нию интеграла
# int — приближения к точному решению интеграла
# i, j — переменные циклов
# c — точки разбиения промежутка интегрирования
# k — количество частичных отрезков
# r, l — правые и левые границы частичных отрезков
# вычисление определенного интеграла (метод Гаусса n=6)
int1:=Gauss6(a,b);
c:=(b-a)/2; # точка разбиения промежутка интегрирования
# вычисление определенного интеграла (метод Гаусса n=11)
int2:=Gauss6(a,c)+Gauss6(c,b);
# количество частичных отрезков, на которое рабивается отрезок
[a,b]
k:=2;
for j from 1 while abs(int1-int2)>eps do
int:=0;
c:=(b-a)/(2*k);
for i from 0 to 2*k-1 do
# вычисление границ каждого частичного отрезка
l:=a+c*i; r:=a+c*(i+1);
# вычисление определенного интеграла (n=6+2*k-1)
15