Решение обыкновенных дифференциальных уравнений

состоит в нахождении решения (интеграла) дифференциального уравнения (системы дифференциальных уравнений),

удовлетворяющего начальным условиям (начальным данным).
В задаче Коши область, в которой должно быть определено искомое решение, заранее не указывается.

или нет) и означаемым h-параметром
Модель:
разобьём ось х на промежутке

(0,1] на М равных частей и получим дискретное пространство

на дискретном пространстве h (на сетке h).

наша схема запишется таким образом:




Это схема Эйлера

шаблоном.
Шаблон, используемый на предыдущем слайде для аппроксимации – двухточечный.
Разностная формула

представляет собой разностную схему для модельной дифференциальной задачи на дискретном пространстве h (или просто на сетке).

в любой точке нашего дискретного пространства находится по формуле:
Значение сеточной

функции в последней точке промежутка:

Значение сеточной функции в произвольной точке промежутка:

уравнения стремиться, при уменьшении шага, к значению в этой точке

функции, полученной в результате решения исходного дифференциального уравнения, то говорят, что имеется
сходимость.

разностная схема имеет первый порядок точности.

на на дискретном пространстве h (на сетке h).
Типы разностей!!!!!

(1)
U(x) - искомое решение ( пусть непрерывное)
L

- заданный дифференциальный оператор
f - заданная функция (правая часть)
ω(h) - совокупность узлов разностной сетки, принадлежащая дискретному пространству D (считаем, что по любому направлению шаг постоянен).
h - шаг сетки
U(h) - сеточная функция
LhU(h)=f(h) , Lh – заданный разностный оператор (2)

Сравнению функций U(x) и U(h) препятствует их разная функциональная природа. Снесём U(x) на сетку:

проектирование решения на сетку

(1) если

в пространстве сеточных функций.
Если ,

то

имеет место сходимость порядка “k”.

где Fh - пространство сеточных функций правых частей, т.е.

и

δf(h) - погрешность аппроксимации или невязка или если существуют такие C1, k1, большие нуля, что

то имеет место аппроксимация порядка

y (tn) в правых частях дифференциальных уравнений f(t,y).

Классический метод

Рунге-Кутты используют для расчета стандартных моделей достаточно часто, так как при небольшом объеме вычислений он обладает точностью метода Ο4(h).

Данный метод не применим для жестких ОДУ.

В зависимости от типа экстраполяции будут получаться алгоритмы различного порядка

аппроксимации.

Методы Адамса-Бэшфорта могут быть как явными так и неявными.

Данный метод не применим для жестких ОДУ.

точную формулу интегрирования: 
Для построения схемы по значениям функции, вычисленным в

n предшествующих узлах, строится интерполяционный полином – Ln-1(x), который используется при интегрировании дифференциального уравнения. Интеграл при этом выражается через квадратурную формулу:

2-м предыдущим точкам
Схема 4-шагового алгоритма получается при использовании кубической экстраполяции

в точке h/2 проходит через состояние системы после двух таких

шагов, в точке h/4 проходит через состояние системы после четырех таких шагов, и т. д., а затем вычисляет значение этой функции в точке h = 0, проводя экстраполяцию. Гладкость правых частей приводит к тому, что вычисленное при помощи экстраполяции состояние системы оказывается очень близко к действительному, а использование рациональной экстраполяции вместо полиномиальной позволяет повысить точность.
После проведения одного шаг принимается решение, следует ли изменять шаг, а если да, то в какую сторону.
Считается, что метод Булирша— Штера часто оказывается более эффективным для поиска гладких решений, чем метод Рунге-Кутты. 

понятием жесткости системы: система будет жесткой, если на всем промежутке

интегрирования с помощью явных схем необходимо вести интегрирование с очень малым шагом.
Жесткость системы проявляется тогда, когда длина промежутка интегрирования, T, удовлетворяет соотношению :

где λmax — наибольшее по абсолютной величине собственное число матрицы системы, записанной в виде: y’=A∙y.

и явным. Однако при пересчет каждого шага требуется:
во-первых, численное

определение производных функции правых частей (в случае системы ОДУ правая часть - это вектор);
во-вторых, решения системы линейных уравнений (поскольку искомые компоненты вектора yn+1 входят в матричное уравнение в линейной комбинации).

алгоритм RADAU5 – алгоритм решения систем ОДУ любого типа (в

том числе и жестких) основанный на неявном методе Рунге-Кутты 5-го порядка с помощью квадратур Радау
адаптированный для жестких систем метод Булишера-Штера
алгоритм BDF – метод обратного дифференцирования

уравнений, линейных относительно старшей производной.
Оdesolve решает для записанных в

общепринятом в математической литературе виде дифференциальных уравнений:
a(x)y(n)+F(x,y,y',.y(n-1))=f(x)
задачу Коши y(x0 )=y0 , y'(x0 )=y0,1 , y''(x0 )=y0,2 , ..., y(n-1)(x0 )=y0,n-1
или простейшую граничную задачу y(k) (a)=ya,k y(m) (b)=yb,k , 0<= k<= n-1, 0<= m<= n-1.

Адамса/метод обратного дифференцирования
Для решения задачи методом Рунге-Кутты с фиксированным шагом

нужно щелкнуть в рабочем документе по имени функции правой кнопкой мыши и пометить во всплывающем меню пункт Fixed.
Для решения задачи методом Рунге-Кутты автоматическим выбором шага пометить во всплывающем меню пункт Adaptive.
Для решения жестких систем во всплывающем меню необходимо выбрать пункт Radau.

уравнений и начальных/граничных условий
Y:=Оdesolve(x,b,step)
Y - имя функции, содержащей значения

найденного решения,
x — переменная интегрирования,
b — конец промежутка интегрирования,
step — шаг, который используется при интегрировании уравнения. Данный параметр необязателен и может быть опущен.

(+).
Для записи производных можно использовать:
оператор дифференцирования из меню математического анализа.
знак

производной, например, вторую производную можно вводить в виде y''(x) (штрих записывается как +).
В любом случае необходимо обязательно записывать аргумент искомой функции.

точке промежутка интегрирования
Правая часть ОДУ
Ключевое слово начала блока решения
Запись ОДУ
Запись

граничного условия

Присвоение результатов решения ОДУ

можно так же с помощью следующих функций:
Adams(y, x1, x2,

npoints, D, {acc}) —решение задачи на отрезке методом Адамса;
rkfixed(y, x1, x2, npoints, D) —решение задачи на отрезке методом Рунге—Кутты с постоянным шагом;
Rkadapt(y, x1, x2, npoints, D) —решение задачи на отрезке методом Рунге—Кутты с автоматическим выбором шага;
Bulstoer(y, x1, x2, npoints, D) —решение задачи на отрезке методом Булирша-Штера .

помощью функций:

BDF(y, x1, x2, npoints, D, {J},{acc}) —решение задачи

на отрезке методом обратного дифференцирования (backward differentiation formula );
Radau(y, x1, x2, npoints, D, {J},{M},{acc}) —решение задачи Radau5;
Stiffr(y, x1, x2, D, AJ) — решение задачи для жестких систем на отрезке с использованием алгоритма Розенброка;
Stiffb(y, x1, x2, D, AJ) —решение задачи для жестких систем на отрезке с использованием алгоритма Булирша—Штера;

x2, npoints, D, {J},{acc})

В функции определяется тип системы

:
жесткая/нежесткая
и выбирается соответственный решатель.

Переключение осуществляется автоматически.

конечная точки отрезка интегрирования системы; для функций, вычисляющих решение в

заданной точке, x1 — начальная точка, x2 — заданная точка;
npoints — число узлов на отрезке [x1, x2]; при решении задачи на отрезке результат содержит npoints+1 строку;
D — имя вектор-функции D(x,y) правых частей, записанное через элементы искомого вектора;

первом столбце которой хранятся выражения частных производных по x правых

частей системы, а в остальных n столбцах содержится матрица Якоби правых частей:







Может быть вычислен с помощью функции Jacob (D(x,y),x)
acc — параметр, контролирующий погрешность решения при автоматическом выборе шага интегрирования (если погрешность решения больше acc, то шаг сетки уменьшается; шаг уменьшается до тех пор, пока его значение не станет меньше save );
AJ — матрица, в которой слиты два Якобиана- по аргументу и по неизвестной функции augment (Jacob (D(x,y),x), Jacob (D(x,y),y));

столбец содержит координаты узлов сетки,
второй столбец — вычисленные приближенные

значения решения y1 (x) в узлах сетки,
(k+1) -й — значения решения yk (x) в узлах сетки.

входящие в систему, должны иметь первый порядок (то есть содержать

только первые производные).
Все уравнения должны быть предварительно разрешены относительно производных и записаны в нормальной форме вида:

порядка уравнений путем замены переменных.
Приведение системы дифференциальных уравнений к

явному виду путем ее решения относительно производных

в Mathcad достигается применением двух встроенных функций. Первая предназначена для

двухточечных задач с краевыми условиями, заданными на концах интервала:
sbval (z,x0,x1, D, load, score)
— поиск вектора недостающих L начальных условий для двухточечной краевой задачи для системы N ОДУ.
z — вектор размера Lx1, присваивающий недостающим начальным условиям (на левой границе интервала) начальные значения;
х0 — левая граница расчетного интервала;
x1 — правая граница расчетного интервала;
load(x0,z) — векторная функция размера Nx1 левых граничных условий, причем недостающие начальные условия поименовываются соответствующими компонентами векторного аргумента z;
score (x1, у) — векторная функция размера Lx1, выражающая L правых граничных условий для векторной функции у в точке x1;
D(x,y) — векторная функция, описывающая систему N ОДУ, размера Nx1 и двух аргументов — скалярного х и векторного у. При этом у — это неизвестная векторная функция аргумента х того же размера Nx1.
 

Необходимо помнить, что число элементов векторов х0 и load равно

количеству уравнений N, а векторов z, score и результата действия функции sbval— количеству правых граничных условий L. Соответственно, левых граничных условий в задаче должно быть (N-L) .
Начальное значение фактически является параметром численного метода и поэтому может решающим образом повлиять на решение краевой задачи.

ВНИМАНИЕ! Метод стрельбы не годится для решения жестких краевых задач. Поэтому алгоритмы решения жестких ОДУ в Mathcad приходится программировать самому

границе интервала начальные значения.
z2 — вектор того же

размера, присваивающий недостающим начальным условиям на правой границе интервала начальные значения.
х0 — левая граница расчетного интервала.
x1 — правая граница расчетного интервала.
xf — точка внутри интервала.
Dо(х,у) — векторная функция, описывающая систему N ОДУ размера Nx1 и двух аргументов — скалярного х и векторного у. При этом у — это неизвестная векторная функция аргумента х того же размера Nx1.
load1(x0,z) — векторная функция размера Nx1 левых граничных условий, причем недостающие начальные условия поименовываются соответствующими компонентами векторного аргумента z.
load2(x1,z) — векторная функция размера Nx1 левых граничных условий, причем недостающие начальные условия поименовываются соответствующими компонентами векторного аргумента z.
score (xf, у) — векторная функция размера Nx1, выражающая внутреннее условие для векторной функции у в точке xf.