«Киев – НТУУ КПИ - 2003 УДК 681.3.06(075.8) ББК 32.973.26-018.2 Я73 Л17 Лазарев Юрий Федорович Л17 Начала программирования в среде MatLAB: Учебное пособие. - К.: НТУУ КПИ, 2003. - 424 с. Изложены основные особенности ...»

Юрий ЛАЗАРЕВ

_

Начала программирования

в среде MatLAB

Учебное пособие

для студентов высших учебных заведений

Киев – НТУУ "КПИ" - 2003

УДК 681.3.06(075.8)

ББК 32.973.26-018.2 Я73

Л17

Лазарев Юрий Федорович

Л17 Начала программирования в среде MatLAB: Учебное пособие. - К.:

НТУУ "КПИ", 2003. - 424 с.

Изложены основные особенности проведения вычислений в среде MatLAB как в режиме калькулятора, так и в программном режиме. Ознакомление с системой рассчитано на начинающего. Приведены сведения об основных командах, операторах, функциях и процедурах MatLAB. Изложение ведется таким образом, чтобы пользователь мог сразу применить полученные знания для проведения вычислений. Пособие содержит много примеров, которые поясняют и иллюстрируют работу по использованию процедур. Рассмотрена работа с некоторыми наиболее важными для инженеров пакетами прикладных программ MatLAB. (Signal Toolbox, Control и SimuLink).

Для студентов высших технических учебных заведений. Может быть полезно научным работникам и инженерам для начального ознакомления с системой MatLAB и приобретения навыков работы с ней.

Табл. 8. Илл. 283. Библиогр. 18 назв.

© Ю. Ф. Лазарев, Содержание Содержание Предисловие Вступление 1. MatLAB как научный калькулятор 1.1. Командное окно 1.2. Операции с числами 1.2.1. Ввод действительных чисел 1.2.2. Простейшие арифметические действия 1.2.3. Ввод комплексных чисел 1.2.4. Элементарные математические функции 1.2.5. Специальные математические функции 1.2.6. Элементарные действия с комплексными числами 1.2.7. Функции комплексного аргумента 1.2.8. Задания 1.2.9. Вопросы 1.3. Простейшие операции с векторами и матрицами 1.3.1. Ввод векторов и матриц 1.3.2. Формирование векторов и матриц 1.3.3. Извлечение и вставка частей матриц 1.3.4. Действия над векторами 1.3.5. Поэлементное преобразование матриц 1.3.6. Матричные действия над матрицами 1.3.7. Матричные функции 1.3.8. Задания 1.3.9. Вопросы 1.4. Функции прикладной численной математики 1.4.1. Операции с полиномами 1.4.2. Обработка данных измерений 1.4.3. Функции линейной алгебры 1.4.4. Аппроксимация и интерполяция данных 1.4.5. Векторная фильтрация и спектральный анализ 1.4.6. Задания 1.4.7. Вопросы 1.5. Построение простейших графиков 1.5.1. Процедура plot 1.5.2. Специальные графики 1.5.3. Дополнительные функции графического окна 1.5.4. Вывод графиков в печать 1.5.5. Задания 1.5.6. Вопросы 1.6. Операторы управления вычислительным процессом 1.6.1. Оператор условного перехода 1.6.2. Оператор переключения 1.6.3. Операторы цикла 1.6.4. Задания 1.6.5. Вопросы 2. Программирование в среде MatLAB 2.1. Функции функций 2.2. Создание М-файлов 2.2.1. Особенности создания М-файлов 2.2.2. Основные особенности оформления М-файлов 2.3. Создание простейших файлов-функций (процедур) 2.3.1. Общие требования к построению 2.3.2. Типовое оформление процедуры-функции 2.3.3. Задания 2.3.4. Вопросы 2.4. Создание Script-файлов 2.4.2. Ввод и вывод информации в диалоговом режиме 2.5.1. Общие требования к представлению графической информации 3.3.5. Некоторые особенности использования системы MatLAB 4.3.3. Создание процедуры символьного представления polynom-объекта 5. Цифровая обработка сигналов (пакет Signal Processing Toolbox) 5.2. Общие средства фильтрации. Формирование случайных процессов 5.3. Процедуры спектрального (частотного) и статистического анализа процессов 6. Исследование линейных стационарных систем (пакет Control Toolbox) 7. Моделирование нелинейных систем (пакет SimuLINK) 7.3.1. Моделирование поведения физического маятника 7.3.2. Моделирование поведения гироскопа в кардановом подвесе 7.4.1. Принципы функционирования блоков системы SimuLink 7.4.3. Передача данных между средой MatLab и S-моделью 7.4.4. Запуск процесса моделирования S-модели из среды MatLab 7.4.5. Образования S-блоков путем использования программ 7.5.3. Моделирование процесса ориентации космического аппарата Указатель операторов, команд, функций и функциональных блоков MatLAB Предисловие В последние годы в университетских и инженерно-технических кругах мира наблюдается интенсивное распространение новой компьютерной системы осуществления математических расчетов - системы MatLAB. В чем причина такой популярности этой системы?

Главные преимущества "языка технических вычислений" MatLAB, которые выгодно отличают его среди других существующих ныне математических систем и пакетов, состоят в следующем:

система MatLAB специально создана для проведения именно инженерных расчетов: математический аппарат, который используется в ней, предельно приближен к современному математическому аппарату инженера и ученого и опирается на вычисления с матрицами, векторами и комплексными числами; графическое представление функциональных зависимостей здесь организовано в форме, которую требует именно инженерная документация;

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

в отличие от большинства математических систем, MatLAB является открытой системой; это означает, что практически все процедуры и функции MatLAB доступны не только для использования, но и для корректировки и модифицирования; MatLAB - система, которая может расширяться пользователем по его желанию созданными им программами и процедурами (подпрограммами); ее легко приспособить к решению нужных очень удобной является возможность использовать практически все вычислительные возможности системы в режиме чрезвычайно мощного научного калькулятора; в то же время можно составлять собственные отдельные программы с целью многоразового их использования для исследований; это делает MatLAB незаменимым средством проведения научных расчетных исследований;

последние версии MatLAB позволяют легко интегрировать ее с текстовым редактором Word, что делает возможным использование при создании текстовых документов вычислительных и графических возможностей MatLAB, например, оформлять инженерные и научные отчеты и статьи с включением в них сложных расчетов и выводом графиков в текст.

Возможности системы огромны, а по скорости выполнения задач она опережает многие другие подобные системы. Все эти особенности делают систему MatLAB весьма привлекательной для использования в учебном процессе высших учебных заведений.

Эта книга является вторым, переработанным и существенно дополненным изданием учебного пособия "Початки програмування в середовищі MatLAB" [11], содержит, в основном, описание MatLAB версии 5.3 и в него добавлены материалы по практическому овладению процедурами пакетов CONTROL (анализа и синтеза линейных систем автоматического управления), SIGNAL (цифровой обработки сигналов), SIMULINK (интерактивного моделирования динамических систем) и некоторых новых важных возможностей MatLAB.

Пособие состоит из семи глав.

В первой главе читатель знакомится с возможностями системы в режиме научного калькулятора. Здесь помещены сведения об основных операторах, командах, функциях и процедурах системы.

Во второй главе описаны правила и примеры составления программ на языке MatLAB. Кроме того, в ней представлены некоторые дополнительные процедуры, которые помогают рационально организовать вычислительный процесс.

Третья глава содержит перечень некоторых процедур и команд общего назначения, которые связывают систему MatLAB с операционной системой компьютера. Здесь же описано использование редактора Word с системой MatLAB.

Важной частью MatLAB, которая позволяет приспосабливать систему к задачам пользователя, является возможность образования новых классов вычислительных объектов. С понятием классов вычислительных объектов в MatLAB и правилами создания новых классов пользователь ознакомится в четвертой главе.

В пятой главе сосредоточены сведения об особенностях использования процедур цифровой обработки сигналов пакета SIGNAL.

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

Седьмая глава знакомит с пакетом SimuLink интерактивного (визуального моделирования динамических систем во временной области).

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

В пособии использованы материалы из изданий, указанных в списке литературы.

Введение Система MatLAB создана фирмой MathWork Inc. (США, г. Нейтик, штат Массачусетс). Хотя впервые эта система начала использоваться в конце 70-х лет, расцвет ее применения начался в конце 80-х, в особенности после появления на рынке версии 4.0. Последние версии MatLAB, - это чрезвычайно развитые системы, которые содержат огромную совокупность процедур и функций, необходимых инженеру и научному работнику для осуществления сложных численных расчетов, моделирования поведения технических и физических систем, оформления результатов этих расчетов в наглядном виде.

Система MatLAB (сокращение от MATrix LABoratory - матричная лаборатория) представляет собой интерактивную компьютерную систему для выполнения инженерных и научных расчетов, ориентированную на работу с массивами данных. Система предполагает возможность обращения к программам, которые написаны на языках FORTRAN, C и C++.

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

Основной объект системы MatLAB - прямоугольный числовой массив (матрица), который допускает комплексные элементы. Использование матриц не требует явного указания их размеров. Система позволяет решать многие вычислительные задачи за значительно меньшее время, чем то, которое необходимо для написания соответствующих программ на языках FORTRAN, BASIC и C.

Система MatLAB выполняет операции с векторами и матрицами даже в режиме непосредственных вычислений без какого-либо программирования. Ею можно пользоваться как мощнейшим калькулятором, в котором наряду с обычными арифметическими и алгебраическими действиями могут использоваться такие сложные операции, как обращение матрицы, вычисление ее собственных значений и векторов, решение систем линейных алгебраических уравнений и многое другое. Тем не менее, характерная основная особенность системы – ее "открытость", то есть легкость ее модификации и адаптации к конкретным задачам пользователя. Пользователь может ввести в систему любую новую команду, оператор или функцию и пользоваться потом ими так же просто, как и встроенными операторами и функциями.

В базовый набор слов системы входят: спецзнаки; знаки арифметических и логических операций; арифметические, тригонометрические и некоторые специальные математические функции; функции быстрого преобразования Фурье и фильтрации; векторные и матричные функции; средства для работы с комплексными числами; операторы построения графиков в декартовой и полярной системах координат, трехмерных поверхностей и т.п. То есть MatLAB предоставляет пользователю большой набор готовых средств (более половины из них - внешние расширения в виде m-файлов).

Система MatLAB имеет собственный язык программирования, который напоминает BASIC. Запись программ в системе является традиционной и потому обычной для большинства пользователей персональных компьютеров. И вдобавок система дает возможность редактировать программы при помощи любого привычного для пользователя текстового редактора.

MatLAB имеет широкие возможности для работы с сигналами, для расчета и проектирования аналоговых и цифровых фильтров, для построения их частотных, импульсных и переходных характеристик. В наличии и средства для спектрального анализа и синтеза, в частности, для реализации прямого и обратного преобразования Фурье. Благодаря этому система довольно удобна для проектирования электронных устройств.

С системой MatLAB поставляются свыше ста m-файлов, которые содержат демонстрационные примеры и определения новых операторов и функций.

Эта библиотека, все файлы которой подробно прокомментированы, - подлинная сокровищница прекрасных примеров программирования на языке системы. Изучение этих примеров и возможность работы в режиме непосредственных вычислений значительно облегчают знакомство с системой серьезных пользователей, заинтересованных в использовании математических расчетов.

Робота в среде MatLab может осуществляться в двух режимах:

в режиме калькулятора, когда вычисления осуществляются сразу после набора очередного оператора или команды MatLab; при этом значение результатов вычисления могут присваиваться некоторым переменным, или результаты получаются непосредственно, без присваивания (как в обычных калькуляторах);

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

В обоих режимах пользователю доступны практически все вычислительные возможности системы, в том числе по выводу информации в графической форме. Программный режим позволяет сохранять разработанные вычислительные алгоритмы и, таким образом, повторять вычисления при других входных данных.

MatLAB имеет черты разных известных языков программирования высокого уровня.

С языком BASIC ее роднит то, что она представляет собой интерпретатор (то есть компилирует и выполняет программу пооператорно, не образовывая отдельного исполняемого файла), незначительное количество операторов и отсутствие необходимости в объявлении типов и размеров переменных, то есть все то, что делает язык программирования весьма удобным в пользовании.

От языка Pascal система MatLAB позаимствовала объектноориентированную направленность, то есть такое построение языка, которое обеспечивает образование новых типов вычислительных объектов по желанию пользователя на основе типов объектов, уже существующих в языке. Эти новые типы объектов (в MatLAB они называются классами) могут иметь собственные процедуры их преобразования (они образуют методы этого класса). Весьма удобно то, что новые процедуры могут быть вызваны с помощью обычных знаков арифметических операций и некоторых специальных знаков, которые применяются в математике.

Принципы сохранения значений переменных в MatLAB более всего приближаются к тем, которые присущи языку FORTRAN, а именно: все переменные являются локальными, то есть действуют лишь в границах той программной единицы (процедуры, функции или главной, управляющей программы), где им присвоены некоторые конкретные значения. При переходе к выполнению другой программной единицы, значения переменных предыдущей программной единицы или совсем теряются (в случае, если выполненная программная единица является процедурой или функцией), или становится недосягаемыми (если выполненная программа - управляющая). В отличие от языков BASIC и Pascal, в языке MatLAB нет глобальных переменных, действие которых распространялась бы на все программные единицы. Но и здесь в MatLAB есть особенность, которая отличает ее от других языков. В отличие от интерпретатора BASIC интерпретатор MatLAB позволяет в одном и том же сеансе работы выполнять несколько самостоятельных программ, причем все переменные, используемые в этих программах, являются общими для этих программ и образуют общее рабочее пространство (Work Space). Это позволяет более рационально организовывать сложные (громоздкие) вычисления по принципу, который напоминает оверлейные структуры.

Эти особенности MatLAB делают ее весьма гибкой и удобной в пользовании вычислительной системой.

1. MatLab как научный калькулятор 1.1. Командное окно После вызова MatLAB из среды Windows на экране возникает изображение так называемого "командного окна" среды MatLAB (рис. 1.1).

Это окно является основным в MatLAB. В нем появляются символы команд, которые набираются пользователем на клавиатуре дисплея, отображаются результаты выполнения этих команд, текст исполняемой программы и информация об ошибках выполнения программы, распознанных системой.

Признаком того, что MatLAB готова к восприятию и выполнению очередной команды, является возникновение в последней строке текстового поля окна знака приглашения '»', после которого расположена мигающая вертикальная черта.

В верхней части окна (под заголовком) размещена строка меню, в которой находятся меню File, Edit, View, Windows, Help. Чтобы открыть какое-либо меню, следует установить на нем указатель мыши и нажать ее левую кнопку. Подробнее функции команд меню описаны далее, в главе 3 "Интерфейс MatLab и команды общего назначения. Написание М-книг".

Здесь отметим лишь, что для выхода из среды MatLAB достаточно открыть меню File и выбрать в нем команду Exit MATLAB, или просто закрыть командное окно, нажав левую клавишу мыши, когда курсор мыши установлен на изображении верхней крайней правой кнопки этого окна (с обозначением косого крестика).

1.2. Операции с числами 1.2.1. Ввод действительных чисел Ввод чисел с клавиатуры осуществляется по общим правилам, принятым для языков программирования высокого уровня:

для отделения дробной части мантиссы числа используется десятичная точка (вместо запятой при обычной записи);

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

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

Если, например, ввести в командном окне MatLAB строку 1. 20357651е -17, то после нажатия клавиши в этом окне появится запись:

Следует отметить, что результат выводится в виде (формате), который определяется предварительно установленным форматом представления чисел. Этот формат может быть установлен с помощью команды Preferences меню File (рис.

1.3). После ее вызова на экране появится одноименное окно (рис. 1.4). Один из участков этого окна имеет название Numeric Format. Он предназначен для установки и изменения формата представления чисел, которые выводятся в командное окно в процессе расчетов. Предусмотрены такие форматы:

Short (default) - краткая запись (применяется по умолчанию);

Hex - запись в виде шестнадцатиричного числа;

Plus - записывается только знак числа;

Short E - краткая запись в формате с плавающей запятой;

Long E - длинная запись в формате с плавающей запятой;

Short G - вторая форма краткой записи в формате с плавающей запятой;

Long G - вторая форма длинной записи в формате с плавающей запятой;

Rational - запись в виде рациональной дроби.

Избирая с помощью мыши нужный вид представления чисел, можно обеспечить в дальнейшем выведение чисел в командное окно именно в этой форме.

Как видно из рис. 1.2, число, которое выведено на экран, не совпадает с введенным. Это обусловлено тем, что установленный по умолчанию формат представления чисел (Short) не позволяет вывести більше 6 значащих цифр. На самом рами. Например, если избрать мышью селекторную кнопку Long E (т. е. установить указанный формат представления чисел), то, повторяя те же действия, получим:

где уже все цифры отображены верно (рис. 1.5).

Следует помнить:

- введенное число и результаты всех вычислений в системе MatLAB сохраняются в памяти ПК с относительной погрешностью около 2.10-16 (т. е. с точными значениями в 15 десятичных разрядах);

- диапазон представления модуля действительных чисел лежит в диапазоне между 10 308 и 10 +308.

1.2.2. Простейшие арифметические действия В арифметических выражениях языка MatLAB используются следующие знаки арифметических операций:

Использование MatLAB в режиме калькулятора может происходить путем простой записи в командную строку последовательности арифметических действий с числами, то есть обычного арифметического выражения, например:

(4.5)2*7.23 - 3.14*10.4.

Если после ввода с клавиатуры этой последовательности нажать клавишу, в командном окне возникнет результат выполнения в виде, представленном на рис. 1.6, т. е. на экран под именем системной переменной ans выводится результат действия последнего выполненного оператора.

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

- если запись оператора не заканчивается символом ';', результат действия этого оператора сразу же выводится в командное окно;

- если оператор заканчивается символом ';', результат його действия не отображается в командном окне;

- если оператор не содержит знака присваивания (=), т. е. является просто записью некоторой последовательности действий над числами и переменными, значение результата присваивается специальной системной - полученное значение переменной ans можно использовать в следующих операторах вычислений, применяя это имя ans; при этом следует помнить, что значение системной переменной ans изменяется после действия очередного оператора без знака присваивания;

- в общем случае форма представления результата в командном окне имеет вид:

Пример. Пусть нужно вычислить выражение (25+17)*7. Это можно сделать таким образом. Сначала набираем последовательность 25+17 и нажимаем.

Получаем на экране результат в виде ans = 42. Теперь записываем последовательность ans*7 и нажимаем. Получаем ans = 294 (рис. 1.7). Чтобы предотвратить выведение промежуточного результата действия 25+17, достаточно после записи этой последовательности добавить символ ' ; '. Тогда будем иметь результаты в виде, представленном на рис. 1.8.

Применяя MatLAB как калькулятор, можно использовать имена переменных для записи промежуточных результатов в память ПК. Для этого служит операция присваивания, которая вводится знаком равенства ' = ' в соответствия со схемой:

Имя переменной может содержать до 30 символов и должно не совпадать с именами функций, процедур системы и системных переменных. При этом система различает большие и маленькие буквы в переменных. Так, имена 'amenu', 'Amenu', 'aMenu' в MatLAB обозначают разные переменные.

Выражение справа от знака присваивания может быть просто числом, арифметическим выражением, строкой символов (тогда эти символы нужно заключить в апострофы) или символьным выражением. Если выражение не заканчивается символом ' ; ', после нажатия клавиши в командном окне возникнет результат выполнения в виде :

Например, если ввести в командное окно строку 'x = 25 + 17', на экране появится запись (рис. 1.9) :

Система MatLAB имеет несколько имен переменных, которые используются самой системой и входят в состав зарезервированных:

i, j - мнимая единица (корень квадратный из -1); pi - число (сохраняется в виде 3.141592653589793); inf - обозначение машинной бесконечности; Na - обозначение неопределенного результата (например, типа 0/0 или inf/inf); eps - погрешность операций над числами с плавающей запятой; ans - результат последней операции без знака присваивания; realmax и realmin – максимально и минимально возможные величины числа, которое может быть использованы.

Эти переменные можно использовать в математических выражениях.

1.2.3. Ввод комплексных чисел Язык системы MatLAB, в отличие от многих языков программирования высокого уровня, содержит в себе очень простую в пользовании встроенную арифметику комплексных чисел. Большинство элементарных математических функций допускают в качестве аргументов комплексные числа, а результаты формируются как комплексные числа. Эта особенность языка делает его очень удобным и полезным для инженеров и научных работников.

Для обозначения мнимой единицы в языке MatLAB зарезервированы два имени i и j. Ввод с клавиатуры значения комплексного числа осуществляется путем записи в командное окно строки вида:

где ДЧ - действительная часть комплексного числа, МЧ - мнимая часть. Например:

Из приведенного примера видно, в каком виде система выводит комплексные числа на экран (и на печать).

1.2.4. Элементарные математические функции Общая форма использования функции в MatLAB такова:

В языке MatLAB предусмотрены следующие элементарные арифметические функции.

Тригонометрические и гиперболические функции asinh(Z) - обратный гиперболический синус;

cosh(Z) - гиперболический косинус;

acosh(Z) - обратный гиперболический косинус;

tanh(Z) - гиперболический тангенс;

atan2(X,Y) - четырехквадрантный арктангенс (угол в диапазоне от - до + между горизонтальным правым лучом и лучом, который проходит через точку с координатами X и Y);

atanh(Z) - обратный гиперболический тангенс;

asech(Z) - обратный гиперболический секанс;

acsch(Z) - обратный гиперболический косеканс;

acoth(Z) - обратный гиперболический котангенс.

log10(Z) - округление к ближайшему целому в сторону отрицательной бесfloor(Z) конечности;

- округление к ближайшему целому в сторону положительной бесконечности;

round(Z) mod(X,Y) rem(X,Y) 1.2.5. Специальные математические функции Кроме элементарных в языке MatLAB предусмотрен целый ряд специальных математических функций. Ниже приведен перечень и краткое содержание этих функций. Правила обращения к ним и использования пользователь может отыскать в описаниях этих функций, которые выводятся на экран, если набрать команду help и указать в той же строке имя функции.

- преобразование декартовых координат в сферические;

cart2sph - преобразование декартовых координат в полярные;

cart2pol - преобразование полярных координат в декартовые;

pol2cart - преобразование сферических координат в декартовые.

sph2cart - модифицированная функция Бесселя первого рода;

- модифицированная функция Бесселя второго рода.

gammainc - неполная гамма-функция;

gammaln - логарифм гамма-функции.

- функция экспоненциального интеграла.

- масштабированная дополнительная функция ошибок;

legendre - представление числа в виде рациональной дроби;

- представление чисел в виде рациональной дроби.

1.2.6. Элементарные действия с комплексными числами Простейшие действия с комплексными числами - сложение, вычитание, умножение, деление и возведение в степень - осуществляются при помощи обычных арифметических знаков +, -, *, /, \ и ^ соответственно.

Примеры использования приведены на рис. 1.11.

Примечание. В приведенном фрагменте использована функция disp (от слова 'дисплей'), которая тоже выводит в командное окно результаты вычислений или некоторый текст. При этом численный результат, как видно, выводится уже без указания имени переменной или ans.

1.2.7. Функции комплексного аргумента Практически все элементарные математические функции, приведенные в п.

1.2.4, вычисляются при комплексных значениях аргумента и получают в результате этого комплексные значения результата.

Благодаря этому, например, функция sqrt вычисляет, в отличие от других языков программирования, квадратный корень из отрицательного аргумента, а функция abs при комплексном значении аргумента вычисляет модуль комплексного числа. Примеры приведены на рис. 1.12.

В MatLAB есть несколько дополнительных функций, рассчитанных только на комплексный аргумент:

real(Z) - выделяет действительную часть комплексного аргумента Z;

imag(Z) - выделяет мнимую часть комплексного аргумента;

angle(Z) - вычисляет значение аргумента комплексного числа Z (в радианах в диапазоне от - до + );

conj(Z) - выдает число, комплексно сопряженное относительно Z.

Примеры приведены на рис. 1.13.

Кроме того, в MatLAB есть специальная функция cplxpair(V), которая осуществляет сортировку заданного вектора V с комплексными элементами таким образом, что комплексно-сопряженные пары этих элементов располагаются в векторе-результате в порядке возрастания их действительных частей, при этом элемент с отрицательной мнимой частью всегда располагается первым. Действительные элементы завершают комплексно-сопряженные пары. Например (в дальнейшем в примерах команды, которые набираются с клавиатуры, будут написаны жирным шрифтом, а результат их выполнения - обычным шрифтом):

» v = [ -1, -1+2i,-5,4,5i,-1-2i,-5i] Columns 1 through Columns 5 through » disp(cplxpair(v)) Columns 1 through Columns 5 through Приспособленность большинства функций MatLAB к оперированию с комплексными числами позволяет значительно проще строить вычисления с действительными числами, результат которых является комплексным, например, находить комплексные корни квадратных уравнений.

1.2.8. Задания Задание 1.1. Вычислите указанное арифметическое выражение. Укажите последовательность нажатия клавиша. Сравните полученный результат с приведенным ответом.

Задание 1.2. Проведите вычисления по заданной формуле при заданных значениях параметров. Укажите необходимую последовательность действий.

Сравните полученный результат с приведенным ответом.

Указание. В системе MatLAB несколько последних команд запоминаются.

Повторный вызов этих команд в командное окно осуществляется нажатием клавиш < > и. Используйте эту возможность для повторного обращения к набранной функции.

ОТВЕТ: а) 1. 5633e+008; б) 5. 0651e-002.

ОТВЕТ: а) 1. 0498e+000; б) 1. 2429e-001.

ОТВЕТ: а) 2. 9464e+000; б) 4. 9445e+000.

ОТВЕТ: а) 5. 4283e-001; б) 8. 9703e-018+ 1. 4650e-001i.

ОТВЕТ: а) 1. 1790e+000; б) 1. 6630e+000.

ОТВЕТ: а) 1. 5880e+004; б) 5. 4516e-001.

ОТВЕТ: а) 3. 5028e-002; б) 1. 4003e+000.

ОТВЕТ: а) 4. 1645e+002; б) 4. 1101e+002.

ОТВЕТ: а) 1. 6193e+003; б) 3. 5238e+003.

21. lg 3x x 9 + ; a) x = e1,648 ; б) x = tg1,21.

ОТВЕТ: а) 6. 1109e+000; б) -5. 1927e-001.

22.

ОТВЕТ: а) 8. 5511e+075 ; б) 4. 0272e+000;

ОТВЕТ: а) -2. 0936e+003; б) 1. 7858e-001.

25.

ОТВЕТ: а) 4. 4466e+000; б) 5. 2145e-001.

Задание 1.3. Выполните такие действия (см. таблицу 1.1):

а) число z1, заданное в алгебраической (экспоненциальной) форме, переведите в экспоненциальную (алгебраическую), проверьте и запишите результат;

б) число z2, заданное в экспоненциальной (алгебраической) форме, переведите в алгебраическую (экспоненциальную), проверьте и запишите результат;

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

Задание 1.4. Найдите корни квадратного уравнения при заданных значениях коэффициентов a, b и c (см. таблицу 1.2).

1.2.9. Вопросы 1. Как представляются действительные числа при вычислениях в системе MatLAB?

2. Как изменить формат представления действительных чисел в командном окне?

3. Каким образом объявляются переменные в языке MatLAB?

4. Как сделать так, чтобы результат действий, записанных в очередной строке а) выводился в командное окно; б) не выводился на экран?

5. Какую роль играет системная переменная ans?

6. Как возвратить в командную строку ранее введенную команду?

7. Как ввести значения комплексного числа, и в каком виде оно выведется на экран?

8. Как на языке MatLAB обеспечить сложение, вычитание, умножение, деление и возведение в степень комплексных чисел?

9. Какие функции работы с комплексными числами предусмотрены в языке MatLAB?

1.3. Простейшие операции с векторами и матрицами MatLAB - система, специально предназначенная для осуществления сложных вычислений с векторами, матрицами и полиномами. Под вектором в MatLAB понимается одномерный массив чисел, а под матрицей - двумерный массив. При этом по умолчанию предполагается, что любая заданная переменная является вектором или матрицей. Например, отдельное заданное число система воспринимает как матрицу размером (1*1), а вектор-строку из N элементов - как матрицу размером (1*N).

1.3.1. Ввод векторов и матриц Начальные значения векторов можно задавать с клавиатуры путем поэлементного ввода. Для этого в строке следует сначала указать имя вектора, потом поставить знак присваивания ' = ', затем, - открывающую квадратную скобку, а за ней ввести заданные значения элементов вектора, отделяя их пробелами или запятыми. Заканчивается строка записью закрывающей квадратной скобки.

Например, запись строки V = [ 1.2 -0.3 1.2e-5] задает вектор V, который содержит три элемента со значениями 1.2, -0.3 и 1.2е-5 (рис. 1.14):

После введения вектора система выводит его на экран. То, что в приведенном примере последний элемент выведен как 0, обусловлено установленным форматом short, в соответствии с которым выводятся данные на экран.

Длинный вектор можно вводить частями, которые потом объединять с помощью операции объединения векторов в строку : v = [ v1 v2 ]. Например:

Язык MatLAB дает пользователю возможность сокращенного введения вектора, значения элементов которого составляют арифметическую прогрессию.

Если обозначить: nz - начальное значение этой прогрессии (значение первого элеОперации с векторами и матрицами мента вектора); kz - конечное значение прогрессии (значение последнего элемента вектора); h - разность прогрессии (шаг), то вектор можно ввести с помощью короткой записи V = nz : h : kz. Например, введение строки V = - 0.1 : 0.3 :

1.4 приведет к такому результату:

Если средний параметр (разность прогрессии) не указан, то он по умолчанию принимается равным единице. Например, команда приводит к формированию такого вектора ans = -2.1000 -1.1000 -0.1000 0.9000 1.9000 2.9000 3.9000 4. Так вводятся векторы-строки. Вектор-столбец вводится аналогично, но значения элементов отделяются знаком ";".

Ввод значений элементов матрицы осуществляется в MatLAB в квадратных скобках, по строкам. При этом элементы строки матрицы один от другого отделяются пробелом или запятой, а строки одна от другой отделяются знаком ";" (рис. 1.17).

1.3.2. Формирование векторов и матриц MatLAB имеет несколько функций, которые позволяют формировать векторы и матрицы некоторого определенного вида. К таким функциям относятся:

zeros(М,N) - создает матрицу размером (М*N) с нулевыми элементами, например:

» zeros(3,5) ones(М,N) - создает матрицу размером (М*N) с единичными элементами, например:

» ones(3,5) eye(М,N) - создает единичную матрицу размером (М*N), т. е. с единицами по главной диагонали и остальными нулевыми элементами, например:

rand(М,N) - создает матрицу размером (М*N) из случайных чисел, равномерно распределенных в диапазоне от 0 до 1, например:

» rand(3,5) 2.1896e-001 6.7930e-001 5.1942e-001 5.3462e-002 7.6982e- 4. 7045e-002 9. 3469e-001 8. 3097e-001 5. 2970e-001 3. 8342e- 6. 7886e-001 3. 8350e-001 3. 4572e-002 6. 7115e-001 6. 6842e- randn(М,N) - создает матрицу размером (М*N) из случайных чисел, распределенных по нормальному (гауссовому) закону с нулевым математическим ожиданием и стандартным (среднеквадратичным) отклонением, равным единице, например:

» randn(3,5) hadamard(N) - создает матрицу Адамара размером (N*N), например:

hilb(N) - создает матрицу Гільберта размером (N*N), например:

invhilb(N) - создает обратную матрицу Гильберта размером (N*N), например:

pascal(N) - создает матрицу Паскаля размером (N*N), например:

В языке MatLAB предусмотрено несколько функций, которые позволяют формировать матрицу на основе другой (заданной) или используя некоторый заданный вектор. К таким функциям принадлежат:

fliplr(A) - формирует матрицу, переставляя столбцы известной матрицы А относительно вертикальной оси, например:

» fliplr(A) flipud(A) - переставляет строки заданной матрицы А относительно горизонтальной оси, например:

» flipud(A) rot90(A) - формирует матрицу путем "поворота" заданной матрицы А на 90 градусов против часовой стрелки:

reshape(A,m,n) - образует матрицу размером (m*n) путем выборки элементов заданной матрицы А по столбцам и последующему распределению этих элементов по 'n' столбцам, каждый из которых содержит 'm' элементов; при этом число элементов матрицы А должно равняться m*n, например:

» reshape(A,2,9) tril(A) - образует нижнюю треугольную матрицу на основе матрицы А путем обнуления ее элементов выше главной диагонали:

triu(A) - образует верхнюю треугольную матрицу на основе матрицы А путем обнуления ее элементов ниже главной диагонали:

hankel(V) - образует квадратную матрицу Ганкеля, первый столбец которой совпадает с заданным вектором V, например:

» hankel(V) Процедура diag(х) - формирует или извлекает диагональ матрицы.

Если х - вектор, то функция diag(х) создает квадратную матрицу с вектором х на главной диагонали:

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

Если х - матрица, то функция diag создает вектор-столбец, который состоит из элементов главной диагонали заданной матрицы х, например, для матрицы А, указанной перед примером применения процедуры fliplr:

Если при этом указать дополнительно номер диагонали, то можно получить вектор-столбец из элементов любой диагонали матрицы х, например:

» diag(A,3) Функция zeros(1,N) формирует (создает) вектор-строку из N нулевых элементов. Аналогично zeros(N,1) создает вектор-столбец из N нулей.

Векторы, значения элементов которых являются случайными равномерно распределенными, формируются таким образом: rand(1,n) - для вектора-строки и rand(m,1) - для вектора-столбца.

1.3.3. Извлечение и вставка частей матриц Прежде всего отметим, что обращение к любому элементу заданной матрицы в MatLAB осуществляется путем указания (в скобках, через запятую) после имени матрицы двух целых положительных чисел, которые определяют соответственно номера строки и столбца матрицы, на пересечении которых расположен этот элемент.

Пусть имеем некоторую матрицу А:

Тогда получить значение элемента этой матрицы, расположенного на пересечении второй строки с третьим столбиком, можно следующим образом:

Если нужно, наоборот, установить на это место некоторое число, например,, то это можно сделать так:

5. 0000 6. 0000 3. 1416 8. 9. 0000 10. 0000 11. 0000 12. Иногда нужно создать меньшую матрицу из большей, формируя ее путем извлечения из последней матрицы элементов ее нескольких строк и столбцов, или, наоборот, вставить меньшую матрицу таким образом, чтобы она стала определенной частью матрицы большего размера. Это в MatLAB делается с помощью знака двоеточия (" : ").

Рассмотрим эти операции на примерах.

Пусть нужно создать вектор V1, состоящий из элементов третьего столбца последней матрицы А. Для этого произведем такие действия:

Чтобы создать вектор V2, состоящий из элементов второй строки матрицы А, поступают так:

Допустим, что необходимо из матрицы А образовать матрицу В размером (2*2), которая состоит из элементов левого нижнего угла матрицы А. Тогда делают так:

Аналогично можно вставить матрицу В в верхнюю середину матрицы А:

>> A(1:2,2:3)=B Как видно, для этого вместо указания номеров элементов матрицы можно указывать диапазон изменения этих номеров путем указания нижней и верхней границ, разделяя их двоеточием.

Примечание. Если верхней границей изменения номеров элементов матрицы является ее размер в этом измерении, вместо него можно использовать служебное слово end. Например:

>> A(2:end,2:end) Эти операции очень удобны для формирования матриц, большинство элементов которых одинаковы, в частности, так называемых разреженных матриц, которые состоят, в основном, из нулей, за исключением отдельных элементов.

Для примера рассмотрим формирование разреженной матрицы размером (5*7) с единичными элементами в ее центре:

>> A = zeros(5,7);

>> B = ones(3,3);

>> A(2:4,3:5)=B "Растянуть" матрицу (А) в единый вектор (V) можно с помощью обычной записи "V = A(:)". При этом создается вектор-столбец с количеством элементов (m*n), в котором столбцы заданной матрицы размещены сверху вниз в порядке самих столбцов:

Наконец, "расширить" матрицу, составляя ее из отдельных заданных матриц ("блоков") можно тоже довольно просто. Если заданы несколько матриц-блоков А1, А2,... АN с одинаковым количеством строк, то из них можно "слепить" единую матрицу А, объединяя блоки в одну "строку" таким образом:

Эту операцию называют горизонтальной конкатенацией (сцеплением) матриц. Вертикальная конкатенация матриц реализуется (при условии, что все составные блоки-матрицы имеют одинаковое количество столбцов) аналогично, путем применения для отделения блоков вместо запятой точки с запятой:

Приведем примеры. Пример горизонтальной конкатенации:

>> A2 = [10;11;12];

Пример вертикальной конкатенации:

1.3.4. Действия над векторами Будем различать две группы действий над векторами:

а) векторные действия - т. е. такие, которые предусмотрены векторным исчислением в математике;

б) действия по преобразованию элементов - это действия, которые преобразуют элементы вектора, но не являются операциями, разрешенными математикой.

Сложение векторов. Как известно, суммироваться могут только векторы одинакового типа (т. е. такие, которые являются или векторами-строками, или векторами-столбцами), имеющие одинаковую длину (т. е. одинаковое количество элементов). Если X и Y - именно такие векторы, то их сумму Z можно получить, введя команду Z = X + Y, например:

Аналогично с помощью арифметического знака ' - ' осуществляется вычитание векторов, имеющих одинаковую структуру (Z = X - Y).

Например:

Транспонирование вектора осуществляется применением знака апострофа, который записывается сразу за записью имени вектора, который нужно транспонировать. Например:

Умножение вектора на число осуществляется в MatLAB с помощью знака арифметического умножения ' * ' таким образом: Z = X*r или Z = r*X, где r - некоторое действительное число.

Умножение двух векторов определено в математике только для векторов одинакового размера (длины) и лишь тогда, когда один из векторов-множителей строка, а второй - столбец. Иначе говоря, если векторы X и Y являются строками, то математическое смысл имеют лишь две формы умножения этих векторов: U = X' * Y и V = X * Y'. При этом в первом случае результатом будет квадратная матрица, а во втором - число. В МatLAB умножение векторов осуществляется применением обычного знака умножения ' * ', который записывается между множителями-векторами.

Для трехкомпонентных векторов в MatLAB существует функция cross, которая позволяет найти векторное произведение двух векторов. Для этого, если заданы два трехкомпонентных вектора v1 и v2, достаточно ввести оператор cross(v1, v2).

» cross(v1,v2) На этом перечень допустимых математических операций с векторами исчерпывается.

В языке MatLAB предусмотрен ряд операций, которые преобразуют заданный вектор в другой того же размера и типа, но не являются операциями с вектором как математическим объектом. К таким операциям относятся, например, все элементарные математические функции, приведенные в разделе 1.2.4 и которые зависят от одного аргумента. В языке MatLAB запись, например, вида Y = sin(X), где X - некоторый известный вектор, приводит к формированию нового вектора Y, имеющего тот же тип и размер, но элементы которого равняются синусам соответствующих элементов вектора-аргумента X. Например:

» x = [ -2,-1,0,1,2];

Кроме этих операций в МаtLAB предусмотрено несколько операций поэлементного преобразования, осуществляемых с помощью знаков обычных арифметических действий. Эти операции применяются к векторам одинакового типа и размера. Результатом их есть вектор того же типа и размера.

Добавление (отнимание) числа к (из) каждому элемента вектора. Осуществляется с помощью знака ' + ' (' - ').

Поэлементное умножение векторов. Проводится с помощью совокупности знаков '.* ', которая записывается между именами перемножаемых векторов.

В результате получается вектор, каждый элемент которого является произведением соответствующих элементов векторов - "сомножителей".

Поэлементное деление векторов. Осуществляется с помощью совокупности знаков './ '. Результат - вектор, каждый элемент которого является частным от деления соответствующего элемента первого вектора на соответствующий элемент второго вектора.

Поэлементное деление векторов в обратном направлении. Осуществляется с помощью совокупности знаков '.\ '. В результате получают вектор, каждый элемент которого является частным от деления соответствующего элемента второго вектора на соответствующий элемент первого вектора.

Поэлементное возведение в степень. Осуществляется с помощью совокупности знаков '.^ '. Результат - вектор, каждый элемент которого является соответствующим элементом первого вектора, возведенным в степень, величина которой равняется значению соответствующего элемента второго вектора. Примеры:

» x = [1,2,3,4,5]; y = [-2,1,4,0,5];

Warning: Divide by zero Вышеуказанные операции позволяют очень просто вычислять (а затем строить графики) сложные математические функции, не используя при этом операторы цикла, т. е. осуществлять построение графиков в режиме калькулятора.

Для этого достаточно задать значение аргумента как арифметическую прогрессию так, как это было показано в п. 1.3.1, а потом записать нужную функцию, используя знаки поэлементного преобразования векторов.

Например, пусть нужно вычислить значения функции:

при значениях аргумента х от 0 до 10 с шагом 1. Вычисление массива значений этой функции в указанных условиях можно осуществить с помощью лишь двух простых операторов :

Columns 1 through 0 1.5311 1.0035 0.0945 -0.3073 -0.2361 -0. Columns 8 through 0. 0595 0. 0544 0. 0137 -0. 1.3.5. Поэлементное преобразование матриц Для поэлементного преобразования матрицы пригодны все указанные ранее в п. 1.2.4 алгебраические функции. Каждая такая функция формирует матрицу того же размера, что и заданная, каждый элемент которой вычисляется как указанная функция от соответствующего элемента заданной матрицы. Кроме этого, в MatLAB определены операции поэлементного умножения матриц одинакового размера (совокупностью знаков '.* ', записываемой между именами перемножаемых матриц), поэлементного деления (совокупности './ ' и '.\'), поэлементного возведения в степень (совокупность '.^' ), когда каждый элемент первой матрицы возводится в степень, равную значению соответствующего элемента второй матрицы.

Приведем несколько примеров:

» A = [1,2,3,4,5; -2, 3, 1, 4, 0] » B = [-1,3,5,-2,1; 1,8,-3,-1,2] Warning: Divide by zero Оригинальной в языке MatLAB является операция прибавления к матрице числа. Она записывается следующим образом: A + x, или х + A (А - матрица, а x число). Такой операции нет в математике. В MatLAB она эквивалентна совокупности операций где Е - обозначение матрицы, которая состоит только из единиц, тех же размеров, что и матрица А. Например:

1.3.6. Матричные действия над матрицами К матричным действиям над матрицами относят такие операции, которые используются в матричном исчислении в математике и не противоречат ему.

Базовые действия с матрицами - сложение, вычитание, транспонирование, умножение матрицы на число, умножение матрицы на матрицу, возведение матрицы в целую степень - осуществляются в языке MatLAB с помощью обычных знаков арифметических операций. При использовании этих операций важно помнить условия, при которых эти операции являются возможными:

при сложении или вычитании матрицы должны иметь одинаковые размеры;

при умножении матриц количество столбцов первой матрицы должно совпадать с количеством строк второй матрицы.

Невыполнение этих условий приведет к появлению в командном окне сообщения об ошибке. Приведем несколько примеров.

Пример сложения и вычитания:

Пример умножения на число:

Пример транспонирования матрицы:

Пример умножения матрицы на матрицу:

Функция обращения матрицы - inv(A) - вычисляет матрицу, обратную заданной матрице А. Исходная матрица А должна быть квадратной, а ее определитель не должен равняться нулю.

Приведем пример:

-2. 6000e-001 1. 0000e- -8. 1739e-002 3. 4783e- Проверим правильность выполнения операции обращения, применяя ее еще раз к полученному результату:

-4. 0000e+001 1. 1500e+ -9. 4000e+001 2. 9900e+ Как видим, мы получили исходную матрицу С, что является признаком правильности выполнения обращения матрицы.

Возведение матрицы в целую степень осуществляется в MatLAB с помощью знака "^": А^n. При этом матрица должна быть квадратной, а n - целым (положительным или отрицательным) числом. Это матричное действие эквивалентно умножению матрицы А на себя n раз (если n - положительно) или умножению обратной матрицы на себя (при n отрицательно).

Приведем пример:

1.5385e-001 -7.6923e-002 3.0769e- 7. 6923e-002 3. 0769e-001 -4. 6154e- 2. 1328e-018 -1. 5385e-001 3. 8462e- Оригинальными в языке MatLAB являются две новые, неопределяемые в математике функции деления матриц. При этом вводятся понятие деления матриц слева направо и деление матриц справа налево. Первая операция записывается с помощью знака ' / ', а вторая - ' \ '.

Операция В / A эквивалентна последовательности действий B * inv(A), где функция inv осуществляет обращение матрицы. Ее удобно использовать для решения матричного уравнения:

Аналогично операция A\B равносильна совокупности операций inv(A)*B, которая представляет собой решение матричного уравнения:

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

В среде MatLAB это можно сделать таким образом:

1.3.7. Матричные функции Вычисление матричной экспоненты (eА) осуществляется с помощью функций expm, expm1, expm2, expm3. Эти функции следует отличать от прежде рассмотренной функции exp(A), которая формирует матрицу, каждый элемент которой равняется е в степени, которая равняется соответствующему элементу матрицы А.

Функция expm является встроенной функцией MatLAB. Функция expm1(A) реализована как М-файл, который вычисляет матричную экспоненту путем использования разложения Паде матрицы А. Функция еxpm2(A) вычисляет матричную экспоненту, используя разложение Тейлора матрицы А. Функция expm3(A) вычисляет матричную экспоненту на основе использования спектрального разложения А.

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

» A = [1,2,3; 0, -1,5;7, -4,1] 131.3648 -9.5601 80. 97. 8030 -7. 1768 59. 123. 0245 -8. 8236 75. 131.3648 -9.5601 80. 97. 8030 -7. 1768 59. 123. 0245 -8. 8236 75. 131.3648 -9.5601 80. 97. 8030 -7. 1768 59. 123. 0245 -8. 8236 75. Функция logm(А) осуществляет обратную операцию - логарифмирование матрицы по натуральному основанию, например:

Функция sqrtm(А) вычисляет такую матрицу Y, что Y*Y = A:

0. 0000 - 0. 0000i 1. 0000 - 0. 0000i 5. 0000 - 0. 0000i 1.3.8. Задания Задание 1.5. Вычислите значения функции f(x) на отрезке [a; b] с шагом h.

1.3.9. Вопросы 1. Как вводятся векторы в языке MatLAB? Какими функциями можно формировать векторы в языке MatLAB?

2. Какие функции MatLAB разрешают преобразовывать вектор поэлементно?

3. С помощью каких средств в MatLAB осуществляются основные операции с векторами?

4. Как вводятся матрицы в системе MatLAB?

5. Какие функции имеются в MatLAB для формирования матриц определенного вида?

6. Как сформировать матрицу: а) по заданным векторам ее строк? б) по заданным векторам ее столбцов? в) по заданным векторам ее диагоналей?

7. Какие функции поэлементного преобразования матрицы есть в MatLAB?

8. Как осуществляются в MatLAB обычные матричные операции?

9. Как решить в MatLAB систему линейных алгебраических уравнений?

1.4. Функции прикладной численной математики 1.4.1. Операции с полиномами В системе MatLAB предусмотрены некоторые дополнительные возможности математического оперирования с полиномами.

Полином (многочлен) как функция определяется выражением:

В системе MatLAB полином задается и сохраняется в виде вектора, элементами которого являются коэффициенты полинома от a n до a 0 в указанном порядке:

Введение полинома в MatLAB осуществляется так же, как и ввод вектора длиной n+1, где n - порядок полинома.

Умножение полиномов. Произведением двух полиномов степеней n и m соответственно, как известно, называют полином степени n+m, коэффициенты которого определяют простым перемножением этих двух полиномов. Фактически операция умножения двух полиномов сводится к построению расширенного вектора коэффициентов по заданным векторам коэффициентов полиномовсомножителей. Эту операцию в математике называют сверткой векторов (а сам вектор, получаемый в результате такой процедуры - вектором-сверткой двух векторов). В MatLAB ее осуществляет функция conv(P1, P2).

Аналогично, функция deconv(P1, P2) осуществляет деление полинома P1 на полином P2, т. е. обратную свертку векторов P1 и P2. Она определяет коэффициенты полинома, который является частным от деления P1 на P2.

» p1 = [1,2,3]; p2 = [1,2,3,4,5,6];

» p = conv(p1,p2) » deconv(p,p1) В общем случае деление двух полиномов приводит к получению двух полиномов - полинома-результата (частного) и полинома-остатка. Чтобы получить оба этого полинома, следует оформить обращение к функции таким образом:

[Q,R] = deconv(B,A).

Тогда результат будет выдан в виде вектора Q с остатком в виде вектора R таким образом, что будет выполненное соотношение Система MatLAB имеет функцию roots(P), которая вычисляет вектор, элементы которого являются корнями заданного полинома Р.

Пусть нужно найти корни полинома:

Ниже показано, как просто это сделать:

» p = [1,8,31,80,94,20];

» disp(roots(p)) -1.0000 + 3.0000i -1.0000 - 3.0000i Обратная операция - построение вектора р коэффициентов полинома по заданному вектору его корней - осуществляется функцией poly:

Здесь r - заданный вектор значений корней, p - вычисленный вектор коэффициентов полинома. Приведем пример:

» p = [1,8,31,80,94,20] -1.0000 + 3.0000i -1.0000 - 3.0000i p1 =8. 0000 31. 0000 80. 0000 94. 0000 20. Заметим, что получаемый вектор не показывает старшего коэффициента, который по умолчанию полагается равным единице.

Эта же функция в случае, если аргументом ее является некоторая квадратная матрица А размером (n*n), строит вектор характеристического полинома этой матрицы. Обращение формирует вектор p коэффициентов характеристического полинома p(s) = det(s*E - A) = p1*sn +... + pn*s + pn+1, где Е - обозначение единичной матрицы размером (n*n).

Рассмотрим пример:

1. 0000 -10. 0000 20. 0000 -36. Для вычисления значения полинома по заданному значению его аргумента в MatLAB предусмотрена функция polyval. Обращение к ней осуществляется по схеме:

y = polyval(p,x), где p - заданный вектор коэффициентов полинома, а x - заданное значение аргумента. Пример:

» y = polyval(p,2) Если в качестве аргумента полинома указана матрица Х, то функция polyval(p,X) вычисляет матрицу Y, каждый элемент которой является значением указанного полинома при значении аргумента, равном соответствующему элементу матрицы Х, например:

» disp(polyval(p,X)) В этом случае функция вычисляет значение полинома для каждого элемента матрицы Х, и поэтому размеры исходной и конечной матриц одинаковы size(Y) = size(X).

Вычисление производной от полинома осуществляется функцией polyder.

Эта функция создает вектор коэффициентов полинома, представляющего собой производную от заданного полинома. Она имеет три вида обращений:

dp = polyder(p) по заданному полиному р вычисляет вектор dp, элементы которого являются коэффициентами полинома-производной от заданного:

» dp = polyder(p) dp = polyder(p1,p2) вычисляет вектор dp, элементы которого являются коэффициентами полинома-производной от произведения двух полиномов р1 и р2:

» p1 = [1,8,31,80,94,20];

» p = conv(p1,p2) » dp = polyder(p) » dp1 = polyder(p1,p2) [q,p] = polyder(p1,p2) вычисляет производную от отношения (p1/p2) двух полиномов р1 и р2 и выдает результат в виде отношения (q/p) полиномов q и p:

» p1 = [1,8,31,80,94,20];

» [q,p] = polyder(p1,p2) » y = deconv(p1,p2) » z1 = polyder(y) 1.4.2. Обработка данных измерений Система MatLAB дает пользователю дополнительные возможности для обработки данных, которые заданы в векторной или матричной форме.

Допустим, что есть некоторая зависимость y(x), заданная рядом точек Ее можно задать в командном окне MatLAB как матрицу xydata, содержащую две строки - значения x и значения y:

>> xydata =[2 4 6 8 10; 5.5 6.3 6.8 8 8.6] На примере этой зависимости рассмотрим основные средства для обработки данных.

Функция size(xydata) предназначена для определения числа строк и столбцов матрицы xydata. Она формирует вектор [n, p], содержащий эти величины:

>> size(xydata) Обращение к ней вида >> [n, p] = size(xydata);

позволяет сохранить в памяти машины и использовать потом при дальнейших вычислениях данные о числе строк n и столбцов p этой матрицы:

С помощью этой функции можно установить длину и тип (строка или столбец) вектора:

v1 = 2. 0000 5. 5000 4. 0000 6. 3000 6. 0000 6. 8000 8. 0000 8. 0000 10. 0000 8. Функция max(V), где V - некоторый вектор, выдает значение максимального элемента этого вектора. Аналогично, функция min(V) извлекает минимальный элемент вектора V. Функции mean(V) и std(V) определяют, соответственно, среднее значение и среднеквадратичное отклонение от него значений элементов вектора V.

Функция сортировки sort(V) формирует вектор, элементы которого расположены в порядке возрастания их значений.

Функция sum(V) вычисляет сумму элементов вектора V.

Функция prod(V) выдает произведение всех элементов вектора V.

Функция cumsum(V) формирует вектор того же типа и размера, любой элемент которого является суммой всех предшествующих элементов вектора V (вектор кумулятивной суммы).

Функция cumprod(V) создает вектор, элементы которого являются произведением всех предшествующих элементов вектора V.

Функция diff(V) создает вектор, который имеет размер на единицу меньший размера вектора V, элементы которого являются разностью между соседними элементами вектора V.

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

» v = [ 1, 0.1, 0.5, 0.1, 0.1,0.4 ];

» disp(size(v)) » disp(max(v)) » disp(min(v)) » disp(mean(v)) » disp(std(v)) » disp(sort(v)) » disp(sum(v)) » disp(prod(v)) » disp(cumsum(v)) » disp(cumprod(v)) » disp(diff(v)) Если указать второй выходной параметр, то можно получить дополнительную информацию об индексе первого элемента, значение которого является максимальным или минимальным:

>> [M,n]=max(v) >> [N,m]=min(v) Интегрирование методом трапеций осуществляет процедура trapz. Обращение к ней вида trapz(x,y) приводит к вычислению площади под графиком функции y(x), в котором соседние точки, заданные векторами х і у, соединены отрезками прямых. Если первый вектор х не указан в обращении, по умолчанию допускается, что шаг интегрирования равняется единице (т. е. вектор х представляет собой вектор из номеров элементов вектора у).

Пример. Вычислим интеграл от функции y = sin(x) в диапазоне от 0 до.

Его точное значение равно 2. Возьмем равномерную сетку из 100 элементов. Тогда вычисления сведутся к совокупности операций:

» disp(trapz(x,y)) Те же функции size, max, min, mean, std, sort, sum, prod, cumsum, cumprod, diff могут быть применены и к матрицам. Основным отличием использования в качестве аргументов этих функций именно матриц является то, что соответствующие описанные выше операции ведутся не по отношению к строкам матриц, а к каждому из столбцов заданной матрицы. Т. е. каждый столбец матрицы А рассматривается как переменная, а каждая строка - как отдельное наблюдение. Так, в результате применения функций max, min, mean, std получаются векторы-строки с количеством элементов, которое равняется количеству столбцов заданной матрицы. Каждый элемент содержит, соответственно, максимальные, минимальное, среднее или среднеквадратичное значения элементов соответствующего столбца заданной матрицы.

Приведем примеры. Пусть имеем 3 величины y1, y2 и y3, измеренные при некоторых пяти значениях аргумента (они не указаны). Тогда данные измерений образуют 3 вектора по 5 элементов:

Сформируем из них матрицу измерений так, чтобы векторы y1, y2 и y3 образовывали столбцы этой матрицы:

5.5000 -1.2000 3. 6.3000 0.5000 5. 6.8000 -0. 8. 6000 0. 1000 10. Применим к этой матрице измерений описанные функции. Получим Если при обращении к функциям max и min указать второй выходной параметр, то он даст информацию о номерах строк, где находятся в соответствующем столбце первые элементы с максимальным (или минимальным) значением. Например:

>> [M,n]=max(A) M= 8.6000 1.0000 10. >> [N,m]=min(A) N = 5.5000 -1. Функция sort сортирует элементы любого из столбцов матрицы. Результатом является матрица того же размера.

Функция sum и prod формируют вектор-строку, каждый элемент которого является суммой или произведением элементов соответствующего столбца исходной матрицы.

Функции cumsum, cumprod образуют матрицы того же размера, элементы каждого столбца которых являются суммой или произведением элементов этого же столбца начальной матрицы, начиная с соответствующего элемента и выше.

Наконец, функция diff создает из заданной матрицы размером (m*n) матрицу размером ((m-1)*n), элементы которой являются разностью между элементами соседних строк начальной матрицы.

Применяя эти процедуры к принятой матрице измерений, получим:

6.3000 -0.6000 3. 6.8000 0.1000 5. ans = 35. 2000 -0. 2000 27. » cumsum(A) » cumprod(A) 0.8000 1.7000 2. 0.6000 -0. 9000 1. Рассмотрим некоторые другие функции, предоставляемые пользователю системой MatLAB.

Функция cov(A) вычисляет матрицу ковариаций измерений. При этом получают квадратную симметричную матрицу с количеством строк и столбцов, равФункции прикладной численной математики ным количеству измеренных величин, т. е. количеству столбцов матрицы измерений. Например, при применении к принятой матрице измерений она дает такой результат:

На диагонали матрицы ковариаций размещены дисперсии измеренных величин, а вне ее - взаимные корреляционные моменты этих величин.

Функция corrcoeff(A) вычисляет матрицу коэффициентов корреляции при тех же условиях. Элементы матрицы S = corrcoef(A) связаны с элементами матрицы ковариаций C=cov(A) таким соотношением:

» corrcoef(A) 1.4.3. Функции линейной алгебры Традиционно к линейной алгебре относят такие задачи, как обращение и псевдообращение матрицы, спектральное и сингулярное разложения матриц, вычисление собственных значений и векторов, сингулярных чисел матриц, вычисление функций от матриц. Ознакомимся с некоторыми основными функциями MatLAB в этой области.

Функция k = cond(A) вычисляет и выдает число обусловленности матрицы относительно операции обращения, которое равняется отношению максимального сингулярного числа матрицы к минимальному.

Функция k = norm(v,p) вычисляет р-норму вектора v по формуле:

k = sum(abs(v). p)(1/p), где р - целое положительное число. Если аргумент р при обращении к функции не указан, вычисляется 2-норма.

Функция k = norm(А,p) вычисляет р-норму матрицы, где р = 1,2, 'fro' или inf. Если аргумент р не указан, вычисляется 2-норма. При этом справедливы такие соотношения:

Функция rd = rcond(А) вычисляет величину, обратную значению числа обусловленности матрицы А относительно 1-нормы. Если матрица А хорошо обуФункции прикладной численной математики словлена, значение rd близко к единице. Если же она плохо обусловленная, rd приближается к нулю.

Функция r = rank(A) вычисляет ранг матрицы, который определяется как количество сингулярных чисел матрицы, превышающие порог max(size(A))*norm(A)*eps.

Приведем примеры применения этих функций:

» disp(cond(A)) » disp(norm(A,1)) » disp(norm(A)) » disp(rcond(A)) » disp(rank(A)) Процедура d = det(A) вычисляет определитель квадратной матрицы на основе треугольного разложения методом исключения Гаусса.

Функция t = trace(A) вычисляет след матрицы А, равный сумме ее диагональных элементов.

Q = null(A) вычисляет ортонормированный базис нуль-пространства матрицы А.

Q = orth(A) выдает ортонормированный базис матрицы А.

Процедура R = rref(A) осуществляет приведение матрицы к треугольному виду на основе метода исключения Гаусса с частичным выбором ведущего элемента.

» disp(det(A)) » disp(trace(A)) » disp(null(A)) » disp(orth(A)) 0.8982 -0. 4082 0. » disp(rref(A)) Функция R=chol(A) осуществляет разложение Холецького для действительных симметричных и комплексных эрмитовых матриц. Например:

» disp(chol(A)) Функция lu(A) осуществляет LU-разложение матрицы А в виде произведения нижней треугольной матрицы L (возможно, с перестановками) и верхней треугольной матрицы U, так что A = L * U.

Обращение к этой функции вида позволяет получить три составляющие этого разложения - нижнюю треугольную матрицу L, верхнюю треугольную U и матрицу перестановок P такие, что Приведем пример:

» disp(lu(A)) 3.0000 8.0000 400. Из него вытекает, что в первом, упрощенном варианте обращения функция выдает комбинацию из матриц L и U.

Обращение матрицы осуществляется с помощью функции inv(А):

» disp(inv(A)) -0. 1806 0. 0910 -0. -0. 0067 -0. 0005 0. Процедура pinv(А) находит матрицу, псевдообратную матрице А, которая имеет размеры матрицы А’ и удовлетворяет условиям -0.0423 0. -0.0423 0. Для квадратных матриц эта операция равнозначна обычному обращению.

Процедура [ Q, R, P ] = qr(A) осуществляет разложение матрицы А на три унитарную матрицу Q, верхнюю треугольную R с диагональными элементами, уменьшающимися по модулю, и матрицу перестановок P - такие что Например:

» [Q,R,P] = qr(A) Q = -0.5547 -0. R = -7.2111 -2.7735 -4.9923 -4.7150 -0. 0 -4.1603 -0.2774 1.9415 -2. Определение характеристического полинома матрицы A можно осуществить с помощью функции poly(A). Обращение к ней вида p=poly(A) дает возможность найти вектор-строку p коэффициентов характеристического полинома p(s) = det(s*E - A) = p1*sn +... + pn*s + pn+1, где Е - обозначение единичной матрицы размером (n*n). Например :

1. 0000 -10. 0000 20. 0000 -36. Вычисление собственных значений и собственных векторов матрицы осуществляет процедура eig(А). Обычное обращение к ней позволяет получить вектор собственных значений матрицы А, т. е. корней характеристического полинома матрицы. Если же обращение имеет вид:

то в результате получают диагональную матрицу D собственных значений и матрицу R правых собственных векторов, которые удовлетворяют условию Эти векторы являются нормированными так, что норма любого из них равна единице. Приведем пример:

» disp(eig(A)) 0.9979 -0.0798 -0. 0.0492 -0.3915 -0. 0.0416 -0.9167 0. Сингулярное разложение матрицы осуществляет процедура svd(А). Упрощенное обращение к ней позволяет получить сингулярные числа матрицы А.

Более сложное обращение вида:

позволяет получить три матрицы - U, которая состоит из ортонормированных собственных векторов, отвечающих наибольшим собственным значениям матрицы А*АT; V - из ортонормированных собственных векторов матрицы АT*А и S диагональную матрицу, которая содержит неотрицательные значения квадратных корней из собственных значений матрицы АT*А (их называют сингулярными числами). Эти матрицы удовлетворяют соотношению:

Рассмотрим пример:

» disp(svd(A)) » [U,S,V] = svd(A) -0.0207 0.1806 -0. -0.0869 0.9795 0. -0.9960 -0.0892 0. -0. 9978 -0. 0453 -0. -0. 0442 0. 9987 -0. Приведение матрицы к форме Гессенберга осуществляется процедурой hess(А). Например:

» disp(hess(A)) Более развернутое обращение [P,H] = hess(A) дает возможность получить, кроме матрицы Н в верхней форме Гессенберга, также унитарную матрицу преобразований Р, которая удовлетворяет условиям:

» [P,H] = hess(A) 1.0000 -3.3340 -1. 5. 0990 25. 5000 96. Процедура schur (А) предназначена для приведения матрицы к форме Шура. Упрощенное обращение к ней приводит к получению матрицы в форме Шура.

Комплексная форма Шура - это верхняя треугольная матрица с собственными значениями на диагонали. Действительная форма Шура сохраняет на диагонали только действительные собственные значения, а комплексные изображаются в виде блоков (2*2), частично занимая нижнюю поддиагональ.

Обращение [U,T] = schur(A) позволяет, кроме матрицы Т Шура, получить также унитарную матрицу U, удовлетворяющую условиям:

Если начальная матрица А является действительной, то результатом будет действительная форма Шура, если же комплексной, то результат выдается в виде комплексной формы Шура.

Приведем пример:

» disp(schur(A)) 1.2234 -6.0905 -4. » [U,T] = hess(A) 1. 0000 -3. 3340 -1. 5. 0990 25. 5000 96. Функция [U,T] = rsf2csf(U,T) преобразует действительную квазитреугольную форму Шура в комплексную треугольную:

» [U,T] = rsf2csf(U,T) -0.9934 -0. -0.0449 0.3892 -0. -0.1055 0.9140 0. Процедура [AA, BB, Q, Z, V] = qz(A,B) приводит пару матриц А и В к обобщенной форме Шура. При этом АА и ВВ являются комплексными верхними треугольными матрицами, Q, Z - матрицами приведения, а V - вектором обобщенных собственных векторов такими, что Обобщенные собственные значения могут быть найдены, исходя из такого условия:

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

Пусть задана система обычных дифференциальных уравнений в неявной форме Коши с одним входом u и одним выходом y такого вида:

причем матрицы Q, R и векторы b, c и d равны соответственно Необходимо вычислить значения полюсов и нулей соответствующей передаточной функции.

Эта задача сводится к отысканию собственных значений, которые удовлетворяют матричным уравнениям:

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

Ниже приведена совокупность операторов, которая приводит к расчету полюсов:

» [AA, BB] = qz(R,-Q) % Приведение матриц к форме Шура 5.5039 + 2.7975i 24.8121 -25.3646i -0. 6457 + 0. 7622i -0. 1337 + 0. 1378i » diag(AA). /diag(BB) % Расчет полюсов Расчет нулей осуществляется таким образом:

31.0963 -0.7169 -36. 0.0000 1.0647 0. Вычисление собственных значений матричного полинома осуществляет процедура polyeig. Обращение [ R, d ] = polyeig(A0, A1,..., Ap ) позволяет решить полную проблему собственных значений для матричного полинома степени р вида Входными переменными этой процедуры являются р+1 квадратные матрицы A0, A1,... Ap порядка n. Исходными переменными - матрица собственных векторов R размером (n*(n*p)) и вектор d собственных значений размером (n*p).

Функция polyvalm предназначена для вычисления матричного полинома вида по заданному значению матрицы Х и вектора p = [pn, pn-1,..., p0] коэффициентов полинома. Для этого достаточно обратиться к этой процедуре по схеме:

» disp(polyvalm(p,X)) Примечание. Следует различать процедуры polyval и polyvalm. Первая вычисляет значение полинома для любого из элементов матрицы аргумента, а вторая при вычислении полинома возводит в соответствующую степень всю матрицу аргумента.

Процедура subspace(А,В) вычисляет угол между двумя подпространствами, которые "натянуты на столбцы" матриц А и В. Если аргументами являются не матрицы, а векторы A и B, вычисляется угол между этими векторами.

1.4.4. Аппроксимация и интерполяция данных Полиномиальная аппроксимация данных измерений, которые сформированы как некоторый вектор Y, при некоторых значениях аргумента, которые образуют вектор Х такой же длины, что и вектор Y, осуществляется процедурой polyfit(X, Y, n). Здесь n - порядок аппроксимирующего полинома. Результатом действия этой процедуры является вектор длиной (n +1) из коэффициентов аппроксимирующего полинома.

Пусть массив значений аргумента имеет вид:

а массив соответствующих значений измеренной величины - вид:

Тогда, применяя указанную функцию при разных значениях порядка аппроксимирующего полинома, получим:

» y = [ -1.1 0.2 0.5 0.8 0.7 0.6 0.4 0.1];

» polyfit(x,y,1) » polyfit(x,y,2) » polyfit(x,y,3) » polyfit(x,y,4) Это означает, що заданную зависимость можно аппроксимировать или прямой или квадратной параболой или кубической параболой или параболой четвертой степени Построим в одном графическом поле графики заданной дискретной функции и графики всех полученных при аппроксимации полиномов:

y = [ -1.1 0.2 0.5 0.8 0.7 0.6 0.4 0.1];

P1=polyfit(x,y,1) ;

P2=polyfit(x,y,2);

P3=polyfit(x,y,3);

P4=polyfit(x,y,4) ;

stem(x,y);

y1=polyval(P1,x1);

y2=polyval(P2,x1);

y3=polyval(P3,x1);

y4=polyval(P4,x1);

plot(x1,y1,x1,y2,x1,y3,x1,y4), grid, set(gca, 'FontName', 'Arial Cyr', 'FontSize', 16), title('Полиномиальная аппроксимация ');

xlabel('Аргумент');

ylabel('Функция') Результат представлен на рис. 1.18.

Функция spline(X,Y,Xi) осуществляет интерполяцию кубическими сплайнами. При обращении Yi = spline(X,Y,Xi) она интерполирует значение вектора Y, заданного при значениях аргумента, представленных в векторе Х, и выдает значение интерполирующей функции в виде вектора Yi при значениях аргумента, заданных вектором Xi. В случае, если вектор Х не указан, по умолчанию принимается, что он имеет длину вектора Y и любой его элемент равен номеру этого элемента.

В качестве примера рассмотрим интерполяцию вектора x = -0.5:0.1:0.2;

y = [-1.1 0.2 0.5 0.8 0.7 0.6 0.4 0.1];

x1 = -0.5:0.01:0.2;

y2 = spline(x,y,x1);

plot (x,y,x1,y2), grid set(gca,'FontName','Arial Cyr','FontSize',16), title('Интерполяция процедурой SPLINE ');

xlabel('Аргумент');

ylabel('Функция') Результат приведен на рис. 1.19.

Одномерную табличную интерполяцию осуществляет процедура interp1.

Обращение к ней в общем случае имеет вид:

Yi = interp1(X,Y,Xi,’’), и позволяет дополнительно указать метод интерполяции в четвертом входном аргументе:

'nearest' - ступенчатая интерполяция;

'linear' - линейная;

‘cubic' - кубическая;

‘spline' - кубическими сплайнами.

Если метод не указан, осуществляется по умолчанию линейная интерполяция. Например, (для одного и того же вектора):

x = -0.5:0.1:0.2;

y = [-1.1 0.2 0.5 0.8 0.7 0.6 0.4 0.1];

x1 = -0.5:0.01:0.2;

y1 = interp1(x,y,x1);

y4 = interp1(x,y,x1,'nearest');

y2 = interp1(x,y,x1,'cubic');

y3 = interp1(x,y,x1,'spline');

plot (x1,y1,x1,y2,x1,y3,x1,y4), grid plot (x1,y1,x1,y2,'.',x1,y3,'--',x1,y4,':'), grid set(gca,'FontName','Arial Cyr','FontSize',8), legend('линейная','кубическая','сплайновая','ступенчатая',0) set(gca,'FontName','Arial Cyr','FontSize',14), set(gcf,'color','white') title('Интерполяция процедурой INTERP1 ');

xlabel('Аргумент');

ylabel('Функция') Результаты приведены на рис.1.20.

1.4.5. Векторная фильтрация и спектральный анализ В системе MatLAB есть несколько функций для проведения цифрового анализа данных наблюдений (измерений).

Так, функция y = filter(b,a,x) обеспечивает формирование вектора y по заданным векторам b, a, x в соответствии с соотношением:

y(k) = b(1)*x(k) + b(2)*x(k-1) +... + b(nb+1)*x(k-nb) a(2)*y(k-1) - a(3)*y(k-3) -... - a(na+1)*y(k-na), (1.1) где вектор b имеет такой состав a вектор а Соотношение (1) можно рассматривать как конечно-разностное уравнение фильтра с дискретной передаточной функцией вида рациональной дроби, коэффициенты числителя которого образовывают вектор b, а знаменателя - вектор а, на вход которого подается сигнал x(t), а на выходе формируется сигнал y(t).

Тогда вектор y будет представлять собой значение исходного сигнала этого фильтра в дискретные моменты времени, соответствующие заданным значениям входного сигнала x(t) (вектор х).

Ниже приведен пример применения функции filter.

» y = filter(b,a,x) Columns 8 through 0.2503 -13.4768 2.8466 56. Функции fft (Fast Fourier Transformation) и ifft (Invers Fast Fourier Transformation) осуществляют преобразование заданного вектора, сооствествующие дискретному прямому и обратному преобразованиям Фурье.

Обращение к этим функциям вида:

приводит к формированию вектора y в первом случае и х - в втором по формулам:

где j - обозначение мнимой единицы; n - число элементов заданного вектора х (оно представляет также размер выходного вектора y).

Приведем пример. Сформируем входной сигнал в виде вектора, элементы которого равняются значениям функции, являющейся суммой двух синусоид с частотами 5 и 12 Гц (рис. 1.21).

t =0:0.001:2;

x = sin(2*pi*5*t) + cos(2*pi*12*t);

plot(t, x); grid set(gcf,'color','white') set(gca,'FontName','Arial Cyr','FontSize',16), title('Входной процесс ');

xlabel('Время (с)');

ylabel('Х(t)') Найдем Фурье-изображение этого сигнала и выведем графическое представление модуля його Фурье-изображения:

plot(a); grid set(gca,'FontName','Arial Cyr','FontSize',16), title('Модуль Фурьє - изображения ');

xlabel('Номер элемента вектора');

ylabel('abs(F(X(t))') Результат отображен на рис. 1.22.

Теперь осуществим обратное преобразование с помощью функции ifft и результат также выведем в форме графика:

plot(t, z); grid set(gca,'FontName','Arial Cyr','FontSize',16), title('Обратное Фурье-преобразование ');

xlabel('Время (с)');

ylabel('Z(t)') На рис. 1.23 изображен результат. Рассматривая его, можно убедиться, что воспроизведенный процесс точно совпадает с исходным.

Внимательно изучая формулу дискретного преобразования Фурье, можно заметить:

а) номер m отвечает моменту времени tm, в который измерен входной сигнал x(m); при этом t1 = 0 ;

б) номер k - это индекс значения частоты fk, которому отвечает найденный элемент y(k) дискретного преобразования Фурье;

в) чтобы перейти от индексов к временной и частотной областям, необходимо знать значение h дискрета (шага) времени, через который измерен входной сигнал x(t) и промежуток T времени, на протяжении которого он измеряется; тогда шаг (дискрет) по частоте в изображении Фурье определится соотношением:

а диапазон изменения частоты - формулой так, в анализируемом примере (h =0. 001, T = 2, n = 21):

г) из (2) следует, что индексу k = 1 отвечает нулевое значение частоты (f0 = 0); иначе говоря, первый элемент вектора y(1) является значением Фурьеизображения при нулевой частоте, т. е. - просто суммой всех заданных значений вектора x; отсюда получаем, что вектор y(k) содержит значение Фурьеизображения, начиная из частоты f0 = 0 (которой отвечает k = 1) до максимальной частоты fmax = F (которой отвечает k = n); таким образом, Фурье-изображение определяется функцией fft только для положительных частот в диапазоне от 0 к F; это неудобно для построения графиков Фурье-изображения от частоты; более удобным и привычным представляется переход к вектору Фурье-изображения, определенному в диапазоне частот от (-F/2) до F/2; частота FN = F / 2 получила название частоты Найквиста;

д) как известно, функция e jz является периодической по z с периодом 2 ; поэтому информация об Фурье-изображении при отрицательных частотах расположена во второй половине вектора y(k).

Сформируем для анализируемого примера массив частот, исходя из вышеупомянутого:

и выведем графік с аргументом-частотой (рис. 1.24):

plot(f,a); grid set(gca,'FontName','Arial Cyr','FontSize',16), title('Модуль Фурье - изображения ');

xlabel('Частота (Гц)');

ylabel('abs(F(X(t))') Как следует из рассмотрения рис. 1.24, по нему непросто распознать те частоты (5 и 12 Гц), с которыми изменяется входной сигнал. Это - следствие того обстоятельства, которое было отмечено в примечании г). Чтобы определить частотный спектр входного сигнала, нужно сначала преобразовать полученный вектор y Фурье-изображения с помощью процедуры fftshift.

Функция fftshift (обращение к ней осуществляется таким образом: z = fftshift(y)) предназначена для формирования нового вектора z из заданного вектора y путем перестановки второй половины вектора y в первую половину вектора z. При этом вторая половина вектора z состоит из элементов первой половины вектора y. Более точно эту операцию можно задать соотношениями:

z(1) = y(n/2+1);..., z(k) = y(n/2+k);..., z(n/2) = y(n); z(n/2+1) = y(1);...

Примечание. Операцию fftshift удобно использовать для преобразования массива Фурье-изображение с целью построения его графика в частотной области. Тем не менее этот массив не может быть использован для обратного преобразования Фурье.

Проиллюстрируем применение этой функции к предыдущему примеру:

f1 = -500 : 0.5 : 500; % Перестройка вектора частот v = fftshift(y); % Перестройка вектора Фурье-изображения a =abs(v); % Отыскание модуля plot(f1(970:1030),a(970:1030)); grid set(gca,'FontName','Arial Cyr','FontSize',16), title('Модуль Фурье - изображения');

xlabel('Частота (Гц)'); ylabel('abs(F(X(t))') Из графика рис. 1.25 уже становится очевидным, что в спектре входного сигнала есть две гармоники - с частотами 5 и 12 Гц.

Остается лишь то неудобство, что из графика спектра невозможно установить амплитуды этих гармоник. Во избежание этого, нужно весь вектор y Фурьеизображения разделить на число его элементов (n), чтобы получить вектор комплексного спектра сигнала:

N=length(y); a =abs(v)/N;

plot(f1(970:1030),a(970:1030)); grid set(gca,'FontName','Arial Cyr','FontSize',16,'Color','white'), title('Модуль комплексного спектра');

xlabel('Частота (Гц)');

ylabel('abs(F(X(t)) / N') Результат приведен на рис. 1.26. Из его рассмотрения вытекает, что "амплитуды" всех составляющих гармоник равны 0.5. При этом нужно принять во внимание, что "амплитуды" распределены между положительными и отрицательными частотами поровну, поэтому они вдвое меньшее действительной амплитуды соответствующей гармоники.

1.4.6. Задания Задание 1.6.

1. Введите произвольную матрицу размером (4*6). Найдите сумму наибольших элементов ее строк.

2. Введите квадратную матрицу (5*5) с одним наименьшим элементом.

Найдите сумму элементов строки, в которой размещен элемент с наименьшим значением.

3. Введите матрицу (6*9), в которой есть единственные наибольший и наименьшие элементы и они расположены в разных строках. Поменяйте местами строку с наибольшим элементом и строку с наименьшим элементом.

4. Введите матрицу (5*6) с разными значениями элементов. В каждой строке выберите элемент с наименьшим значением, из полученных чисел выберите наибольшее. Найдите индексы полученных элементов.

5. Введите матрицу (5*6). Найдите вектор, элементами которого являются наибольшие элементы соответствующей строки матрицы.

6. Введите матрицу (5*6). Постройте вектор, элементами которого являются суммы наибольшего и наименьшего элементов соответствующей строки матрицы.

7. Введите матрицу (5*6). Постройте вектор, элементами которого являются средние значения элементов соответствующей строки матрицы.

8. Введите матрицу (5*6). Постройте вектор, элементами которого являются среднеквадратичные отклонения элементов соответствующей строки матрицы от их среднего значения.

9. Введите матрицу (5*6). Постройте вектор, элементами которого являются средние арифметические наибольшего и наименьшего элементов соответствующей строки матрицы.

10. Введите матрицу (6*5). Постройте вектор, элементами которого являются суммы квадратов элементов соответствующего столбца матрицы.

11. Введите матрицу (5*5). Постройте векторы, элементами которых являются суммы элементов столбцов матрицы, произведения элементов столбцов и наименьшие значения элементов столбцов.

12. Введите матрицу (5*6). Найдите среднее арифметическое наибольших и наименьших ее элементов.

13. Введите матрицу (5*5). Постройте вектор, элементами которого являются элементы главной диагонали матрицы. Найдите след матрицы.

14. Введите две матрицы (4*4). Постройте новую матрицу размером (4*8), включая в первые 4 столбца строки первой матрицы, а в другие - столбцы второй матрицы.

15. Найдите сумму всех элементов матрицы размером (4*3).

Задание 1.7. Вычислите векторы:

а) модуля частотной передаточной функции (ЧПФ);

б) аргумента ЧПФ;

в) действительной части ЧПФ;

г) мнимой части ЧПФ по заданным числителю и знаменателю передаточной функции (таблица 1.4).

Предварительно найдите корни знаменателя передаточной функции, определите наибольшую собственную частоту max системы. Обеспечьте вычисление ЧПФ при 100 значениях частоты в диапазоне от 0 до 5max.

4. 61p2+1. 82p+67.56 p4+3. 65p3+45p2+7. 04p+ p3+4. 61p2+1. 82p+67.56 p4+2. 65p3+68p2+5p+ p3+4. 61p2+1. 82p+67.56 p3+3. 09p2+70p+ p3+4. 61p2+1. 82p+67.56 7p2+7p+ p3+4. 61p2+1. 82p+67.56 p4+5p3+30p2+7p+ Указание. Частотной передаточной функцией называют передаточную функцию системы при мнимых значениях її аргумента ( p = j ).

Собственные частоты системы - это значения модулей мнимых частей корней характеристического уравнения системы (которое получается приравниванием нулю знаменателя передаточной функции).

Задача 1.8. Введите произвольную матрицу размером (5*5). Найдите:

1) определитель матрицы; в случае, если определитель окажется равным нулю, или слишком малым, измените некоторые элементы матрицы и повторите вычисления;

2) обратную матрицу; проверьте правильность путем обращения обратной 3) характеристический полином матрицы;

4) корни характеристического полинома матрицы; рассортируйте корни по комплексно-спряженным парам и в порядке возрастания величин;

5) собственные значения матрицы; сравните с ранее найденными корнями характеристического полинома;

6) LU-разложение матрицы; проверьте его правильность;

7) QR-разложение матрицы; проверьте его правильность;

8) сингулярные числа матрицы; сравните их с получаемыми при svdразложении;

9) след матрицы;

10) число обусловленности матрицы;

11) экспоненту от матрицы;

12) логарифм от экспоненты матрицы; сравните с исходной матрицей.

1.4.7. Вопросы 1. Какой объект в MatLAB называется полиномом?

2. Как в MatLAB осуществляется перемножение и деление полиномов?

3. При помощи каких функций можно найти корни заданного полинома, значение полинома по известному значению аргумента?

4. Какие функции позволяют найти производную от полинома?

5. Как найти характеристический полином матрицы?

1.5. Построение простейших графиков 1.5.1. Процедура plot Вывод графиков в системе MatLAB - настолько простая и удобная операция, что ее можно использовать даже в режиме калькулятора.

Основной функцией, обеспечивающей построение графиков на экране дисплея, является функция plot. Общая форма обращения к ней такова:

plot(x1, y1, s1, x2, y2, s2,...).

Здесь x1, y1 - заданные векторы, элементами которых являются массивы значений аргумента (х1) и функции (y1), отвечающие первой кривой графика; x2, y2 - массивы значений аргумента и функции второй кривой и т.д. При этом предполагается, что значения аргумента откладываются вдоль горизонтальной оси графика, а значения функции - вдоль вертикальной оси. Переменные s1, s2,... являются символьными (их указание не является обязательным). Любая из них может содержать до трех специальных символов, определяющих соответственно: а) тип линии, которая соединяет отдельные точки графика; б) тип точки графика; в) цвет линии. Если переменные s не указаны, то тип линии по умолчанию - отрезок прямой, тип точки - пиксел, а цвет устанавливается (в версии 5.3) в такой очередности:

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

Например, обращение вида plot(x1,y1,x2,y2,...) приведет к построению графика, в котором первая кривая будет линией из отрезков прямых синего цвета, вторая кривая - такого же типа зеленой линией и т.д.

Графики в MatLAB всегда выводятся в отдельное графическое окно, которое называют фигурой.

Приведем пример. Пусть нужно вывести график функции на промежутке от -3 до +3 с шагом /100.

Сначала надо сформировать массив значений аргумента х:

потом вычислить массив соответствующих значений функции:

и, наконец, построить график зависимости y(х).

В командном окне последовательность операций будет выглядеть так:

» x = -3*pi:pi/100:3*pi;

» y = 3*sin(x+pi/3);

» plot(x,y) В результате на экране появится окно с графиком (см. рис. 1.27).

Если вектор аргумента при обращении к функции plot не указан явно, то система по умолчанию принимает в качестве аргумента номера элементов вектора функции. Например, если ввести команду то результатом будет появление графика в виде, приведенном на рис. 1.28.

Графики, приведенные на рис. 1.27, 1.28, имеют несколько недостатков:

на них не нанесена сетка из координатных линий, что затрудняет "чтение" графиков;

нет общей информации о кривой графика (заголовка);

неизвестно, какие величины отложены по осям графика.

Первый недостаток устраняется с помощью функции grid. Если к этой функции обратиться сразу после обращения к функции plot:

» x = -3*pi:pi/100:3*pi;

» y = 3*sin(x+pi/3);

» plot(x,y), grid, то график будет снабжен координатной сеткой (рис. 1.29).

Ценной особенностью графиков, построенных в системе MatLAB, является то, что сетка координат всегда соответствует "целым" шагам изменения, что делает графики "читабельными", т. е. такими, что по графику можно отсчитывать значение функции при любом заданном значении аргумента и наоборот.

Заголовок графика выводится с помощью процедуры title. Если после обращения к процедуре plot вызвать title таким образом:

то над графиком появится текст, записанный между апострофами в скобках. При этом следует помнить, что текст всегда должен помещаться в апострофы.

Аналогично можно вывести объяснения к графику, которые размещаются вдоль горизонтальной оси (функция xlabel) и вдоль вертикальной оси (функция ylabel).

Например, совокупность операторов » x = -3*pi : pi/100 : 3*pi;

» y = 3*sin(x+pi/3);

» plot(x,y), grid » title('Функция y = 3*sin(x+pi/3)');

» xlabel('x'); ylabel('y');

приведет к оформлению поля фигуры в виде, представленном на рис. 1.30. Очевидно, такая форма уже целиком удовлетворяет требованиям, предъявляемым к инженерным графикам.

Не сложнее вывод в среде MatLAB графиков функций, заданных параметрически. Пусть, например, необходимо построить график функции y(х), заданной формулами:

Параметрическая функция x=4*exp(-0.05t)*sin(t); y= 0.2*exp(-0.1t)*sin(2t) Выберем диапазон изменения параметра t от 0 до 50 с шагом 0.1. Тогда, набирая совокупность операторов x = 4*exp(-0.05*t).*sin(t);

y = 0.2*exp(-0.1*t).*sin(2*t);

set(gcf,'color','white') title('Параметрическая функция x=4*exp(-0.05t)*sin(t); y= 0.2*exp(-0.1t)*sin(2t) ') grid, получим график рис. 1.31.

1.5.2. Специальные графики Большим удобством, предоставляемым системой MatLAB, является указанная ранее возможность не указывать аргумент функции при построении ее графика. В этом случае в качестве аргумента система принимает номер элемента вектора, график которого строится. Пользуясь этим, например, можно построить "график вектора":

» title('График вектора Х') » ylabel('Значение элементов') » xlabel(' Номер элемента').

Результат представлен на рис. 1.32.

Еще более наглядным является представление вектора в виде столбцовой диаграммы с помощью функции bar (см. рис. 1.33):