МАТЕМАТИЧЕСКАЯ ОБРАБОТКА ИНФОРМАЦИИ Методические указания к выполнению лабораторных работ
для студентов специальности
5В074600 – Космическая техника и технологии
Алматы 2016
Некоммерческое акционерное общество
Кафедра системы управления
аэрокосмической техникой
АЛМАТЫ ЭНЕРГЕТИКА ЖӘНЕ БАЙЛАНЫС УНИВЕРСИТЕТІ
АЛМАТИНСКИЙ УНИВЕРСИТЕТ ЭНЕРГЕТИКИ И СВЯЗИ
2
СОСТАВИТЕЛИ: М.К. Шимырбаев, Д.Т. Тулекенова. Математическая обработка информации. Методические указания к выполнению лабораторных работ для студентов специальности 5В074600 – Космическая техника и технологии. – Алматы: АУЭС, 2016. – 41 с.
Данная работа представляет собой методические указания для студентов, изучающих дисциплину «Математическая обработка информации». Методические указания включают теоретический материал и задания для лабораторных и контрольных работ, позволяющие освоить основные принципы работы с пакетом MATLAB: арифметические вычисления, работа с массивами, построение графиков и поверхностей, решение уравнений и систем уравнений, вычисление интегралов, пределов, производных функций, решение дифференциальных уравнений.
Ил.24, табл.5, библиогр. – 6 назв.
Рецензент канд.экон.наук, доцент каф. ИС Ж.Г. Аренбаева
Печатается по плану издания некоммерческого акционерного общества
«Алматинский университет энергетики и связи» на 2016 г.
©НАО “Алматинский университет энергетики и связи”, 2016 г.
3
Содержание
Введение ... 4
1 Лабораторная работа №1. Построение двумерных графиков. ... 5
2 Лабораторная работа №2.Построение трёхмерных графиков ... 9
3 Лабораторная работа №3. Решение системы линейных алгебраических уравнений ... 12
4 Лабораторная работа №4. Решение задач алгебры и анализа ... 14
5 Лабораторная работа №5. Решение обыкновенных дифференциальных уравнений ... 24
6 Лабораторная работа №6. Решение ОДУ с помощью m-file ... 28
7 Лабораторная работа №7. Моделирование в Simulink ... 35
Список литературы ... 40
4 Введение
MATLAB – одна из старейших, тщательно проработанных и проверенных временем систем автоматизации математических расчетов, построенная на расширенном представлении и применении матричных операций. Однако синтаксис языка программирования системы настолько продуман, что эта ориентация почти не ощущается теми пользователями, которых не интересуют непосредственно матричные вычисления. Матрицы широко применяются в сложных математических расчетах, например, при решении задач линейной алгебры и математического моделирования статических и динамических систем и объектов. Они являются основой автоматического составления и решения уравнений состояния динамических объектов и систем.
В MATLAB реализованы современные численные методы компьютерной математики, используются мощные средства графической визуализации и анимационной графики. Возможности системы весьма обширны, а по скорости выполнения задач она нередко превосходит своих конкурентов и применима для расчетов практически в любой области науки и техники.
5
1 Лабораторная работа №1. Построение двумерных графиков 1.1 Основные теоретические сведения
Функции одной переменной у(х) находят широкое применение в практике математических и других расчетов, а также в технике компьютерного математического моделирования. Для отображения таких функций используются графики в декартовой (прямоугольной) системе координат. При этом обычно строятся две оси – горизонтальная X и вертикальная Y, и задаются координаты х и у, определяющие узловые точки функции у(х). Эти точки соединяются друг с другом отрезками прямых, т.е.
при построении графика осуществляется линейная интерполяция для промежуточных точек. Поскольку MATLAB–матричная система, совокупность точек у(х) задается векторами X и Y одинакового размера.
При построении графиков появляется графическое окно. Иногда оно бывает скрыто ранее имеющимися окнами как системы MATLAB, так и других приложений. Если вы не увидели окна графика, то поищите его в списке открытых окон (приложений) на панели задач или с помощью клавиш [Alt]+[Tab].
Порядок работы при построении графика функции следующий:
1) Задать значения аргумента функции.
2) Задать функцию.
3) Построить график.
4) Отформатировать график.
5) Добавить на график дополнительные элементы.
Для построения графиков функций в MATLAB служит команда plot, имеющая несколько вариантов записи (x – аргумент функции, y – функция):
- plot(x,y) – строит график одной функции;
- plot(x,y,s) – строит график функции с заданным типом и цветом линии и точек (s – строковая константа);
- plot(x,y1,x,y2,…) – строит графики нескольких функций в одной системе координат;
- plot(x,y1,s1,x,y2,s2,…) – строит графики нескольких функций в одной системе координат c заданным типом и цветом линии и точек.
С помощью строковой константы s можно изменять цвет линии, представлять узловые точки различными отметками (точка, окружность, крест и т. д.) и менять тип линии графика. Значения строковой константы представлены в таблице 1-3.
6 Таблица 1 - Цвет линии
Код Описание Код Описание
Y Желтый G Зеленый
M Фиолетовый B Синий
C Голубой W Белый
R Красный K Черный
Таблица 2 - Тип точки
Код Описание Код Описание
. Точка D Ромб
0 Окружность V Треугольник
Х Крест < Треугольник
+ Плюс > Треугольник
* Звездочка P Пятиугольник
S Квадрат H Шестиугольник
Таблица 3 - Тип линии
Код Описание
- Сплошная
-. Штрихпунктир
-- Штриховая
Пример. Построить график функции y=sinx на отрезке [–4; 4], шаг 0,2.
Порядок ввода:
>> x=–4:0.2:4;
>> y=sin(x);
>> plot(x,y).
Пример. Построить график этой же функции штриховой линией фиолетового цвета, отметив точки ромбами.
Порядок ввода:
>> x=–4:0.2:4;
>> y=sin(x);
>> plot(x,y,'dm--').
В результате каждого построения получатся графики, рисунок 1.
7
Рисунок 1 – Графики функции y=sinx стандартного вида и с заданными параметрами
Пример. Построить в одной системе координат графики функций y=sinx и z=cosx на отрезке [–5; 5] с шагом 0,2.
Порядок ввода:
>> x=–5:0.2:5;
>> y=sin(x);
>> z=cos(x);
>> plot(x,y,'-.+r',x,z,'--ok').
Рисунок 2 – Графики функций y=sinx и z=cosx 1.2 Добавление объектов на график
Дополнительно на график можно нанести различные объекты с помощью соответствующих команд пункта меню Insert:
8 X Label – подпись оси Х;
Y Label – подпись оси У;
Title – заголовок графика;
Legend – легенда;
Line – линия;
Arrow – стрелка;
Text Arrow – стрелка с текстом;
Double Arrow – двойная стрелка;
TextBox – текстовое поле;
Rectangle – прямоугольник;
Ellipse – эллипс;
Axes – оси.
Пример. Построить и отформатировать графики функций y=2𝑥2 + 3𝑥2 и z=3|3x–5| на отрезке [–5; 5] с шагом 0,1.
Порядок ввода:
>> x=–5:0.1:5;
>> y=2*x.^3+3*x.^2;
>> z=3*abs(3*x–5);
>> g=plot(x,y,x,z)
>> set(g,'linewidth',3,'linestyle','--')
>> title('Plots')
>> xlabel('x')
>> ylabel('y, z')
>> grid
>> legend('y','z')
>> text(–2.5,80,'z')
>> line([–2.4 –2.2],[70 50])
>> text(–1.5,–50,'y')
>> line([–1.4 –1.2],[–40 5])
Рисунок 3 – Графики функций y=2𝑥2 + 3𝑥2 и z=3|3x–5|
9 1.3 Содержание отчета
1. Цель работы.
2. Пример расчета и вывода данных.
1.4 Контрольные вопросы
1. Перечислите основные команды MATLAB для построения графиков.
2. Напишите команды пункта меню Insert.
1.5 Варианты заданий
Построить графики функций в одной системе координат, отформатировав их с помощью окна свойств графика по образцу:
1) f= lnx + 𝑥2, x∈[1; 7], шаг 0,4;
2) f= 𝑥2, y=sinx, x∈ [–5; 5], шаг 0,5;
3) f= sin𝑥2 – cosx, y= 𝑥2 – 3, x∈ [–4; 4], шаг 0,3;
4) f= sin𝑥2 + cosx, y= 𝑥2 – 4, z= x ∈ 0,8, x∈ [–4; 4], шаг 0,4.
5) x(t)=t∙cost, y(t)=t∙sint, t ∈ [0; 10], шаг 1/10;
6) b=3, x(t)=b∙cos3t, y(t)=b∙sin3t, t ∈ [0; 20], шаг 5/12;
7) a=4, x(t)=a∙(𝑡2 – 1)/( 𝑡2 + 1), y(t)=a∙t∙( 𝑡2– 1), t∈ [–10; 10], шаг 0,5 8) x(t)=t∙cost, y(t)=t∙sint, t ∈ [0; 10], шаг 1/10;
9) f=ln|x+2,5|+1, х∈ [–5; 5], шаг 1;
10) f= sin𝑥3 + cosx, y= 𝑥2 – 6, z= x ∈ 0,1, x∈ [–5; 5], шаг 0,01.
2 Лабораторная работа №2. Построение трёхмерных графиков 2.1 Основные теоретические сведения
Построение трехмерных графиков (поверхностей) во многом похоже на построение двумерных графиков. Для этого используется команда plot3, которая имеет несколько вариантов записи:
- plot3(X,Y,Z) – строит поверхность по точкам, координаты которых берутся из матриц X, Y, Z;
- plot3(X1,Y1,Z1,Х2,Y2,Z2,...) – строит несколько поверхностей Z1, Z2 и т. д.;
- plot3 (X,Y,Z,S) – строит поверхность заданным типом и цветом линии и точек (S – строковая константа, задающая тип и
цвет линии и точек);
- plot3(X1,Y1,Z1,S1,Х2,Y2,Z2,S2,Х3,Y3,Z3,S3,...) – строит несколько поверхностей Z1, Z2 и т.д. заданным типом и цветом линии и точек (S1, S2 и т.д.).
Значения строковой константы S приведены ранее (таблицы 4-6).
Особенно наглядное представление о поверхностях дают сетчатые графики, использующие функциональную закраску ячеек. Например, цвет окраски поверхности z(x, у) может быть поставлен в соответствие с высотой z
10
поверхности с выбором для малых высот темных тонов, а для больших – светлых. Для построения таких поверхностей используются команды класса surf:
- surf(X,Y,Z) – строит цветную параметрическую поверхность по данным матриц X, Y и Z;
- surfc(X,Y,Z) – строит цветную параметрическую поверхность по данным матриц X, Y, Z и проекцию фигуры на опорную плоскость.
Иногда бывают полезны графики трехмерных слоеных поверхностей, как бы состоящие из тонких пластинок – слоев.
Такие поверхности строит функция waterfall:
waterfall(X,Y,Z) – строит поверхность, состоящую из тонких пластинок – слоев.
Порядок построения трехмерных графиков следующий:
1) Задать матрицы X и Y на основе диапазонов значений переменных x и y с помощью команды преобразования диапазонов значений переменных в соответствующие матрицы:
[X,Y]=meshgrid(диапазон1, дипазон2);
если диапазоны одинаковые, то [X,Y]=meshgrid(диапазон);
2) Задать функцию Z(X,Y).
3) Построить поверхность нужного вида с помощью соответствующей команды.
4) Отформатировать график.
Отформатировать трехмерный график можно с помощью окна свойств графика или специальных команд (см. построение двумерных графиков).
Кроме того, для форматирования цветных поверхностей есть дополнительные команды:
- colormap(gray) – задает окраску тонами серого цвета;
- shading interp – устраняет изображения линий и задает интерполяцию для оттенков цвета поверхности;
- colorbar – выводит на экран цветовую шкалу.
Пример. Построить график функции z(х,у)=𝑥2+ 𝑦2 на отрезке [–3; 3] с шагом 0,15.
Порядок ввода:
>> [X,Y]=meshgrid([–3:0.15:3]);
>> Z=X.^ 2+Y.^2;
>> plot3(X,Y,Z);
>> grid.
Пример. Построить цветной график функции g(х,у)=5x∙sinx–3,5𝑦2 и его проекцию на отрезке [–3; 3] с шагом 0,15.
Порядок ввода:
>> [X,Y]=meshgrid([–3:0.15:3]);
>> G=5∙X.∙sin(Y)–3.5∙Y.^2;
11
>> surfc(X,Y,G).
Рисунок 4 – Графики функций z(х,у)= 𝑥2+ 𝑦2 и g(х,у)=5x∙sinx–3,5 𝑦2 2.2 Содержание отчета
1. Цель занятия.
2. Листинг программ и результаты выполнения программ.
2.3 Контрольные вопросы
1. С помощью какой команды осуществляется построение графиков декартовой системах координат?
2. Для чего служит команда mesh?
3. Как осуществляется задание надписей?
4. Для чего используется команда grid?
5. Как осуществляется разбивка окна на меньшие окна?
6. Для чего используется команда hold?
2.4 Варианты заданий
Построить цветные поверхности функции на заданных отрезках и отформатировать их по образцу:
1) z=2xsinx+3ycosy на отрезке [–2; 2], шагом 0,2;
2) z=2xsinx+3ycosy на отрезке [–5; 5], шагом 0,5;
3) z=𝑥2+𝑦2 на отрезке [–2; 2] с шагом 0,2;
4) g = 2𝑥2+𝑦
3 ∙ (3 cos 3𝑦+sin 𝑦) на отрезке [–4; 4] с шагом 0,2;
5) 𝑥2
10+ 𝑦2
5= 2𝑧 на отрезке [–4; 4] с шагом 0,2;
6) x2
8 − 𝑦2
8= 2𝑧 на отрезке [–4; 4] с шагом 0,2;
7) 𝑥
2
50 + 𝑦2
50+𝑧2
50= 1 на отрезке [–5; 5], шагом 0,5;
12 8) f = 2𝑦2+𝑥
5 ∙ (5 cos 5𝑦+sin 𝑦) на отрезке [–5; 5] с шагом 0,5;
9) z=5xsinx+2ycosy на отрезке [–3; 3], шагом 0,3;
10) 𝑥
2
8 + 𝑦2
4= 4𝑧 на отрезке [–5; 5] с шагом 0,5.
3 Лабораторная работа №3. Решение системы линейных алгебраических уравнений
3.1 Основные теоретические сведения. Матричный способ решения систем линейных уравнений
Решение систем линейных уравнений относится к самой массовой области применения матричных методов. Как известно, обычная система линейных уравнений имеет вид:
{
𝑎11𝑥1+ 𝑎12𝑥2 + ⋯ + 𝑎1𝑛𝑥𝑛 = 𝑏1; 𝑎21𝑥1+ 𝑎22𝑥2 + ⋯ + 𝑎2𝑛𝑥𝑛 = 𝑏2;
… … … . .
𝑎𝑛1𝑥1 + 𝑎𝑛2𝑥2 + ⋯ + 𝑎𝑛𝑛𝑥𝑛 = 𝑏𝑛.
(1)
Рассмотрим различные способы решения систем уравнений.
Пусть А – матрица коэффициентов при неизвестных;
В – вектор-столбец свободных членов;
Х – вектор-столбец неизвестных, т.е.
A=
(
𝑎𝑎1121 𝑎 𝑎1222 … 𝑎 … 𝑎1𝑛2𝑛… … … … 𝑎𝑛1 𝑎𝑛2 … 𝑎𝑛𝑛
)
, B=(
𝑏𝑏1…2 𝑏𝑛
)
, X=(
𝑥𝑥…12𝑥𝑛
).
(2)Тогда систему уравнений можно записать в матричном виде А∙Х=В.
Решение системы уравнений имеет вид Х=А^.–1∙В, где А^–1 – матрица, обратная матрице А.
Пример. Решить систему уравнений { 2𝑥1− 𝑥2 = 4;
5𝑥1 + 2𝑥2 = 3;
Порядок ввода:
>> A=[2 –1;5 2]
>> В=[4; 3]
>> X=A^–1∙B
В результате получим х1=1,222, х2=–1,556.
13
3.2 Решение систем уравнений с помощью функции solve
Для решения систем линейных или нелинейных уравнений в MATLAB существует специальная функция solve. Чтобы решить систему уравнений с помощью этой функции, нужно:
1) Определить символьные переменные (неизвестные системы уравнений).
2) Вычислить неизвестные по формуле [х1,х2,…]=
solve('уравнение1','уравнение2',…).
3) Вывести найденное решение с заданной точностью, используя функцию vpa(переменная, число знаков).
Пример. Решить систему уравнений { 2𝑥1− 𝑥2 = 4;
5𝑥1+ 2𝑥2 = 3.
Порядок ввода:
>> syms x1 x2
>> [x1,x2]=solve('2*x1–x2=4','5*x1+2*x2=3');
>> vpa(x1,4)
>> vpa(x2,4)
В результате получим х1=1,222, х2=–1,556.
3.3 Содержание отчета
1. Перечислите способы решения СЛАУ.
2. Листинг программы для вывода графиков функций.
3.4 Контрольные вопросы
1. С помощью какой команды осуществляется построение графиков декартовой системах координат?
2. Для чего служит функция vpa?
3. Для чего используется команда solve?
3.5 Варианты заданий
1. Решить системы линейных уравнений матричным способом, методом Крамера и с помощью функции solve:
1) {
2𝑥1 + 𝑥2 + 𝑥4 = 8;
𝑥1−3𝑥2+ 2𝑥3 + 4𝑥4 = 9;
−5𝑥1− 𝑥3− 7𝑥4 = −5;
𝑥1−6𝑥2+ 2𝑥3+ 6𝑥4 = 0.
3) {
2𝑥1+ 3𝑥2− 4𝑥3 + 𝑥4 = 3,1;
0,1𝑥1−2𝑥2 − 5𝑥3 + 𝑥4 = 2;
0,15𝑥1− 3𝑥2 + 𝑥3− 4𝑥4 = 1;
10𝑥1 + 2𝑥2 − 𝑥3+ 2,1𝑥4 = −4,7.
14 2) {
3𝑥1− 𝑥2 = 5;
−2𝑥1+ 𝑥2− 𝑥3 = 0;
2𝑥1 − 𝑥2+ 4𝑥3 = 15.
4) {
4𝑥1+ 0,24𝑥2 − 0,08𝑥3 = 8;
0,09𝑥1 + 3𝑥2− 0,15𝑥3 = 9;
0,04𝑥1 − 0,08𝑥2+ 4𝑥3 = 20.
5) {2𝑥𝑦 = 16;
𝑥2 + 𝑦 = 8. 6) {𝑥2− 𝑦 = 3;
𝑥 + 2𝑦 = 4.
7) {𝑥3 − 2𝑦 = 1;
𝑥2+ 3𝑦 = 2. 8) {3𝑥3+ 2𝑥2 + 𝑦 = 5;
2𝑥2 − 3𝑥 + 𝑦 = 4.
2. Решить системы нелинейных уравнений с помощью функции solve:
1) {𝑥3− 2𝑦 = 1;
3𝑥 + 3𝑦2 = 2. 3) { √𝑥 + 𝑦 + √𝑥 − 𝑦4 = 8;
√𝑥3+ 𝑥2𝑦 − 𝑥𝑦2− 𝑦3
4 = 12.
2) {
𝑢𝑣𝑥2 = 8;
𝑣𝑥2𝑤 = 24;
𝑥3𝑤𝑢 = 12;
𝑢 + 𝑣 + 𝑤 = 𝑥 + 4.
4) {
𝑥3− 2𝑥2+ 3𝑦 + 𝑧 = 2;
2𝑥2− 𝑥 + 𝑦2 − 𝑧2 = 1;
𝑦3 − 2𝑧2+ 3𝑧 = 1.
4 Лабораторная работа №4. Решение задач алгебры и анализа 4.1 Основные теоретические сведения
В системе MATLAB для решения систем линейных уравнений предусмотрены знаки операций / и \. Чтобы решить систему линейных уравнений вида
A⋅y =B, (1) где A – заданная квадратная матрица размером NxN;
B – заданный вектор–столбец длины N, достаточно применить операцию
\ и вычислить выражение A\ B .
Операция \ называется левым делением матриц и, будучи примененная к матрицам A и B в виде A\B , примерно эквивалентна вычислению выражения inv(A)∙B. Здесь под inv(A) понимается вычисление матрицы, обратной к матрице A . Операцию / называют правым делением матрицы. Выражение A/B примерно соответствует вычислению выражения B∙inv(A). Значит, эта операция позволяет решать системы линейных уравнений вида y⋅A= B .
Решение уравнения F(x)=0, или нахождение нулей функции, осуществляется с помощью функции fzero(name, x0) . В качестве первого аргумента ей передается имя функции, задающей исходное уравнение, вторым
15
аргументом служит начальное приближение к корню. Определим, в качестве примера, нули функции cos(x) на отрезке от 0 до pi. В качестве начального приближения примем x0=1.
>> x=fzero('cos', 1);
x = 1.5708.
Если требуется найти корень функции, отличной от стандартной (встроенной в систему MATLAB) и тем самым не имеющей в рамках системы MATLAB фиксированного имени, то нужно приписать некоторое имя выражению, вычисляющему функцию.
Пусть, например, требуется найти корни уравнения cos(x)=x, что эквивалентно нахождению нулей функции, вычисляемой по формуле y=cos(x)−x, не имеющей в рамках системы MATLAB фиксированного имени.
В этом случае нужно создать Mat–функцию вида function y=MyFunction1(x);
y=cos(x)-x;
После этого можно воспользоваться функцией fzero:
>> x=fzero('MyFunction1',pi/2);
x = 0.7391.
Если найдено абсолютно точное значение корня, то значение функции в этой точке равно нулю. Таким образом, величина функции в приближенно найденном нуле косвенно характеризует погрешность результата. Чтобы управлять погрешностью, нужно осуществлять вызов функции fzero с тремя аргументами fzero(name, x0, tol ), где tol задает требуемую величину погрешности (ошибки). Необходимо отметить, что функция fzero находит нули только вещественнозначных функций одной вещественной переменной.
Однако часто бывает необходимо найти комплексные корни вещественнозначных функций, особенно в случае многочленов. Для этой цели в системе MATLAB существует специальная функция roots(p), которой в качестве аргумента передается массив коэффициентов многочлена (p).
Например, для многочлена px4 3x3 3x2 3x2 нужно сначала сформировать массив его коэффициентов (расположив по порядку убывания степени х):
>> p = [ 1 -3 3 -3 2 ]:
>> r = roots( p);
r = 2.0000
-0.0000 +1.0000i -0.0000 –1.0000i 1.0000.
В задаче о нахождении нулей функции сложным моментом является нахождение начального приближения к нулю функции, а также априорная оценка их количества. Поэтому важно параллельно с применением функций
16
типа roots или fzero визуализировать поведение искомых функций на том или ином отрезке значений аргумента.
В системе MATLAB имеются специальные функции для поиска минимума заданных функций. При этом возможен поиск минимума как для функции одной вещественной переменной, так и для функций многих переменных. Для функций одной переменной их минимумы разыскивает функция xmin =fmin(name, x0, x1).
Здесь name представляет имя функции, у которой находятся минимумы, а x0 и x1 задают отрезок поиска. Для поиска минимума функции нескольких переменных применяется функция fmins: xmins = fmins( name, x0). Здесь name является именем функции нескольких переменных, для которой ищется минимум, а x0 – это вектор ее аргументов, с которого начинается поиск. Для иллюстративного примера создадим простую функцию двух переменных имеющую минимумом точку (0,0).
function у = MyFunc2(х) у = x(1)^2+ x(2)^2
После этого можно вызвать функцию fmins, которая приближенно находит вектор xmin координат точки минимума:
>> xmin = fmins( 'MyFunc2', [1,1] );
>> xmin(1);
ans =-2.1023e-005;
>> xmin(2);
ans =2.5484е-005.
Для вычисления интегралов методом трапеций в системе MATLAB предусмотрена функция trapz: Integ = trapz(х, у). Точность вычисления интеграла зависит от величины шага интегрирования: чем меньше этот шаг, тем больше точность.
Вычислим простой интеграл
0
)
cos(x dx методом трапеций с шагом интегрирования рi/10.
>> dx = pi/10;
>>х = 0:dx:pi;
>> y=cos(x);
>>I1 = trapz(x,y);
I1 = 5.5511e-017.
Обычно для достижения высокой точности требуется выполнять интегрирование с очень малыми шагами, а контроль достигнутой точности осуществлять путем сравнения последовательных результатов. При одном и том же шаге интегрирования методы более высоких порядков точности достигают более точных результатов. В системе MATLAB методы интегрирования более высоких порядков точности реализуются функциями:
17
quad (метод Симпcона) и quad8 (метод Ньюто-на-Котеса 8-го порядка точности).
Двойные интегралы в системе MATLAB вычисляются с помощью специальной функции dblquad.
Вычислим интеграл вида
1 0 2
1
)) sin(
) sin(
(x y y x dxdy.
Оформим подынтегральную mat–функцию и вызовем функцию dblquad:
function z = integ(x, у) z = x.*sin(y) + y.*sin(x);
>>J = dblquad( 'integ', 0, 1, 1, 2 );
J = 1.1678.
Возможности встроенного пакета символьных вычислений и операции Symbolic Math Toolbox достаточно обширны, рассмотрим лишь некоторые его возможности. Символьный объект создается при помощи функции syms.
Команда создает три символьные переменные х, а и b
» syms х a b.
Конструирование символьных функций от переменных класса sym object производится с использованием обычных арифметических операций и обозначений для встроенных математических функций. Запись формулы для выражения в одну строку не всегда удобна, более естественный вид выражения выводит в командное окно функция pretty:
Рисунок 5 – Демонстрация работы функции pretty Упростим выражение
y x
y x
2
2
, используя функцию simplify.
18
Рисунок 6 – Демонстрация работы функции simplify
Символьную функцию можно создать без предварительного объявления переменных при помощи sym, входным аргументом которой является строка с выражением, заключенная в апострофы:
Рисунок 7 – Демонстрация создания символьной переменной без предварительного объявления
Упрощение тригонометрических, логарифмических, экспоненциальных функций и полиномов осуществляется функцией expand, формат обращения к которой имеет следующий вид: rez=expand(S), где S – символьное выражение, которое нужно упростить, rez – результат упрощения.
Например:
>> syms x y;
>> rezl=expand(sin(x+y));
rezl = sin (x) ∙cos (у) +cos (x) ∙sin (y).
19
С помощью функции factor можно раскладывать многочлены на простые множители, а целые числа – в произведение простых чисел
>> factor(sym('x^5 – 1') ) ans =
(х-1) ∙ (х^4+х^3+х^2+х+1)
Функция subs осуществляет подстановку новых выражений для указанных символьных переменных:
Рисунок 8 – Демонстрация работы функции subs
Symbolic Math Toolbox позволяет работать как с неопределенными интегралами, так и с определенными. Неопределенные интегралы от символьных функций вычисляются при помощи функции int, в качестве входных аргументов указываются символьная функция и переменная, по которой происходит интегрирование (рисунок 5).
Рисунок 9 – Демонстрация работы функции int
Для нахождения определенного интеграла в символьном виде следует задать нижний и верхний пределы интегрирования, соответственно, в третьем и четвертом аргументах int:
20
Рисунок 10 – Демонстрация работы функции int для нахождения определенного интеграла
Перечислим еще несколько функций, часто используемых при символьных вычислениях:
inv – вычисляет обратную матрицу;
limit –вычисляет пределы;
taylor – осуществляет разложение функций в ряд Тейлора;
solve – решает алгебраическое уравнение и систему алгебраических уравнений.
Для работы с символьными данными предусмотрена оболочка funtool.
Она представляет собой интерактивный графический калькулятор, позволяющий быстро построить графики двух функций f(x) и g(x). Интерфейс данного приложения представлен на рисунке 10. При запуске выводятся три автономных окна: два графически и управляющее.
Управляющее окно содержит два поля для ввода функций, поле ввода пределом изменения переменной x, поле ввода масштабирующего коэффициента.
21
Рисунок 11 – Интерфейс приложения funtool
4.2 Содержание отчета 1. Цель занятия.
2. Листинг программы и результаты выполнения.
4.3 Контрольные вопросы
1.Для чего служит функция fmin?
2. Для чего служит функция fzero?
4. Перечислите основные внешние расширения системы MATLAB для поиска экстремумов функций.
5. Какая функция служит для вычисления определенного интеграла?
6. Назовите функции работы с символьными переменными.
4.4 Варианты задании
1. Составить и отладить программы для нахождения корней уравнения 𝑓1(x) = 0 и 𝑓2 (x) = 0 и вывести графики функции на основании задания из таблицы 1.
22
2. Найти определенный интеграл для подынтегральной функции, заданной в таблице 2.
3. Найти определенный интеграл для подынтегральной функции, заданной в таблице 2 с использованием пакета символьных вычислений.
Таблица 1 – Варианты заданий для нахождения корней уравнения
№ варианта
)
1(x
f – полином 3-ей степени с
коэффициентами a f2(x)
1 2 3 4 5 6
1 0 -1 4 -1 0.2ex 20
2 0 2 -2 -15 40 cos(x)
3 0 1 4 -1 10log(x5.5)
4 0 9 -8 -70 100 sin(x)
5 0 -4 4 50 70cos(x)
6 1 -5 4 40 60e.0.1x 100
7 2 -3 2 30 20sin(2x)
8 3 -6 1 50 e.x sin(2x)
9 4 -9 1 70 e.x cos(3x)
10 5 -7 5 60 60 cos(x)
11 -1 -4 9 60 15log(x1.5)
12 -2 -6 -7 55 50lg(x1.5)
13 -3 -9 -8 75 100 cos(x)
14 -4 7 8 -75 20sin(x/2)
15 -5 1 4 -1 40cos(x/2)
23
Таблица 2 – Варианты функций для нахождения значения интеграла
№ варианта
Функция Интервал интегрирования начало
интервала
конец интервала
1 2 3 4
1 f(x)x2 ex -2 2
2 f(x)xcos(x) -2 2
3 f(x)sin(x)2cos(x) -2 2 4 f(x)cos(x)(1x2)1 -2 2 5 f(x)(x2)2 ln(x) -0.5 4.5
6 f(x)2xtg(x) -1.4 1.4
7 f(x)x2 20cos(x1) -5 5 8 f(x)20ln(x2 1)0.1x3 -5 15
9 f(x)x2 ex -4 2
10 f(x)2xcos(x) -2 2
11 f(x)3sin(x)cos(x) -2 2 12 f(x)x2 20cos(x1) -2 5 13 f(x)20ln(x2 1)0.1x3 -1 1
14 f(x)xcos(x) -2 2
15 f(x)(x2)2 ln(x) -2 2
24
5 Лабораторная работа №5. Решение обыкновенных дифференциальных уравнений
5.1 Основные теоретические сведения. Решение ОДУ с помощью команды dsolve
Для решения обыкновенных дифференциальных уравнений в символьном виде в MATLAB существует команда dsolve. Она может быть использована, если решение существует в аналитическом виде.
Практически это означает, что командой dsolve можно пользоваться только при поиске решения линейного дифференциального уравнения (или системы линейных уравнений).
Пример. Решить дифференциальное уравнение x dx/dt=-0,5 x с начальным условием x(0)=10. Построить график решения в интервале
[–0,5; 7].
Порядок ввода:
>> x=dsolve('Dx=–0.5*x','x(0)=10');
>> ezplot(x,[–0.5,7]);
>> grid.
В результате получим функцию х=10e–1/2t и график
Рисунок 12 – График функции-решения уравнения
Пример. Решить дифференциальное уравнение 2-го порядка с нулевыми начальными
2,5𝑑2𝑥
𝑑𝑡2 + 3𝑑𝑥
𝑑𝑡 + 5𝑥 = 12
25
условиями и построить график решения в интервале [–0,2; 9].
Порядок ввода:
>> x=dsolve('2.5D2x+3Dx+5x=12','Dx(0)=0','x(0)=0');
>> ezplot(x,[–0.2 9]);
>> grid.
Рисунок 13 – График функции-решения уравнения
Пример. Построить график решения дифференциального уравнения 3-го порядка с нулевыми начальными условиями в интервале [–0,2; 21]:
1.5𝑑3𝑥
𝑑𝑡3 + 4𝑑2𝑥
𝑑𝑡2 + 3𝑑𝑥
𝑑𝑡 = 12.
Порядок ввода:
>> x=dsolve('1.5D3x+4D2x+3Dx+5x=12','D2x(0)=0','Dx(0)=0','x(0)=0');
>> ezplot(x,[–0.2 21]);
>> grid.
Рисунок 14 – График решения уравнения
26
Пример. Решить неоднородную систему дифференциальных уравнений третьего порядка с нулевыми начальными условиями:
{
𝑑𝑥1
𝑑𝑡 = −𝑥1 + 10;
𝑑𝑥2
𝑑𝑡 = 2 ∙ 𝑥1 − 𝑥3;
𝑑𝑥3
𝑑𝑡 = 2,5 ∙ 𝑥1 + 3 ∙ 𝑥2 − 2 ∙ 𝑥3.
Порядок ввода:
>> [x1,x2,x3]=dsolve('Dx1= –x1+10', 'Dx2=2∙x1–x3', 'Dx3=2.5∙x1–3∙x2–
2∙x3','x1(0)=0','x2(0)=0','x3(0)=0').
5.2 Содержание отчета 1. Цель занятия.
2. Листинг программы и результаты выполнения.
5.3 Контрольные вопросы
1. Какие существуют методы решения систем дифференциальных уравнений?
2. Для чего служит функция ode?
5.4 Вариант задании
1. Решить дифференциальные уравнения старшего порядка при заданных начальных условиях:
1) 𝑑
2𝑥
𝑑𝑡2 + 4𝑥 = 0, 𝑥(0) = 0, 𝑑𝑥
𝑑𝑡 (0) = 2.
2) 𝑑
2𝑥
𝑑𝑡2 + 3𝑑𝑥
𝑑𝑡 + 2𝑥 = 0, 𝑥(0) = 1. 𝑑𝑥
𝑑𝑡 (0) = −1.
3) 𝑑
2𝑥
𝑑𝑡2 + 2𝑑𝑥
𝑑𝑡 = 0, 𝑥(0) = 1, 𝑑𝑥
𝑑𝑡 (0) = 0.
4) 𝑡𝑑2𝑥
𝑑𝑡2 = 𝑑𝑥
𝑑𝑡, 𝑥(0) = 0, 𝑑𝑥
𝑑𝑡 = 0.
5) (1 − 𝑡2)𝑑2𝑥
𝑑𝑡2 − 2𝑡𝑑𝑥
𝑑𝑡 = 0, 𝑥(0) = 0, 𝑑𝑥
𝑑𝑡 (0) = 3.
6) 𝑑
2𝑥
𝑑𝑡2(1 − 3𝑙𝑛𝑡) +1𝑑𝑥
𝑡𝑑𝑡 = 3 + 𝑙𝑛𝑡, 𝑥(1) = 0.8, 𝑑𝑥
𝑑𝑡 = −1.
27 7) 𝑑2𝑥
𝑑𝑡2 − 6𝑥 = 0, 𝑥(1) = 0, 𝑑𝑥
𝑑𝑡(0) = 1.
8) 5𝑡𝑑2𝑥
𝑑𝑡2 =𝑑𝑥
𝑑𝑡 , 𝑥(0) = 0, 𝑑𝑥
𝑑𝑡 = −1.
9) (1 − 𝑡2)𝑑2𝑥
𝑑𝑡2 − 6𝑡𝑑𝑥
𝑑𝑡 = 0, 𝑥(1) = 0, 𝑑𝑥
𝑑𝑡 (0) = −1.
10) 𝑑2𝑥
𝑑𝑡2(1 − 𝑙𝑛𝑡) +1𝑑𝑥
𝑡𝑑𝑡 = 3 + 𝑙𝑛𝑡, 𝑥(0) = 0.5, 𝑑𝑥
𝑑𝑡 = −1.
2. Решить системы дифференциальных уравнений при заданных начальных условиях и построить графики решения:
1) {
𝑑𝑥1
𝑑𝑡 = 2 ∙ 𝑥1 − 1;
𝑑𝑥2
𝑑𝑡 = 3 ∙ 𝑥1 + 2 ∙ 𝑥2 + 1; 𝑥1(0) = 1; 𝑥(2) = 1.
2) {
𝑑𝑥1
𝑑𝑡 = 𝑥1 + 𝑥2;
𝑑𝑥2
𝑑𝑡 = 2 ∙ 𝑥1 + 2; 𝑥1(1) = 0; 𝑥2(1) = 1.
6 Лабораторная работа №6. Решение ОДУ с помощью m-file 6.1 Основные теоретические сведения. Методы решения ОДУ
Для решения систем ОДУ в MatLAB реализованы различные методы.
Их реализации названы решателями ОДУ. Решатели реализуют следующие методы решения систем дифференциальных уравнений:
- ode45 – одношаговые явные методы Рунге-Кутта 4-го и 5-го порядка.
Это классический метод, рекомендуемый для начальной пробы решения. Во многих случаях он дает хорошие результаты;
- ode23 – одношаговые явные методы Рунге-Кутта 2-го и 4-го порядка.
При умеренной жесткости системы ОДУ и низких требованиях к точности этот метод может дать выигрыш в скорости решения;
- ode133 – многошаговый метод Адамса-Башворта-Мултона переменного порядка. Это адаптивный метод, который может обеспечить высокую точность решения;
- ode15s – многошаговый метод переменного порядка (от 1-го до 5-го, по умолчанию 5), использующий формулы численного дифференцирования.
Это адаптивный метод, его стоит применять, если решатель ode45 не обеспечивает решения;
- ode23s – одношаговый метод, использующий модифицированную формулу Розенброка 2-го порядка. Может обеспечить высокую скорость вычислений при низкой точности;
28
- ode23t – метод трапеций с интерполяцией. Этот метод дает хорошие результаты при решении задач, описывающих осцилляторы с почти гармоническим выходным сигналом;
- ode23tb – неявный метод Рунге-Кутта в начале решения и метод, использующий формулы обратного дифференцирования 2-го порядка в последующем. При низкой точности этот метод может оказаться более эффективным, чем ode15s.
Все решатели могут решать системы уравнений явного вида:
y’ = F(t, y). Решатели ode15s, ode23s, ode23t, ode23tb могут решать уравнения неявного вида My’ = F(t, y).
В описанных далее функциях для решения систем дифференциальных уравнений приняты следующие обозначения и правила:
- F – название ODE-файла, то есть функции от t и y, которая возвращает вектор-столбец;
- tspan – вектор, определяющий интервал интегрирования [to tfinal].
Для получения решений в конкретные моменты времени to, t1, …, tfinal (расположенные в порядке уменьшения или увеличения), нужно использовать tspan = [t0 t1 … tfinal];
- y0 –вектор начальных условий;
- options – аргумент, создаваемый функцией odeset;
- p1, p2, … - произвольные параметры, передаваемые в функцию F;
- T, Y – матрица решений Y, где каждая строка соответствует времени, возвращенном в векторе-столбце T.
Перейдем к описанию функций для решения систем дифференциальных уравнений (будем обозначать понятием solver один из возможных численных методов решения ОДУ: ode45, ode23, ode113, ode15s, ode23s, ode23t, ode23tb):
- [T, Y] = solver(‘F’, tspan, y0) интегрирует систему дифференциальных уравнений вида y’ = F(t, y) на интервале tspan с начальными условиями y0. ‘F’
– строка, содержащая имя ODE-файла. Функция F(t, y) должна возвращать вектор-столбец. Каждая строка в массиве решений Y соответствует времени, возвращаемом в векторе-столбце T;
- [T, Y] = solver(‘F’, tspan, y0, options) дает решение, подобное описанному выше, но с параметрами, определяемыми значениями аргумента options, созданного функцией odeset. Обычно используемые параметры включают допустимое значение относительной погрешности RelTol (по умолчанию 1е-3) и вектор допустимых значений абсолютной погрешности AbsTol (все компоненты по умолчанию равны 1е-6);
- [T, Y] = solver(‘F’, tspan, y0, options, p1, p2, …) дает решение, подобное описанному выше, передавая дополнительные параметры р1, р2, … в m-файл F всякий раз, когда он вызывается. Используйте options = [], если никакие опции не задаются;
- [T, Y, TE, YE, IE] = solver(‘F’, tspan, y0, options) в дополнение к описанному решению содержит свойства Events, установленные в структуре options в положение on. ODE-файл должен быть создан так, чтобы вызов F(t,
29
y, ‘events’) возвращал соответствующую информацию. Выходной аргумент ТЕ – вектор-столбец времен, в которые происходят события, строки YE являются соответствующими решениями, а индексы в векторе IE определяют, какое событие произошло;
- [T, X, Y] = solver(‘model’, tspan, y0, options, ut, p1, p2, …) использует модель Simulink, вызывая соответствующий решатель из нее: [T, X, Y] = sim(solver, ‘model’, …).
Решатель систем ОДУ дает возможность получать решения систем из n уравнений. Система ОДУ может быть как однородной, так и неоднородной.
Решение сводится к следующему:
1) Создание m-файла. Независимо от вида системы он имеет вид:
function dy = solverDE(t, y);
dy = zeros(n, 1);
dy(1) = f1 (t, y(1), y(2), …, y(n));
dy(2) = f2 (t, y(1), y(2), …, y(n));
………
dy(n) = fn (t, y(1), y(2), …, y(n)).
2) Получение решения и сопровождающий его график:
>> [T, Y] = solver(‘solverDE’, [t0 tfinal], [y10 y20 … yn0];
>> plot(T, Y).
6.2 Линейные однородные уравнения с постоянными коэффициентами
Дифференциальное уравнение
𝑎𝑛𝑦(𝑛) + 𝑎𝑛−1𝑦(𝑛−1)+ . . . + 𝑎2𝑦′′+ 𝑎1𝑦′ + 𝑎0𝑦 = 0 , (1)
где 𝑎𝑛, 𝑎𝑛−1, … , 𝑎2, 𝑎1, 𝑎0 - постоянные величины, называется линейным однородным уравнением n - го порядка с постоянными коэффициентами.
𝑦 = 𝐶1𝑦2 + 𝐶2𝑦2+ 𝐶3𝑦3+ . . . +𝐶𝑛𝑦𝑛, (2) где 𝑦1, … , 𝑦𝑛 - его линейно независимые частные решения. Последние находятся с помощью характеристического уравнения
𝑎𝑛𝑟𝑛 + 𝑎𝑛−1𝑟𝑛−1+ . . . +𝑎2𝑟2+ 𝑎1𝑟 + 𝑎0 = 0. (3) Если характеристическое уравнение имеет n действительных различных корней 𝑟𝑚(𝑚 = 1,2, … , 𝑛), то каждому из них соответствует частное решение
𝑦𝑚 = 𝑒𝑟𝑛𝑥 (4) и общее решение уравнения (3) принимает вид: