A.3. Метод Трапеций
Метод
Трапеций применяется для вычисления
определенного интеграла
у которого подынтегральная функция
линейная или содержит элементарные
функции.
Разобьем интервал интегрирования [a;b] на n равных отрезков и восстановим из полученных точек a, , , …, b перпендикуляры до пересечения с графиком функции. Соединив последовательно точки пересечения, представим площадь полученной криволинейной трапеции как сумму прямолинейных трапеций, площади которых легко подсчитать. Заменив подынтегральную функцию f(x) в пределах элементарного отрезка [ ; ] рисунок 6.
Рисунок 6 – Иллюстрация метода Трапеций
следующие формулы для элементарных площадей:
Тогда общая площадь равна:
I=
Отсюда получаем формулу трапеций:
Для
достижения заданной точности, при
вычислении определенного интеграла
используется метод двойного просчета.
Он основан на двукратном вычислении
значения интеграла вначале с шагом h
, а затем с шагом h/2.
Полученные значения интегралов
и
могут быть применены для оценки
погрешности интегрирования по формуле:
где k = 2 для формулы Трапеций.
Если полученная погрешность не удовлетворяет требуемой точности, то вычисляется значение интеграла при h=h/4 и снова оценивается погрешность, и т.д. до тех пор, пока не окажется, что погрешность стала меньше заданной точности.
Приложение b. Тестирование функций, реализующих выбранные численные методы
Рисунок 7 – Форма с данными тест
B.1.Метод Эйлера с автоматическим выбором шага
Для
тестирования функций, реализующих
Эйлера с автоматическим выбором шага,
найдем решение дифференциального
уравнения
на отрезке [a;b]
с табличным шагом
.
Исходные данные:
Аналитическое решение ОДУ.
+
x + C
y=
y=
Из
начальных условий
получим
.
c=1
Получаем аналитическое решение уравнения
Табл.
4
Значение аргумента x |
Аналитическое решение |
Решения, полученные программой |
Решения, полученные Scilab |
0 |
1 |
1 |
1 |
0.05 |
1.05259 |
1.0519187 |
1.0525860 |
0.1 |
1.11071 |
1.1092670 |
1.1107106 |
0.15 |
1.17498 |
1.1726304 |
1.1749786 |
0.2 |
1.24608 |
1.2426722 |
1.2460767 |
0.25 |
1.32478 |
1.3201456 |
1.3247847 |
0.3 |
1.41199 |
1.4064718 |
1.4119899 |
0.35 |
1.5087 |
1.5021636 |
1.5087024 |
0.4 |
1.61607 |
1.6083484 |
1.6160743 |
0.45 |
1.73542 |
1.7263123 |
1.7354209 |
0.5 |
1.86825 |
1.8576247 |
1.8682459 |
Исходя из результатов, представленных на таблице 4, можно сделать вывод, что процедура работают правильно.
Форма с результатами, полученными при помощи тестируемых функций, представлена на рисунке 8.
Рисунок 8 – Форма с решением ОДУ, тест
Решение тестового дифференциального уравнения при помощи математического пакета Scilab 6.0.2 представлено ниже:
--> //Вычисление решения ОДУ на отрезке [0;0.5]с шагом 0.05
--> deff('yd=f(y,x)','yd=y*(x+1)')
--> //Задание начальных условий
--> y0=1;x0=0;x=0:0.05:0.5;
--> //Описание функции f(x)
--> y=ode(y0,x0,x,f);
--> yr=[x;y];//Создание и вывод матрицы решений
--> yr'
ans =
0. 1.
0.05 1.0525860
0.1 1.1107106
0.15 1.1749786
0.2 1.2460767
0.25 1.3247847
0.3 1.4119899
0.35 1.5087024
0.4 1.6160743
0.45 1.7354209
0.5 1.8682459
Схема алгоритма функции, реализующей метод автоматического выбора шага, представлена на рисунках с 9 по 11.
Рисунок 9 – Схема алгоритма процедуры Runge_step
Рисунок 10 – Схема алгоритма функции r
Рисунок 11 – Схема алгоритма функции test_f_odu
Код тестируемого модуля.
Option Strict On
Option Explicit On
Imports System.Math
Module Module2
'Процедуры, реализующей метод автоматического выбора шага
Sub Runge_step(ByVal x0 As Double, ByVal y0 As Double, ByVal h0 As Double, ByVal n As Integer, ByVal eps As Double, ByRef xm() As Double, ByRef ym() As Double)
Dim h, y, y1, x As Double
Dim m As Integer 'коэффициент дробления шага
Dim p As Integer
Dim i As Integer = 0
Form3.ListBox1.Items.Add(i)
Form3.ListBox4.Items.Add(vbNewLine)
Form3.ListBox5.Items.Add(vbNewLine)
p = 1 'метод Эйлер - это формула Рунге-Кутта 1 порядка отсюда следует ,что p=1
xm(0) = x0
ym(0) = y0
For i = 1 To n
h = h0
m = 1
y = r(x0, y0, h, m)
Do
y1 = y : h = h / 2
x = x0 : y = y0
m = 2 * m
y = r(x, y, h, m)
Loop Until Abs(y - y1) / (2 ^ p - 1) < eps
x0 = x0 + h0
y0 = y
xm(i) = x0
ym(i) = y0
Form3.ListBox1.Items.Add(i)
Form3.ListBox4.Items.Add(h)
Form3.ListBox5.Items.Add(m)
Next i
printLBx(xm, Form2.ListBox1)
printLBy(ym, Form2.ListBox2)
printLBx(xm, Form3.ListBox2)
printLBy(ym, Form3.ListBox3)
End Sub
'Процедура решения ОДУ методом Эйлера
Function r(ByVal x As Double, ByVal y As Double, ByVal h As Double, ByVal m As Double) As Double
For j = 1 To m
y = y + h * test_f_odu(x, y)
x = x + h
Next j
Return y
End Function
'Функция для теста
Function test_f_odu(ByVal x As Double, ByVal y As Double) As Double
test_f_odu = y * (x + 1)
End Function
End Module