--PAGE_BREAK--2 Решение дифференциальных уравнений с помощью Mathcad
Mathcad имеет ряд встроенных функций, предназначенных для решения обыкновенных дифференциальных уравнений (ОДУ). При решении ОДУ искомой величиной является функция. При использовании любых методов численного интегрирования необходимо, чтобы были заданы по крайней мере следующие величины:
начальные условия;
набор точек в которых нужно найти решение;
само дифференциальное уравнение, записанное в некотором специальном виде, который будет описан ниже.
Один из наиболее эффективных алгоритмов интегрирования ОДУ основан на численном методе Рунге-Кутты четвертого порядка. Функция, реализующая этот метод, имеет вид rkfixed (y,x1,x2, npoints,D)
Здесь:
y-вектор начальных условий размерности n, где n- порядок дифференциального уравнения или число уравнений в системе (если решается система уравнений);
x1, x2– граничные точки интервала, на котором ищется решение дифференциального уравнения. Начальные условия, заданные в векторе y,- это значение решения в точке x1;
npoints— число точек (не считая начальной точки), в которых ищется приближенное решение. При помощи этого аргумента определяется число строк (1+npoints) в матрице, возвращаемой функцией rkfixed;
D(x,y) – функция, возвращающая значение в виде вектора n элементов, содержащих первые производные неизвестных функций.
2.1 Метод Эйлера
Если задачу об отыскании всех решений дифференциального уравнения удается свести к конечному числу алгебраических операций, операций интегрирования и дифференцирования известных функций, то говорят, что уравнениеинтегрируется в квадратурах. В приложениях крайне редко встречаются уравнения, интегрируемые в квадратурах. Поэтому для исследования дифференциальных уравнений широко используются приближенные, численные методы их решения.
Численное решениена отрезке [a, b] задачи Коши
y' =f(x, y), y(a) =y0
состоит в построении таблицы приближенных значений
y,y1, ...,yi, ...yN
решенияy(x)в узлах сетки
a=xx1xixN=b, y(xi)@yi.
Еслиxi=a+ i h, h=(b-a)/ N,то сетка называетсяравномерной.
Численный метод решения задачи Коши называетсяодношаговым, если для вычисления решения в точкеx0+hиспользуется информация о решении только в точкеx0.
Простейший одношаговый метод численного решения задачи Коши -метод Эйлера. В методе Эйлера величиныyiвычисляются по формуле
yi+1 =yi+hf(xi,yi), i= 0, 1
Найдем методом Эйлера на [0, 1] с шагом h=0.2 приближенное решение задачи Коши
Для того чтобы изменить стиль изображения, щелкните дважды по полю графиков и установите соответствующие параметры
Определим правую часть уравнения
Расчетные формулы метода Эйлера для решения этой задачи имеют вид
x0=0, y0= 1, xi+1 = xi + 0.2, yi+1 = yi + 0.2(sinxi — cosyi), i =0, 1, ..., 4.
Изобразим приближенное решение графически.
y' = sin x – cos y, y(0)=1.
Определим диапазон изменения номера точки i=0,1, ..., 4
Знак присваивания можно ввести щелчком по соответствующей позиции в панели Evaluation.
Определим начальное условие — решение в начальной точке
Для того чтобы ввести символ диапазона изменения индекса , щелкните по соответствующей позиции в панели Matrix или введите с клавиатуры символ («точка с запятой»)
Определим шаг формулы Эйлера — шаг интегрирования
Для того чтобы ввести нижний индекс переменной, щелкните по соответствующей позиции в панели Matrix или в панели Calculator
Определим по формулам Эйлера значения приближенного решения в узлах сетки
Выведем в рабочий документ вычисленные значения решения
Построим график найденного решения y(x)
Для того чтобы вывести значение переменной в рабочий документ, введите имя переменной, знак равенства и щелкните по рабочему документу вне выделяющей рамки
Для того чтобы построить график приближенного решения, щелкните в панели Graph по пиктограмме декартова графика, введите в помеченной позиции возле оси абсцисс обозначение компонент вектора, содержащего значения узлов сетки, а в позиции возле оси ординат — обозначение компонент вектора, содержащего значения приближенного решения в узлах сетки; затем щелкните по свободному месту в рабочем документе вне поля графиков.
2.2 Метод Эйлера с шагом h/2.
Метод Эйлера допускает простуюгеометрическую интерпретацию. Пусть известна точка (xi,yi) интегральной кривой уравненияy'=f(x, y).
Касательная к интегральной кривой уравнения, проходящая через эту точку, определяется уравнением
y=yi+f(xi,yi)(x-xi).
Следовательно, вычисленная методом Эйлера точка (xi+1 ,yi+1 ),
Гдеxi+1=xi+h,yi+1=yi+h f(xi,yi), лежит на этой касательной.
Найдем методом Эйлера на [0, 1] с шагом h=0.2 и с шагом h=0.1 приближенное решение задачи Коши
y' = sin x – cosy,y(0)=1.
Изобразим приближенные решения графически.
Расчетные формулы метода Эйлера для решения этой задачи имеют вид
x0=0, y0= 1, xi+1 = xi + 0.2, yi+1 = yi + 0.2(sinxi — cosyi), i =0, 1, ..., 4
xi+1 = xi + 0.2, yi+1 = yi + 0.2(sinxi — cosyi), i =0, 1, ..., 9
Определим правую часть уравнения
Знак присваивания можно ввести щелчком по соответствующей позиции в панели Evaluation.
Определим диапазон изменения номера точки i=0,1, ..., 4 для вычислений с шагом h=0.2
Для того чтобы ввести символ диапазона изменения индекса , щелкните по соответствующей позиции в панели Matrix или введите с клавиатуры символ («точка с запятой»)
При решении задачи с шагом h=0.2 назовем шаг h1, аргумент — x1, а решение — y1.
Определим начальное условие
Для того чтобы ввести нижний индекс переменной, щелкните по соответствующей позиции в панели Matrix или в панели Calculator
Определим шаг формулы Эйлера — шаг интегрирования
Определим по формулам Эйлера значения приближенного решения в узлах сетки
Выведем в рабочий документ вычисленные значения решения
Для того чтобы вывести значение переменной в рабочий документ, введите имя переменной, знак равенства и щелкните по рабочему документу вне выделяющей рамки
Построим график найденного решения y1(x1)
Для того чтобы построить график приближенного решения, щелкните в панели Graph по пиктограмме декартова графика, введите в помеченной позиции возле оси абсцисс обозначение компонент вектора, содержащего значения узлов сетки, а в позиции возле оси ординат — обозначение компонент вектора, содержащего значения приближенного решения в узлах сетки; затем щелкните по свободному месту в рабочем документе вне поля графиков.
Определим диапазон изменения номера точки i=0,1, ..., 9 для вычислений с шагом h=0.1
Для того чтобы ввести символ диапазона изменения индекса , щелкните по соответствующей позиции в панели Matrix или введите с клавиатуры символ («точка с запятой»)
При решении задачи с шагом h=0.1 назовем шаг h2, аргумент — x2, а решение — y2.
Определим начальное условие
Для того чтобы ввести нижний индекс переменной, щелкните по соответствующей позиции в панели Matrix или в панели Calculator
Определим шаг формулы Эйлера — шаг интегрирования
Определим по формулам Эйлера значения приближенного решения в узлах сетки
Выведем в рабочий документ вычисленные значения решения. Для сравнения рядом выведены значения решения, вычисленные с большим шагом
Построим график решения y2(x2)
Построим на одном графике оба приближенные решения
Для того чтобы одновременно построить графики нескольких функций от разных аргументов, щелкните в панели Graph по пиктограмме декартова графика, введите в помеченной позиции у оси абсцисс имя первого аргумента, запятую, имя второго аргумента, и т.д., разделяя имена аргументов запятой.
Аналогично, в позиции возле оси ординат введите имя функции первого аргумента, запятую, имя функции второго аргумента и т.д.разделяя имена функций запятой.
Когда функции определены, щелкните по рабочему документу вне поля графиков.
2.3 Метод Рунге – Кутты
Методом Рунге-Кутты четвертого порядкаточности называют одношаговый метод, относящийся к широкому классу методов Рунге-Кутты. В этом методе величиныyi+1 вычисляются по следующим формулам:
yi+1 =yi+h(k1+ 2k2+ 2k3+k4)/6, i= 0, 1, ...
k1=f(xi,yi),
k2=f(xi+h/2,yi+hk1/2),
k3=f(xi+h/2,yi+hk2/2),
k4=f(xi+h,yi+hk3).
Найдем на [0, 1]приближенноерешение задачи Кошиy' = sinx– cosy,y(0)=1методом Рунге-Кутты 4-го порядка с шагом h=0.2 и методом Эйлера с тем же шагом.Изобразим оба приближенные решения графически
Для решения задачиметодом Рунге-Кутты воспользуемся функциейrkfixed
Определим начальное условие — решение в начальной точке
продолжение
--PAGE_BREAK--