Добавил:
Upload Опубликованный материал нарушает ваши авторские права? Сообщите нам.
Вуз: Предмет: Файл:
лабы2 2 семmath_pack.doc
Скачиваний:
5
Добавлен:
01.04.2025
Размер:
2 Мб
Скачать
☆

5.1. Теоретическое введение

ОПЕРАЦИИ С ВЕКТОРАМИ И МАТРИЦАМИ

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

eye(m, n) — возвращают единичную матрицу размера m x n с единицами по диагонали и нулями в остальных ячейках;

ones(m, n) — возвращают матрицу размера m x n, состоящую из единиц;

zeros(m, n) — возвращают матрицу размера m x n, состоящую из нулей;

linspace(a, b, n) — генерирует n точек, равномерно распределенных в интервале от а до b;

rand(m, n) — генерирует матрицу случайных чисел размера m x n из промежутка (0, 1), используя равномерный закон распределения; rand () без аргументов возвращает одно случайное число;

В приведенных выше функциях (кроме linspace( n)) в качестве аргумента допустимо использовать size(A), что позволит формировать результирующую матрицу того же размера и размерности, что и А.

cat (dim, А, В) — реализует конкатенацию (объединение) массивов А и В в соответствии со спецификацией размерности dim и возвращает объединенный массив. При этом, при dim = 1 происходит горизонтальная конкатенация, при dim = 2 — вертикальная.

X = diag(v) — помещает вектор v на главную диагональ квадратной матрицы X и нули в остальные ячейки;

v = diag(X) — возвращает главную диагональ матрицы X;

repmat (А, m, n) — возвращает матрицу В размера m x n, в которой каждый элемент заменяется на копию матрицы А); если А — скаляр, то функция возвращает матрицу размера m x n со значениями элементов, заданных А;

tril(X, k) — возвращает матрицу, все элементы которой, расположенные выше k-й диагонали, заменены нулями, а нижняя треугольная часть матрицы X, включая указанную диагональ, остается неизменной. При k=0 это главная диагональ, при k>0 — одна из верхних диагоналей, при k<0 — одна из нижних диагоналей. По умолчанию к = 0.

triu(X. k) — возвращает неизменной верхнюю треугольную часть матрицы X, начиная с k-й диагонали, остальные элементы заменяются нулями. При k=0 — это главная диагональ, при k>0 — одна из верхних диагоналей, при k<0 — одна из нижних диагоналей;

magic(n) – возвращает так называемый магический квадрат – квадратную матрицу размера n x n, в которой сумма значений элементов по строкам, столбцам, а также сумма элементов главной и побочной диагонали равны.

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

prod (A, dim) — возвращает произведение элементов массива А по столбцам при dim=l или по строкам при dim=2;

sum(A, dim) — возвращает сумму элементов массива A по столбцам (dim=1) или по строкам (dim=2).

cross(U, V) — возвращает векторное произведение векторов U и V в трехмерном пространстве;

dot(U, V) — возвращает скалярное произведение векторов U и V;

transpose(X) – транспонирование вектора/матрицы Х;

trace(X) – возвращает след матрицы Х, т. е. сумму ее диагональных элементов;

isequal(U ,V) – проверка равенства векторов/матриц U и V;

norm(X, dim) – возвращает норму вектора Х, при dim=1 – первую норму ( ), при dim=2 – евклидову (вторую) норму ( ), при dim=3 – бесконечную норму( );

det(X) – возвращает определитель матрицы X;

inv(X) – возвращает матрицу, обратную к матрице X;

rank (X) – возвращает ранг матрицы X.

Весьма представителен в MATLAB набор матричных функций. Перечислим некоторые из них:

ехрт(Х) — возвращает от матрицы X;

sqrtm(X) – возвращает логарифм матрицы;

logm(X) – возвращает квадратный корень из X, соответствующий неотрицательным действительным частям собственных значений X.

ЧИСЛЕННЫЕ МЕТОДЫ

Элементарные средства решения систем линейных алгебраических уравнений

При решении инженерных задач достаточно часто приходится сталкиваться с решением систем линейных алгебраических уравнений (СЛАУ). Как известно, обычная СЛАУ может быть представлена в матричном виде как АХ=В, где А — матрица коэффициентов уравнений размера m x n, X — искомый вектор неизвестных с n компонентами и В — вектор свободных членов, также содержащий n компонент. Решением системы называется совокупность чисел , при подстановке которых все уравнения системы обращаются в тождества.

СЛАУ называется совместной, если имеет хотя бы одно решение, и несовместной, если не имеет ни одного решения. Таким образом, решить СЛАУ – значит определить, является ли она совместной или нет.

Существуют точные и приближенные (итерационные) методы решения СЛАУ. К точным относятся методы, которые дают точное решение задачи после конечного числа арифметических и логических операций. К указанным методам относятся, например, матричный метод, метод Крамера, метод Гаусса и т. д. К приближенным – метод Якоби, метод Гаусса – Зейделя и др.

Для определения совместности системы можно использовать теорему Кронекера-Капелли, согласно которой для того, чтобы СЛАУ была совместной, необходимо и достаточно, чтобы ранг матрицы А системы был равен рангу ее расширенной матрицы. Напомним, что расширенная матрица получается после дополнения матрицы коэффициентов столбцом свободных членов. Если ранг совместной системы равен числу неизвестных, то система имеет единственное решение. Если ранг совместной системы меньше числа неизвестных, то система имеет бесчисленное множество решение.

На рис. 5.1 представлен пример исследования на совместность следующей системы уравнений:

Рис. 5.1. Пример исследования на совместность заданной СЛАУ.

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

Рассмотрим точные методы решения СЛАУ в системе MATLAB. Одним из таких методов является матричный метод решения СЛАУ, который предусматривает выполнение следующих шагов:

  1. Убедиться в том, что система не является вырожденной (определитель ее не равен нулю), т. е. существует обратная матрица.

  2. Найти обратную матрицу .

  3. Найти вектор – столбец X, вычислив произведение .

На рис. 5.2 приводится пример решения представленной выше системы уравнений с помощью матричного метода.

Рис. 5.2. Пример решения СЛАУ с помощью матричного метода.

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

В основе решения СЛАУ с помощью формул Крамера лежит теорема Крамера: система n уравнений с n неизвестными, определитель которой отличен от нуля всегда имеет решение, и притом единственное. Значение каждого из неизвестных в этом случае равно дроби, знаменателем которой является определитель системы, а числитель равен определителю системы, в которой вместо столбца коэффициентов при неизвестном располагается столбец свободных членов.

Пусть дана система n линейных уравнений с n неизвестными:

Обозначим определитель системы: . Заменим поочередно столбцы коэффициентов при на столбец свободных членов и получим n определителей:

.

Формулы Крамера для решения системы уравнений имеют вид:

.

Решение СЛАУ из предыдущего примера с помощью формул Крамера представлено на рис. 5.3.

Рис. 5.3. Пример решения СЛАУ с помощью формул Крамера.

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

Существо метода Гаусса состоит в последовательном исключении неизвестных. Систему из n уравнений сначала приводят к равносильной ей системе с треугольной матрицей коэффициентов (прямой ход). Затем выполняют обратный ход, вычисляя значения и последовательно подставляя их в уравнения для получения значений . Для уменьшения влияния ошибок округления на решение СЛАУ указанным методом применяют его модификацию – метод Гаусса с частичным выбором главного элемента.

Этот метод в MATLAB реализован в функции rref(), которая на его основании осуществляет приведение расширенной матрицы к треугольной форме. Пример ее применения для решения СЛАУ представлен на рис. 5.4.

Рис. 5.4. Пример решения СЛАУ с помощью метод Гаусса с частичным выбором главного элемента.

Как было сказано выше, помимо точных методов решения СЛАУ существует ряд итерационных методов. Одним их них является метод наименьших квадратов, реализованный в системе MATLAB функцией lsqr(). При вызове функции lsqr(A, B) итерации производятся либо до сходимости к решению, либо до появления ошибки, либо до достижения максимального числа итераций (по умолчанию равного min(20, m, n), где m – число уравнений, n – число неизвестных). Сходимость достигается, когда отношение вторых норм векторов norm(B-Ax)/norm(B) меньше или равно погрешности метода tol (по умолчанию 1е-6). Функция lsqr(A, B, tol) — возвращает решение с заданной погрешностью tol.

На рис. 5.5 приведен пример решения СЛАУ с использованием функции lsqr(A, B). В данном примере процесс итераций сходится на четвертом шаге с относительным остатком (отношением вторых норм векторов невязки и свободных членов) 1.7∙10-11.

Рис. 5.5. Пример решения СЛАУ методом наименьших квадратов.

Помимо рассмотренных методов решения СЛАУ в MATLAB реализован еще ряд методов решения с помощью функций, расположенных в ядре системы. Описание этих функций можно посмотреть с помощью команды help имя_функции. Функция mldivide() или оператор «\» самостоятельно выбирают лучший метод решения заданной системы уравнений (см. рис. 5.6).

Рис. 5.6. Пример решения СЛАУ.

Решение нелинейных уравнений и их систем

Хотя в ядро MATLAB последних версий встроено около 900 операторов и функций, пользователю всегда может потребоваться та или иная функция, отсутствующая в ядре. Язык программирования системы MATLAB предоставляет ряд возможностей для задания функции пользователя. Одна из таких возможностей – задание своих функций в виде m– файлов, для создания которых может использоваться как встроенный редактор Editor/ Debugger (открывается командой New из меню File), так и любой текстовый редактор, поддерживающий формат ASCII. Подготовленный и записанный на диск m-файл становится частью системы, и его можно вызывать как из командной строки, так и из другого m-файла.

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

  • ключевые слова языка программирования — синий цвет;

  • операторы, константы и переменные — черный цвет;

  • комментарии после знака % — зеленый цвет;

  • символьные переменные (в апострофах) — зеленый цвет;

  • синтаксические ошибки — красный цвет.

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

М-файлы, создаваемые редактором/отладчиком, делятся на два класса: файлы-сценарии или Script-файлы, не имеющие входных параметров, и файлы-функции, имеющие входные параметры. Script-файл использует глобальные переменные, значения которых могут быть изменены в любой момент сеанса работы и в любом месте программы. Для запуска файла-сценария из командной строки MATLAB достаточно просто указать его имя в этой строке. Файл-функция отличается от файла-сценария прежде всего тем, что созданная им функция имеет входные и выходные параметры.

Строка заголовка функции начинается с объявления function, после которого указывается имя выходного параметра, заключенное в квадратные скобки. Затем следует знак равенства, имя самой функции и заключенный в круглые скобки список ее входных параметров. После заголовка идет тело функции, в котором задаются инструкции необходимые для вычисления значений выходного аргумента на основе входных данных. Например,

function [f] = fun(x)

f = 0.25 * x +cos(x)-1;

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

>>fun(0)

>>x = fun(pi)

Еще одна возможность для создания функции пользователя заключается в применении объектов класса дескрипторов функций, задаваемых с помощью символа @, например: f=@exp. Численные значения функций, заданных дескрипторами, вычисляются с помощью функции feval(), например, feval(f,1,0).

Довольно часто возникает задача решения нелинейного уравнения вида или . Последнее, однако, можно свести к виду . Таким образом, данная задача сводится к нахождению значений аргумента х функции f(x) одной переменной, при котором значение функции равно нулю. В системе MATLAB имеется ряд функций для решения данной задачи, например, функция fzero(). В зависимости от формы задания функции fzero() реализуются следующие хорошо известные численные методы поиска нуля функции: деления отрезка пополам, секущей и обратной квадратичной интерполяции. Приведем некоторые из них:

fzero( @fun, x0, options) — возвращает уточненное значение начального приближения х0, при котором достигается нуль функции fun, используя параметры вычислений из вектора options (не является обязательным параметром): tolX – конечное допустимое отклонение по значению х, maxfunevals – максимально число допустимых расчетов функции, maxiter – максимальное число допустимых итераций, display – уровень отображения ('off' отображение не производится, 'iter' отображение проводится на каждой итерации, 'final' (принимается по умолчанию) отображение только конечной информации). Перечисленные параметры следует предварительно установить при помощи команды optimset(). Например, options=optimset('Display','iter'). Возвращенное значение близко к точке, где функция меняет знак, или равно NaN, если такая точка не найдена;

fzero(@fun,[xl x2], options) — возвращает значение х, при котором fun(x)=0, с заданием интервала поиска с помощью вектора x=[xl х2], такого, что знак fun(xl) отличается от знака fun(x2). Если это не так, выдается сообщение об ошибке. Вызов функции fzero с интервалом гарантирует, что fzero возвратит значение, близкое к точке, где fun изменяет знак.

Приведенный ниже пример показывает приближенное вычисление корня уравнения cos(x)=0 с представлением косинуса дескриптором:

» х= fzero(@cos, [1 3])

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

Рассмотрим еще один пример решения нелинейного уравнения . На первом этапе решения создадим m-файл fun.m, листинг которого представлен ниже:

function [ f]=fun(x)

f=exp(x)-2+x;

Далее построим график функции fun, используя следующую последовательность команд:

>>х=-5:0.1:5;

>>plot(x,funl(x));grid on;

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

Рис. 5.7. График функции f=exp(x)-2+x.

Из рисунка нетрудно заметить, что значения корня заключено в интервале [-1, 1]. Найдем его, используя функцию fzero(), предварительно задав точность нахождения корня:

>> options=optimset('Display','final', 'tolX', 10^(-6));

>>x=fzero( @fun, [-1 1], options)

В результате получим корень уравнения х = 0.4429.

Для решения систем нелинейных уравнений предназначена функция fsolve() из пакета Optimization Toolbox, которая использует для решения метод наименьших квадратов и, в отличие от функции fzero(), ищет не только точки пересечения, но и точки касания графиков функций. Синтаксис функции:

fsolve( @fun, x0, options) — возвращает уточненное значение начального приближения х0, при котором достигаются нули функций описанных в fun, используя параметры вычислений из вектора options;

fsolve (@fun,[xl x2], options) — возвращает значение х, при котором достигаются нули функций описанных в fun, с заданием интервала поиска с помощью вектора x=[xl х2].

Применим функцию fsolve() для решения следующей системы нелинейных уравнений: .

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

>> x=-1:0.1:1;

>> y1=acos(sqrt(2/3).*cos(x));

>> y2=asin(sqrt(2).*sin(x));

>> plot(x,y1,x,y2)

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

Рис. 5.8. Графическое решение системы уравнений.

Далее создадим m-файл fun.m, листинг которого представлен ниже:

function [ f]=fun1_(x)

f=[acos(sqrt(2/3)*cos(x(1)))-x(2); asin(sqrt(2)*sin(x(1)))-x(2)];

Найдем решение заданной системы нелинейных уравнений с точностью , используя функцию fsolve():

x=fsolve (@fun,[0 0],optimset('tolX', 10^(-6)'))

В результате получим следующее решение:

x=[0.5236 0.7854]