Проведение численных расчетов

МИНИСТЕРСТВО ОБРАЗОВАНИЯ И НАУКИ РОССИЙСКОЙ ФЕДЕРАЦИИ

ФЕДЕРАЛЬНОЕ ГОСУДАРСТВЕННОЕ ОБРАЗОВАТЕЛЬНОЕ УЧРЕЖДЕНИЕ ВЫСШЕГО ПРОФЕССИОНАЛЬНОГО ОБРАЗОВАНИЯ

ТЮМЕНСКИЙ ГОСУДАРСТВЕННЫЙ НЕФТЕГАЗОВЫЙ УНИВЕРСИТЕТ

ИНСТИТУТ ГЕОЛОГИИ И НЕФТЕГАЗОДОБЫЧИ

 

Кафедра «Автоматизации и управления»

 

МЕТОДИЧЕСКИЕ УКАЗАНИЯ

по дисциплине

«Идентификация и диагностика систем»

к лабораторной работе

«ОСНОВЫПРОГРАММИРОВАНИЯ В СРЕДЕ MATLAB.»

 

 

Составители: доц. Говорков Д.А.

 

Тюмень 2011

Содержание

 

Основные операторы среды 3

Проведение численных расчетов 6

Построение графиков в среде 9

Задания к лабораторной работе 12

 

Основные операторы

Представление данных в виде матриц (массивов)

	вектор-строка	вектор-столбец	матрица
операция	A=[1,2,3]	B=[1;2;3]	C=[1,2,3;4,5,6;7,8,9]
результат	A = 1 2 3	B =	C = 1 2 3 4 5 6 7 8 9

 

Обращение к элементу массива осуществляется через индексы: C(1,1) - первый элемент (в первом столбце, первой строке, индексы в MATLAB начинаются с 1, первый индекс – номер строки, второй индекс – номер столбца).

Операции с матрицами соответствуют операциям, принятым в матричном исчислении (при операциях с матрицами необходимо корректно учитывать их размерность):

	умножение матрицы на число (скаляр)	умножение вектора-строки на вектор-столбец	умножение вектора-столбца на вектор-строку
операция	C*2	A*B	B*A
результат	ans = 2 4 6 8 10 12 14 16 18	ans =	ans = 1 2 3 2 4 6 3 6 9

 

Определитель и ранг матрицы, транспонирование и обращение матриц (обратная матрица существует только при ненулевом определителе исходной матрицы).

Определитель матрицы	Ранг матрицы	Транспонирование матрицы	Обращение матрицы
det(C)	rank(C)	C'	inv(C)
ans =   6.6613e-016	ans =	ans = 1 4 7 2 5 8 3 6 9	ans = 1.0e+016 * -0.4504 0.9007 -0.4504 0.9007 -1.8014 0.9007 -0.4504 0.9007 -0.4504

(запись 1.0e+016 * в результате означает, что все элементы матрицы умножаются на это число)

Решение систем алгебраических уравнений вида , где - вектор неизвестных (искомых) параметров, - заданы.

операция	A=[1,2,3;4,5,6;7,8,9]	b=[22;25;35]	c=A\b
результат	A= 1 3 5 3 2 6 4 5 7	b =	c =

Задание матриц нужной размерности, состоящей из нулей – используются для начальной инициации матриц, что ускоряет расчет программ

A=zeros(2,2) – матрица размера 2 на 2, заполненная нулями:

A =

0 0

0 0

Задание матрицы можно также осуществлять в цикле, путем обращения к индексам (программы желательно писать в файлах - скриптах MATLAB):

A=zeros(2,2);

k=0;

for i=1:2 % внешний цикл перебора индекса строки

for j=1:2 % внутренний цикл перебора индекса столбца

A(i,j)=k; % задание элемента матрицы

k=k+1;

end

end

результат:

A =

 

0 1

2 3

 

Обращение к векторам, столбцам и строкам матрицы. Для того чтобы задать значения нескольких элементов вектора в одной операции используется операнд типа k ₁: k ₂, где k ₁, k ₂ – номера элементов, которые необходимо задать:

x=[1,2,3,4];

x(2:3)=[3,1]; - заменить 2-й и 3-й элементы значениями 3 и 1.

результат:

x =

1 3 1 4

 

Для того, чтобы заменить все элементы вектора используется операнд:

x=[1,2,3,4];

x(:)=[4,3,2,1];

результат:

x =

4 3 2 1

 

Аналогично для столбцов и строк матриц:

B=[1,2,3;4,5,6];

B(1:2,2)=[10;11]; - заменить 1-й и 2-й элементы 2-го столбца

результат:

B =

 

1 10 3

4 11 6

B=[1,2,3;4,5,6];

B(2,2:3)=[10,11]; - заменить 2-й и 3-й элементы 2-й строки

результат:

B =

1 2 3

4 10 11

 

N.B. Данный механизм обращения является ключевым при расчетах выражений заданных в матрично-векторном виде.

 

 

Проведение численных расчетов

Отметим ряд рекомендаций при проведении численных расчетов в программной среде Matlab:

а) код программы необходимо записывать в отдельных файлах (формат M-script), сохраняемых на жесткий диск. Запуск программы также необходимо осуществлять посредством запуска M-script файла (либо из окна Editor, либо непосредственно с диска). Крайне не рекомендуется копировать, писать или запускать код в основном рабочем окне (Command Window).

б) при любом (пере)объявлении переменной в программе, эта переменная и ее значение хранится в памяти Matlab (отображается в списке переменных – окно Workspace), в том числе и при новом «прогоне» программы, что может влиять на корректность счета. Рекомендуется в начале кода программы вставить команду clear all, которая удаляет все переменные из памяти системы.

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

Анализируемый период времени разбивается на точки с временным шагом (периодом дискретизации) , количество точек - , номера точек - , значение времени в каждой из точек - , значение функции в каждой из точек - .

Численный аналог производной непрерывной функции - дискретная функция вида: (смещенная оценка). Т.к. для ее расчета используется значение исходной функции на -м шаге, то интервал, на котором будет определенна производная будет на 1 меньше исходного (т.е. равным ). Аналогично для второй производной: с интервалом .

Схема численного интегрирования Эйлера позволяет находить численное решение систем обыкновенных дифференциальных уравнений. Например, пусть система обыкновенных линейных однородных дифференциальных уравнений (ЛОДУ) 2-го порядка задана в форме Коши:

здесь - некоторые функции (переменные состояния), - заданные параметры (коэффициенты), - производные функций.

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

Оставляя в левой части значения функций на -м шаге, получим:

Данная запись является основной расчетной формулой схемы Эйлера. Для ее расчета (т.е. расчета значений функций ) на анализируемом периоде (с учетом наличия 1-й производной – это период ) необходимо задать начальные значений функций (начальные условия).

В дальнейшем крайне рекомендуется записывать схему Эйлера в матрично-векторном виде по типу:

– из основной формулы Эйлера:

– по правилу умножения векторов:

– окончательно

где - вектор состояний, - матрица параметров, - единичная матрица.

пример построения схемы Эйлера:

clear all;

a1=1;

a2=2;

a3=3;

a4=4;

 

dt=0.01;

km=200;

 

x=zeros(2,km);

x(1,1)=10;

x(2,1)=-10;

 

A=[a1,a2;a3,a4];

I=[1,0;0,1];

 

for k=1:km-1

x(:,k+1)=(dt*A+I)*x(:,k);

end

 

Построение графиков

Базовое построение графиков. Функция plot(…) – основной инструмент построения графиков в MATLAB. Базовая форма задания функции: plot(x,y), где x – вектор значений по оси OX, y – вектор значений по оси OY (размерность этих векторов должна совпадать).

пример: построить график y=x², для

x0=0; % левый конец интервала вариации x

xm=100; % правый конец интервала вариации x

y=zeros(xm,1); % инициализация матрицы значений y

for x=x0:1:xm % цикл перебора переменной x

y(x+1)=x^2; % задание значений y

end

x=0:100; % переопределение вектора переменных x для ф-ии plot

plot(x,y(x+1)) % рисуем график

результат:

N.B. Т.к. индексы в MATLABе не могут быть меньше 1, ф-ия y в цикле определяется с индексом (k+1). Удобно индекс в цикле задавать вспомогательной переменной, начинающейся с 1 и возрастающей на каждом шаге цикла. Отметим, что задание переменной x в заголовке цикла (for x=x0:1:xm) приводит к ее переопределению, т.е. в конце цикла ее значение будет равно xm. Поэтому переменную x необходимо переопределять для использования в ф-ии plot.

Расширенные возможности по построению графиков

а) дополнительные возможности функции plot() – задание свойств (атрибутов) линий графика

plot(x,y,'Имя атрибута №1', Значение атрибута №1,'Имя атрибута №2', 'Значение атрибута №2',...)

N.B. Иногда имя и значение атрибута совмещаются (см. пример – тип линии). Если значение атрибута не является числом, то оно должно быть заключено в одинарные кавычки ''. Иногда значение атрибута может быть вектором и т.д.

тип свойства	имя свойства	значение (пример)	эффект
цвет линии	'color'	[1,0.5,0.2]	задает цвет в формате RGB
тип линии	'--'	пунктирная линия
цвет и тип линии (совместно)	'-k'	сплошная линия черного цвета
толщина линии	'LineWidth'		задает толщину линии

б)функция figure('Имя атрибута №1’, Значение атрибута №1,...) – задание свойств (атрибутов) окна графика

тип свойства	имя свойства	значение (пример)	эффект
создание нового окна графиков	figure();	создает новое окно для отображения графиков
задание типа графика при создании	figure('name','x^2');	создание нового окна с названием «x^2»

 

в) функция legend('имя 1,'имя 2',...) – создание легенды графика. Порядок подписей (имен) в легенде должен соответствовать порядку линий на графике

тип свойства	имя свойства	значение (пример)	эффект
создание легенды	legend('x^2','exp(x)')	создает новую легенду для графика

 

г) функция set(gca, 'Имя атрибута №1', Значение атрибута №1, …) – универсальная функция изменения свойств объектов в MATLAB, в данном случае – свойств осей окна графика (текущих активных осей текущего активного окна графика - gca)

тип свойства	имя свойства	значение (пример)	эффект
изменение пределов на оси y графика	'YLim'	[0 100]	задает пределы по оси y от 0 до 100

 

д) функция xlabel('подпись, 'Имя атрибута №1', Значение атрибута №1, …) – функция отображения подписей по осям окна графика

тип свойства	имя свойства	значение (пример)	эффект
собственно подпись к оси	'x'	подпись к оси – x
тип шрифта	'fontname'	'arial'	задает тип шрифта
размер шрифта	'fontsize'		задает размер шрифта

N.B. При построении более чем одного графика на одних осях необходимо после первого заданного графика (первой записанной функции plot()) вставить команду hold on (иначе будет отображаться только последний заданный график)

N.B.2 Для удобства чтения графиков на них рекомендуется наносить сетку – команда grid on

N.B.3 Полный список атрибутов к любому из объектов нужно искать в справке.

пример: построить графики y = x² и y = e^x для

x0=0;

xm=5;

y=zeros(xm,2);

for x=x0:1:xm

y(x+1,1)=x^2;

y(x+1,2)=exp(x);

end

x=0:5;

plot(x,y(x+1,1),'-k','LineWidth',2);

hold on;

grid on;

plot(x,y(x+1,2),'--','LineWidth',2,'color',[1,0.5,0.2]);

legend('x^2','exp(x)');

set(gca,'YLim',[0 100]);

xlabel('x','fontname','arial');

ylabel('y','fontsize',23)

результат: