Решение дифференциальных уравнений с помощью степенных рядов

ИНТЕРПОЛЯЦИОННЫЙ ПОЛИНОМ ЛАГРАНЖА

РАЗЛОЖЕНИЕ РАЦИОНАЛЬНОЙ ДРОБИ НА СУММУ ПРОСТЕЙШИХ ДРОБЕЙ

При нахождении (вручную) неопределенного интеграла

dx

возникает задача разложения подынтегральной дроби

на сумму простейших дробей. Получить такое разложение непосредственно в MATLAB нельзя. Следует использовать команду convert (с формой разложения parfrac) системы Maple:

>> maple('convert((x^2+2)^2*(x^3+3)/(x+1)/(x^2+1)^2,parfrac,x)')

ans =

x^2-x+3+9/2/(x+1)-1/2*(-7+9*x)/(x^2+1)-(-1+2*x)/(x^2+1)^2

>> pretty(sym(ans))

2 1 -7 + 9 x -1 + 2 x

x - x + 3 + 9/2 ----- - 1/2 ---------- - -----------

1 + x 2 2 2

1 + x (1 + x)

т.е.

= x² - x+3+ - -.

Предположим, что некоторая функция f(x) задана таблицей своих значений:

x	x0	x1	…	xn
y	y0	y1	…	yn

Требуется найти интерполяционный полином Лагранжа – многочлен L_n(x) степени не выше n, значения которого в точках x_k совпадают со значениями данной функции в этих точках, т.е. L_n(x_k) = y_k, k = 0,…, n.

Для нахождения интерполяционного полинома Лагранжа в Maple служит команда interp.

Пример:

Найти интерполяционный полином Лагранжа функции f(x)

x
y

Решение:

>> syms x

>> maple('interp([0,1,3,7],[5,4,2,1],x)')

ans =

1/56*x^3-1/14*x^2-53/56*x+5

Построим на рисунке 6 узлы интерполяции (команда stem) и график найденного интерполяционного полинома Лагранжа

L₃(x) = x³ - x² - x+5:

>> stem([0 1 3 7],[5 4 2 1],'fill')

>> hold

Current plot held

>> ezplot(ans,-1, 8)

>>grid

Рисунок 6

Как видим из рисунка 6, график найденного интерполяционного полинома Лагранжа проходит через узлы интерполирования.

Пример:

Решить дифференциальное уравнение y' = cosx+e^y с начальным условием y(0) = 1.

Решение:

Обращение к dsolve возвращает сообщение о том, что решение не найдено:

>> dsolve('Dy=cos(x)+exp(y)','y(0)=1','x')

ans =

[ empty sym ]

Команда dsolve не нашла аналитического решения в MATLAB. Известно, что решения этого дифференциального уравнения в аналитическом виде не существует. Найти разложение решения в степенной ряд (до 6 -й степени по умолчанию) можно с помощью команды dsolve системы Maple.

>> maple('dsolve({diff(y(x),x)=cos(x)+exp(y(x)),y(0)=1},y(x),series)')

ans =

y(x) = series(1+(1+exp(1))*x+(1/2*exp(1)*(1+exp(1)))*x^2+(-1/6+1/3*exp(1)*(3/2*exp(1)+exp(1)^2+1/2))*x^3+(1/2*exp(1)^3+1/4*exp(1)^4+7/24*exp(1)^2)*x^4+(1/2*exp(1)^4+1/5*exp(1)^5+5/12*exp(1)^3+1/120-1/40*exp(1)+1/12*exp(1)^2)*x^5+O(x^6),x,6)

Имеется возможность управлять порядком разложения. Найдем разложение решения в степенной ряд до 3 -й степени:

>> maple('Order:=3;dsolve({diff(y(x),x)=cos(x)+exp(y(x)),y(0)=1},y(x),series)')

ans =

Order:= 3y(x) = series(1+(1+exp(1))*x+(1/2*exp(1)*(1+exp(1)))*x^2+O(x^3),x,3)

Дата добавления: 2014-12-07