- •4.Численные методы решения задач
- •Общая методика метода Гаусса. Система представлена в общем виде как:
- •4.2.Решение нелинейных уравнений
- •4.3. Методы построения аппроксимационных полиномов
- •Определить аппроксимационный многочлен
- •4.4.Численное дифференцирование и интегрирование
- •Ошибка вычисления интеграла определяется как
- •4.5.Методы решения дифференциальных уравнений. Метод Рунге-Кутта
- •4.6.Методы решения дифференциальных уравнений с частными производными
4.4.Численное дифференцирование и интегрирование
К численному дифференцированию приходиться прибегать в том случае, когда функция f(x), для которой нужно найти производную, задана таблично или же функциональная зависимость f(x) имеет сложное аналитическое выражение. В этих случаях функция разбивается одномсерной сеткой и используются приближенные формулы. Формула первой производной для двух узлов:
(4.4.1), , h-шаг изменения аргумента.
Формула первой производной по трем узлам для второй точки:
(4.4.2) , где
Формула второй производной для трех узлов:
(4.4.3)
Численное интегрирование применяют в случаях когда нельзя найти формулу первообразной в элементарных функциях. Общий метод численного интегрирования сводится к замене интегрируемой площади подынтегральной функции на элементарные площади, получаемые разбиением с заданным шагом , где – квадратурная сумма, -коэффициенты, -узлы квадратурной функции.
Формула прямоугольников.
(4.4.4), где (k- число разбиений).
Иллюстрация метода прямоугольников.
Формула трапеций.
(4.4.5)
где h-шаг разбиения
Ошибка вычисления интеграла определяется как
, где
Иллюстрация метода трапеций.
Формула Симпсона.
(5.4.6 ), где (n-число разбиений)
Ошибка вычисления интеграла
где , четвертая производная.
Рассмотрим численный метод вычисления интеграла с шагом 0.01 по формуле Симпсона.
Интеграл будет иметь вид :
=
Пример. Написать алгоритм и программу вычисления интеграла на основе метода прямоугольников и найти первую производную по двум узлам и вторую производную выражения по двум узлам 2 и 4.
Обозначим в схеме выражения (4.4.2), (4.4.3), (4.4.4) соответственно (a), (b), (c)
var
delta:real;
x0,x1,n,proiz1,proiz2,integ:real;
i : integer;
function fun(x:real):real;
begin
fun:=12*x*x+exp(x);
end;
Begin
write('Введите x0 и x1:');
readln(x0,x1);
n:=x0;
proiz1:=fun(x1)-fun(x0)/(2*abs(x1-x0)); {вычисление первой производной}
proiz2:=(fun(x0)-2*fun(x0+abs(x1-x0)/2)+fun(x1))/(abs(x1-x0)*abs(x1-x0));
{вычисление второй производной}
delta:=abs(x1-x0)/10;
integ:=0;
for i:=1 to 10 do begin
integ:=integ+delta*(fun(x0)+fun(x0+delta))/2; {вычисление интеграла}
x0:=x0+delta;
end;
writeln('Первая производная=',proiz1:6:2,’ в точке x=’,(x1-n)/2);
writeln('Вторая производная=',proiz2:6:2,’ в точке x=’,(x1-n)/2);
writeln('Интеграл=',integ:6:2);
readln;
End.
Результаты задачи 4 :
Введите x0 и x1 : 2 4
Первая производная=232.75 в точке x=3
Вторая производная=11.45 в точке x=3
Интеграл=271.53
4.5.Методы решения дифференциальных уравнений. Метод Рунге-Кутта
Дифференциальным уравнением 1-го порядка называется уравнение вида:
(4.5.1)
Решить дифференциальное уравнение это значит определить функцию y=F(x), которая при подстановке в уравнение (4.5.1) приводит к тождеству. Поиск решения y=F(x) при реализации аналитических методов сводится к интегрированию уравнения (4.5.1), что не всегда приводит к успеху из-за сложности вида f(x,y). Поэтому для решения задач (4.5.1) используют численные методы.
Метод Эйлера. При реализации метода Эйлера для уравнения (4.5.1) первую производную заменяют приближенным соотношением:
, (4.5.2)
которое является первым приближением формулы Тейлора:
Таким образом, уравнение (4.5.2) преобразуется к виду:
(4.5.3)
где .
Для решения уравнения (4.5.1) с помощью приближения (4.5.3) вводят начальное условие Коши .
Поясним применимость формулы (4.5.3) для решения дифференциального укравнения на отрезке [1,3] при начальном условии y(x=1)=2.
Выберем шаг итерации 0.5:
y(1)=2;
Отсюда
y(1)=2;
y(1.5)=2+0.5*(2*2+e)=5.39;
y(2)=5.39+0.5*(2*5.39+ )=13.1;
y(2.5)=13.1+0.5*(2*13.1+ )=43.17;
y(3)=43.17+0.5*(2*43.17+ )=92.79.
Результаты вычислений представлены в таблице
k |
|
|
|
|
0 |
1 |
2 |
- |
- |
1 |
1.5 |
5.39 |
3.39 |
2 |
2 |
2 |
13.1 |
7.71 |
5.39 |
3 |
2.5 |
43.17 |
30.07 |
13.1 |
4 |
3 |
92.79 |
49.62 |
43.17 |
Метод Рунге-Кутта четвертого порядка.
Итерационная формула в этом методе имеет вид:
(4.5.4),
где ; ; ; .
Поясним применимость метода Рунге-Кутта четвертого порядка для решения уравнения на отрезке [1,3] с начальными условиями y(x=1)=2.
Выберем шаг 1:
;
;
;
;
;
для следующего шага:
;
;
;
;
;
Обозначим . Результаты вычислений представлены в таблице:
k |
|
|
|
|
0 |
1 |
2 |
- |
- |
1 |
2 |
26.27 |
24.27 |
2 |
2 |
3 |
201.3 |
175.03 |
26.27 |
Для решения дифференциального уравнения более высокого порядка используется следующий прием. Любое уравнение n-го порядка можно свести к системе n уравнений первого порядка способом замены переменных. Для этого вводят новые переменные , в результате чего получают эквивалентную систему:
Указанный метод позволяет свести решение дифференциального уравнения n-го порядка к решению системы n уравнений первого порядка. В свою очередь методы решения одного уравнения первого порядка распространяются на систему таких уравнений.
Рассмотрим данный метод на примере дифференциального уравнения с начальными условиями , на отрезке [1,3].
Выполним замену переменной и приходим к системе из двух уравнений:
Таким образом любое уравнение n-порядка можно свести к системе уравнений первого порядка.
Пример. Написать алгоритм и программу решения дифференциального уравнения на отрезке [0,2] методом Рунге-Кутта четвертого порядка. Начальное условие Коши y(0)=1.
Обозначим выражение (4.5.4) в схеме алгоритма как (a).
program Runge_Cutta;
uses crt;
var y,x: array [1..11] of real;
k1,k2,k3,k4,y0,h,f0,fn:real;
i:integer;
function fun(x,y:real):real;
begin
fun:=2*y+exp(x);
end;
begin {н.п.}
clrscr;
f0:=0;
fn:=2;
h:=abs(fn-f0)/10; {шаг}
y0:=1;k1:=0;k2:=0;k3:=0;k4:=0;
k1:=fun(f0,y0);
k2:=fun(f0+h/2,y0+h*k1/2);
k3:=fun(f0+h/2,y0+h*k2/2);
k4:=fun(f0+h/2,y0+h*k3/2);
y[1]:=y0+(h/6)*(k1+2*k2+2*k3+k4); {вычисление функции в первой
точке}
x[1]:=f0;
f0:=f0+h;
for i:=2 to 11 do begin
k1:=fun(f0,y[i-1]);
k2:=fun(f0+h/2,y[i-1]+h*k1/2);
k3:=fun(f0+h/2,y[i-1]+h*k2/2);
k4:=fun(f0+h/2,y[i-1]+h*k3/2);
y[i]:=y[i-1]+(h/6)*(k1+2*k2+2*k3+k4); {вычисление функции в
остальных точках}
x[i]:=f0;
f0:=f0+h;
end;
for i:=1 to 11 do
writeln('x[',i,']=',x[i]:2:1,' y[',i,']=',y[i]:5:3);
readln;
end.{к.п.}
Результаты задачи 5:
x[1]=0.0 y[1]=1.733 x[2]=0.2 y[2]=2.870
x[3]=0.4 y[3]=4.618 x[4]=0.6 y[4]=7.282
x[5]=0.8 y[5]=11.315 x[6]=1.0 y[6]=17.391
x[7]=1.2 y[7]=26.510 x[8]=1.4 y[8]=40.151
x[9]=1.6 y[9]=60.505 x[10]=1.8 y[10]=90.814
x[11]=2.0 y[11]=135.871