Матлаб. К. Ю. Петрова введение в matlab учебное пособие
 Скачать 2.57 Mb.
|
|
sin, cos должен указываться в радианах, если задавать его в градусах (degree), то надо использовать команды sind, cosd (впервые они появились в MATLAB 7). Например, sin(pi/6) и sind(30) дадут один и тот же результат ans=0.5000. Аргументом каждой из функций может быть число или вектор (набор чисел, массив). Например, набрав sind([0, 30, 90]) получим ans=0 0.5000 1.0000. Для оперативного получения справок об этих и других командах используется команда help. Например, набрав в MATLAB 7 >> help sind, получим справку: SIND Sine of argument in degrees. SIND(X) is the sine of the elements of X, expressed in degrees. For integers n, sind(n*180) is exactly zero, whereas sin(n*pi) reflects the accuracy of the floating point value of pi. See also asind, sin. 1.2.2. Ввод числовых данных Перечислим несколько простых команд для ввода числовых данных в виде векторов и матриц. Самый простой способ формирования векторов и матриц в MATLAB заключается в непосредственном вводе их элементов с клавиатуры. Например, набирая на клавиатуре данные Х = [1 -2 3 8 5 6], получаем одномерный массив (вектор-строку) Х из шести элементов. Формирование вектора-строки из равноотстоящих значений аргумента выполняется с помощью команды x = x0:h:xn . По умолчанию шаг h принимается равным 1. Например, команда x=0:10 дает целые числа от 0 до 10, а x=0:0.1:10 задает набор значений аргумента от нуля до 10 с шагом 0.1. Двумерные массивы задаются в виде матриц, при этом строки разделяются символом «точка с запятой». Элементы одной и той же строки могут разделяться как пробелами, так и запятыми: >> a=[1 2 3; 4 5 6; 7 8 9 a = 1 2 3 4 5 6 7 8 9 >> a=[1 ,2, 3; 4, 5, 6; 7, 8, 9] a = 1 2 3 4 5 6 7 8 9 Для доступа к элементам массива используются круглые скобки: >> a(1,1) ans = 1 >> a(3,3) ans = 9 >> b=[1 2 3 4 5]; b(4) ans = 4 Для получения строки или столбца матрицы используется символ «двоеточие»: >> a(:,1) ans = 1 4 7 >> a(2,:) ans = 4 5 6 >> b(1:3) ans = 1 2 3 >> b(3:end) ans = 3 4 5 10 В MATLAB имеется ряд команд, облегчающих формирование векторов и матриц (табл.2). Таблица 2 linspace zeros eye logspace ones diag Команда linspace используется для получения набора равноотстоящих значений аргумента. Например, х=linspace(xmin,xmax) создаст массив из 100 значений аргумента между точками хmin и xmax . Если требуется иметь другое число точек, например, N точек на интервале от 0 до 10, следует записать x = linspace(0,10,N). Аналогичные модификации имеет команда logspaсe, обеспечивающая логарифмическое расположение точек массива. Для создания некоторых распространенных матриц имеются команды zeros, ones, eye и diag. Команды zeros и ones, создают матрицы, заполненные нулями и единицами. Так, по командaм A=zeros (3) и B=ones (3) будут сформированы матрицы A= 0 0 0 B = 1 1 1 0 0 0 1 1 1 0 0 0 1 1 1 , а по командам C = zeros(1,4) и D = ones(1,4) – векторы-строки C = 0 0 0 0 и D = 1 1 1 1. Команда eye предназначена для создания (она часто обозначается буквой I, название которой произносится как слово «eye») или ее фрагмента. >> E=eye(3) E = 1 0 0 0 1 0 0 0 1 >> G=eye(2,7) G = 1 0 0 0 0 0 0 0 1 0 0 0 0 0 Команда eye(size(A)) даст единичную матрицу того же размера, что и матрица А. Команда diag имеет два различных назначения. Если ее аргументом является матрица, то результатом будет вектор, содержащий элементы, стоящие на диагонали. Если же аргументом является вектор, то команда diag создает матрицу с заданной диагональю и остальными нулевыми элементами. Например, применяя команду diag к приведенной выше матрице В, получим строку [1 1 1], а повторное применение команды diag даст единичную матрицу Е. 1.3. Построение графиков Основное средство для построения графиков в MATLAB – это команда plot и различные ее модификации. Она может вызываться с одним или несколькими входными аргументами. Стандартный вариант ее вызова – это plot(x,y), где x и y – два массива чисел, содержащие абсциссы и ординаты точек графика функции y = f(x). Выше был приведен пример построения графика синусоиды, аналогично строятся графики любых других функций. При этом вычерчивание осей и выбор масштабов по ним производится автоматически. В случае, если вызов команды plot производится с одним аргументом в формате plot(y), координатами x служат индексы массива y Для того чтобы снабдить рисунок координатной сеткой, используется команда grid. Вызов ее без параметров осуществляет переключение режимов «с сеткой»/«без сетки», а задание grid on и grid off явно указывает, следует включить сетку или отключить. Иногда на одном графике требуется нарисовать несколько кривых. В этом случае в команде plot указывают несколько пар аргументов (по числу функций) plot(х1, у1, х2, у2, ..., хn, уn), где х1, у1; 11 х2, у2 и т.д. – пары векторов. Каждой паре х, у будет соответствовать свой график, при этом они могут быть заданы векторами разной длины. Пример. Пусть требуется построить графики затухающих колебаний t e t x t s i n ) ( 2 , 0 , t e t y t cos ) ( 2 , 0 , причем аргумент t изменяется от 0 до 10 с шагом 0,1. Это делается с помощью следующей группы команд: >>t=0:.1:10; x=exp( –.2*t).*sin(t); y=exp(–.2*t).*cos(t); plot(t, x, t, y), grid. Результат показан на рис. 1.3. Использование точки перед знаком * (умножение) при вычислении переменных x, y, указывает на поэлементное перемножение массивов чисел (каждая из функций sin t, cos t, e -0,2t , представлена вектором из 101 точек). Добавляя команду plot(x,y), grid, получим график логарифмической спирали, показанный на рис. 1.4. Рис. 1.3 -0.5 0 0.5 1 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1 Рис.1.4 После того как график выведен на экран, его можно озаглавить, обозначить оси, сделать текстовую разметку, для чего используются команды title, xlabel, ylabel, text. Например, чтобы нанести обозначения осей на последний график, надо набрать xlabel('x'), ylabel('y') В команде plot в одиночных кавычках можно использовать дополнительный аргумент, указывающий тип символов, используемых для построения графика. Так, plot(X,Y,'x') вычерчивает точечный график, используя символы x (крестики), тогда как plot(X1,Y1,':',X2,Y2,'+') использует символ двоеточия для первой кривой и символ + для втоpой. Цвет линий также может задаваться пользователем. Например, команды plot (X,Y,'r') и plot (X,Y,'+g') используют красную линию для получения первого графика и зеленые + метки для второго. Справку о возможных вариантах типов линий, точек и цветов можно получить, набрав help plot. Команда plot строит графики на плоскости. MATLAB позволяет также наглядно изображать линии и поверхности в трехмерном пространстве. Для изображения линий в пространстве используется команда plot3. Получим, например, график винтовой линии, которая задается уравнениями x=sint; y=cost; z=t. Возьмем диапазон изменения параметра 10 0 t с шагом 0,02 : >>t = 0:pi/50:10*pi; plot3(sin(t),cos(t),t); Результат показан на рис. 1.5. 0 2 4 6 8 10 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1 12 Три стандартные поверхности – сфера, эллипсоид и цилиндр – строятся с помощью команд sphere, ellipsoid, cylinder соответственно. Результат выполнения первой из них показан на рис. 1.6. -1 0 1 -1 0 1 0 10 20 30 40 -1 0 1 -1 0 1 -1 -0.5 0 0.5 1 Рис. 1.5 Рис. 1.6 Дополнительные сведения о графических возможностях MATLAB можно найти в разд. 1.7 и 4.1. 1.4. Матричные операции MATLAB имеет обширный арсенал матричных операций. К простейшим из них относятся сложение и умножение, вычисление ранга и определителя, а также обращение матрицы. Элементарные операции над матрицами перечислены в табл. 3 Таблица 3 A±B A' trace inv A*B det rank pinv Для сложения и вычитания матриц одинакового размера используются знаки + и –, например, С = А + В. Умножение матриц обозначается звездочкой: С=А * В. Оно допускается, если число строк матрицы А равно числу столбцов матрицы В. При этом в общем случае A B B A * * Напомним, что матричное умножение вычисляется по известному правилу «строка на столбец». В частности, произведение матриц второго порядка находится по формуле 22 22 12 22 21 22 11 21 22 12 12 11 21 12 11 11 22 21 12 11 22 21 12 11 b a b a b a b a b a b a b a b a b b b b a a a a Для MATLAB такое умножение – элементарная операция. Приведем простой пример: >>A=[1 2; 3 4] >>B=[5 6; 7 8] >>C=A*B >>CT=B'*A' 1 2 3 4 5 6 7 8 19 22 43 50 19 43 22 50 13 Простейшей матричной операцией является транспонирование (строки заменяются столбцами). Например, результатом транспонирования вектора-строки будет вектор-столбец. В MATLAB транспонирование обозначается штрихом. В последней колонке приведена транспонированная матрица С, она равна произведению транспонированных матриц A и B, взятых в обратном порядке. По команде trace(A) вычисляется след матрицы А, т.е. сумма ее диагональных элементов. Полезно помнить, что он равняется также сумме ее собственных чисел. По команде rank находится ранг матрицы. Он определяется максимальным размером ее минора с ненулевым определителем и одновременно указывает на число линейно независимых строк (столбцов) матрицы. Следующая матричная операция – вычисление определителя, осуществляется командой det. Как известно, определитель матрицы d c b a A размера 2x2 равен bc ad A det . Формула для определителя третьего порядка имеет вид: det 3 1 2 2 3 1 1 2 3 2 1 3 1 3 2 3 2 1 3 2 1 3 2 1 3 2 1 c b a c b a c b a c b a c b a c b a c c c b b b a a a Обращение квадратной матрицы А производится по команде inv(A). Она является основной при решении системы линейных алгебраических уравнений. Напомним формулу для вычисления обратной матрицы , det 1 * 1 A A A где А * – присоединенная матрица, составленная из алгебраических дополнений A ij матрицы А. В соответствии с этой формулой процедура вычисления обратной матрицы содержит 2 шага. Шаг 1. Каждый элемент a ij матрицы А заменяется его алгебраическим дополнением A ij , т.е. определителем матрицы, получаемой вычеркиванием соответствующей строки и столбца. Если сумма индексов i+j нечетна, определитель берется с минусом. Шаг 2. Полученная матрица транспонируется и делится на определитель исходной матрицы. В частности, для матрицы второго порядка получаем a c b d bc ad d c b a 1 1 Пример. Табл. 4 содержит примеры выполнения упомянутых операций. Таблица 4 >> A=[1 2 3; 0 1 4; 0 0 1] >> A' >> det(A) >> inv(A) >> inv(A') 1 2 3 0 1 4 0 0 1 1 0 0 2 1 0 3 4 1 1 1 -2 5 0 1 -4 0 0 1 1 0 0 -2 1 0 5 -4 1 В первой строке приведены команды, набираемые в рабочем окне, во второй – ответы, которые даст MATLAB. В частности видим, что при транспонировании матрицы обратная матрица тоже транспонируется. 14 Типичная задача линейной алгебры – решение системы линейных уравнений. На матричном языке она сводится к тому, чтобы найти вектор X , удовлетворяющий равенству B AX , где матрица A и вектор-столбец B заданы. Решение этого уравнения имеет вид B A X 1 Для получения такого решения в пакете MATLAB достаточно набрать >>X=inv(A)*B либо >>X=(A^-1)*B. Разумеется, для выполнения этих операций матрица A должна быть невырожденной. В случае если система уравнений B AX недоопределена или переопределена, решение получают с помощью команды pinv, выполняющей псевдообращение матрицы А. Можно также использовать знак деления – косую черту (slash и backslash). Пример. Требуется решить систему трех уравнений с тремя неизвестными 5 2 2 2 2 2 2 z y x z y x z y x В данном случае матрицы А и В имеют вид 5 2 2 , 2 2 1 1 2 2 1 1 1 B A Вводя их в MATLAB и набирая >>X=inv(A)*B, получим: >>A=[1 1 1; 2 2 1; 1 2 2] >>B=[2; 2; 5] >>X=inv(A)*B A=1 1 1 2 2 1 1 2 2 B=2 2 5 X=-1 1 2 Следовательно, решение системы имеет вид x=-1, y=1, z=2. 1.5. Работа с полиномами В состав MATLAB входит ряд команд, позволяющих выполнять различные операции с полиномами от одной переменной, включая поиск корней, умножение и деление полиномов, построение полинома, проходящего через заданные точки и др. Полином описывается строкой своих коэффициентов в порядке от старшего к младшему. Так, полином 5 2 3 x x будет представлен вектором 5 2 0 1 Перечень основных команд MATLAB для работы с полиномами приведен в табл. 5. Таблица 5 roots conv polyval poly deconv residue Дадим краткие пояснения к ней. Команда roots предназначена для отыскания корней полинома. Например, чтобы решить квадратное уравнение 0 6 5 2 x x следует набрать r=roots([1 5 6]) , результатом будут значения корней r=[-2; -3]. Функция poly выполняет обратную операцию – строит полином по заданным корням. Так p=poly([-2 -3]) даст p=[1 5 6]. Если в качестве входного аргумента функции poly фигурирует квадратная матрица, то результатом будет ее характеристический полином. 15 Следует иметь в виду, что вычисление корней полиномов высокой степени с помощью команды roots может сопровождаться заметными погрешностями, поэтому пользоваться ей следует с осторожностью. Пример ( демонстрация ошибок при вычислении корней). Все корни полинома 8 ) 1 ( x вещественны и равны – 1. Сформируем этот полином командой poly и найдем его корни командой roots. >>p=poly(-ones(1, 8)), r=roots(p), p = 1 8 28 56 70 56 28 8 1 r = -1.0203; -1.0142 0.0144i; -0.9998 0.0201i; -0.9858 0.0140i; -0.9801 В результате получили комплексные числа, причем их модуль отличается от единицы на 2%. Если учесть, что MATLAB работает с 32-разрядной сеткой, это очень большая погрешность. В наглядной форме ошибки вычисления корней показаны на рис. 1.7. Он получен путем пятикратного применения пары команд roots/poly, когда по корням восстанавливался полином, снова искались его корни и т.д. Кружок в центре рис. 1.7 соответствует истинным корням полинома 8 ) 1 ( x , а звездочками помечены значения корней, вычисленные в процессе итераций. % Корневой тест n=8; r=-ones(1,n); clg, hold on for i=1:5, p=poly(r), r=roots(p), plot(r,'*','LineWidth',2) end plot(-1,0,'o','LineWidth',2),grid on; hold off Избежать подобных ошибок можно, переходя к символьным вычислениям (для этого необходимо наличие тулбокса SYMBOLIC). >> r=-ones(1,8); p=poly(r) % численное представление полинома p = 1 8 28 56 70 56 28 8 1 >> P=poly2sym(p), % символьное представление полинома P =x^8+8*x^7+28*x^6+56*x^5+70*x^4+56*x^3+28*x^2+8*x+1 >>R=solve(P); % решение уравнения 0 ) ( x P R'=[ -1, -1, -1, -1, -1, -1, -1, -1] % найденные корни Как видим, символьный решатель уравнений solve дал точный ответ. К сожалению, его возможности ограничены только уравнениями, допускающими аналитическое решение. Для того чтобы построить график полинома, надо предварительно вычислить его значения в точках заданного интервала. Для этой цели служит функция polyval (сокращение от polynomial value). Например, чтобы построить график полинома 6 5 2 x x y на интервале 5 5 x следует набрать: >> x= -5:0.1:5; p=[1 5 6]; y=polyval(p, x); plot(x, y), grid. В результате будет получен график, который пересекает ось абсцисс в точках х 1 = -3, х 2 = -2 (это найденные выше корни полинома). -1.04 -1.02 -1 -0.98 -0.96 -0.025 -0.02 -0.015 -0.01 -0.005 0 0.005 0.01 0.015 0.02 0.025 (x -1) 8 Рис. 1.7 |