Некоторые функции для работы с массивами.

Вообще, существует множество различных функций, специально предназначенных для работы с массивами. Это, в частности, уже рассмотренные ранее функции ones, zeros, sum, cat, size, ndims, reshape и многие другие. Часть перечисленных функций сообщает служебную информацию о массивах, другая группа функций обеспечивает контролируемое (управляемое) изменение их структуры, третья группа предназначена для генерации массивов с заданными свойствами и т.д. К последней группе можно отнести, в частности, рассмотренные в параграфе 2.15 функции eye и rand.

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

>> A=[2 2 2; 3 3 3; 5 5 5]

A =

2 2 2

3 3 3

5 5 5

>> prod(A)

ans =

30 30 30

Как отмечалось, по умолчанию суммирование и перемножение элементов массива с помощью функций sum и prod выполняется по столбцам. Для того, чтобы выполнить соответствующие операции по строкам, указанные функции следует вызвать с двумя входными аргументами. Первым из аргументов задается имя матрицы, а вторым – номер индекса, по которому требуется выполнить суммирование или перемножение: 1 – по столбцам (по умолчанию), 2 – по строкам. Рассмотрим пример:

>> A=[1 2; 3 4]

A =

1 2

3 4

>> sum(A,2)

ans =

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

>> x=[3 7 5 4.8 8.95 2 6];

>> min(x)

ans =

>> max(x)

ans =

8.9500

>> A=[1 2 3; 4 5 6; 7 8 9]

A =

1 2 3

4 5 6

7 8 9

>> max(A)

ans =

7 8 9

>> min(A)

ans =

1 2 3

Если требуется определить максимальные и минимальные элементы не в столбцах, а в строках заданной матрицы, форма вызова функций max и min должна быть следующей:

>> max(A,[],2)

ans =

>> min(A,[],2)

ans =

Заметим, что результат представлен в виде вектор-столбца.

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

>> b=[3 2 4 8 5 7 3 1 0]

b =

3 2 4 8 5 7 3 1 0

>> sort(b)

ans =

0 1 2 3 3 4 5 7 8

>> B=[3 7 9; 4 8 1; 2 9 5]

B =

3 7 9

4 8 1

2 9 5

>> sort(B,2)

ans =

3 7 9

1 4 8

2 5 9

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

>> [B,ind]=sort(B,2)

B =

3 7 9

1 4 8

2 5 9

ind =

1 2 3

3 1 2

1 3 2

Часть 4.

СТАНДАРТНЫЕ СРЕДСТВА РЕШЕНИЯ НЕКОТОРЫХ ТИПОВЫХ ЗАДАЧ ЛИНЕЙНОЙ АЛГЕБРЫ И МАТЕМАТИЧЕСКОГО АНАЛИЗА

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

Как уже отмечалось, в MATLAB для решения систем линейных алгебраических уравнений (СЛАУ) предусмотрены знаки операций «/» и «\». Второй из этих знаков уже использовался в предыдущем параграфе именно для решения СЛАУ.

В общем случае операция «\» называется «левым делением матриц» и, будучи примененная к матрицам A и B в виде A\B, подобна вычислению выражения inv(A)*B, где под inv(A) понимается вычисление матрицы, обратной к матрице A.

Операция «/» называется «правым делением матриц». Выражение A/B подобно вычислению выражения B*inv(A), т.е. данная операция позволяет, в частности, решать СЛАУ вида Y*A=B.

Решение некоторых задач линейной алгебры.

Матричные функции

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

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

cond(A) – возвращает число обусловленности матрицы как отношение наибольшего сингулярного числа матрицы A к наименьшему (если соответствующий результат близок к единице, то матрица хорошо обусловлена);

condeig(A) – возвращает вектор, содержащий числа обусловленности собственных значений матрицы A;

rcond(A) – возвращает обратную величину числа обусловленности матрицы A (значения, близкие к единице, свидетельствуют о хорошей обусловленности, а близкие к нулю – о плохой);

det(A) – возвращает значение определителя (детерминанта) матрицы A;

rank(A) – возвращает значение ранга матрицы A;

norm(A) – возвращает норму матрицы A (наибольшее по абсолютной

величине сингулярное число);

trace(A) – возвращает след матрицы, который равен сумме значений элементов матрицы, расположенных на главной диагонали;

orth(A) – возвращает ортонормированный базис для разложения матрицы A [20,21,27,40,45,46,51,95,96,117-120,134,135,238];

rref(A) – приводит расширенную матрицу СЛАУ n-го порядка к виду [E v], где E – единичная матрица n-го порядка, v – вектор-столбец, содержащий решение СЛАУ;

rrefmovie(A) – результат аналогичен результату, возвращаемому функцией rref, но при этом выводится информация о каждой выполняемой операции и получаемом при этом результате;

inv(A) – возвращает обратную матрицу (для квадратной невырожденной матрицы (определитель такой матрицы не равен нулю); матрица называется обратной для матрицы A, если при умножении матрицы на матрицу A справа и слева получаем единичную матрицу);

pinv(A) – возвращает псевдообратную матрицу (для заданной матрицы A псевдообратная матрица P имеет размеры (транспонированной матрицы A) и удовлетворяет условиям A*P*A=A и P*A*P=P (заметим, что псевдообратная матрица строится на основе сингулярного разложения));

chol(A) – возвращает разложение (разложение матрицы – это в общем случае представление матрицы в виде суммы или произведения нескольких матриц определенного вида) Холецкого для симметричной положительно определенной матрицы A, т.е. верхнюю треугольную матрицу T такую, что T’*T=A;

lu(A) – возвращает результат LU-разложения матрицы A (разложение на нижнюю треугольную и верхнюю треугольную матрицы; если к этой же функции обратиться в виде [L,U]=lu(A), то она вернет верхнюю треугольную матрицу (U), нижнюю треугольную матрицу (L), причем будет справедливо равенство A=L*U; обращение вида [L,U,P]=lu(A) позволит также получить матрицу перестановки строк, такую что P*A=L*U);

qr(A) – возвращает QR-разложение матрицы A (разложение на унитарную и верхнетреугольную матрицы);

d=eig(A) – возвращает вектор d, координаты (элементы) которого есть собственные значения (собственные числа) матрицы A;

[V,D]=eig(A) – возвращает матрицу собственных векторов V (собственные векторы расположены по столбцам) и диагональную матрицу D собственных значений (матрица Жордана (элементы, расположенные на ее главной диагонали есть собственные значения матрицы A); каноническая форма матрицы A), т.е., по сути, определяется разложение Жордана (но при условия отсутствия в матрице Жордана жордановых клеток порядка, большего единицы) и справедливо равенство A*V=V*D;

s=svd(A) – возвращает вектор, элементы которого являются искомыми сингулярными числами (сингулярными значениями);

[U,S,V]=svd(A) – возвращает диагональную матрицу S (элементами матрицы S являются и унитарные матрицы U (элементами матрицы U являются ортонормированные собственные векторы, соответствующие наибольшим собственным значениям матрицы A*A') и V (элементами матрицы V являются ортонормированные собственные векторы, соответствующие наибольшим собственным значениям матрицы A*A'), такие, что имеет место равенство A=U*S*V;

[U,T]=shur(A) – возвращает ортогональную (унитарную) матрицу U и верхнюю треугольную матрицу T (по диагонали матрицы T расположены собственные значения матрицы A), такие, что имеет место равенство A=U*T*U';

[V,J]=jordan(A) – возвращает жорданово разложение матрицы A, где J – жорданова форма матрицы A; V – матрица, столбцы которой образуют жорданов базис.

В системе MATLAB также имеются функции cdf2rdf, hess, qz и rsf2csf осуществляющие приведение матриц к некоторым специальным формам. Для вычисления функций от матриц используются функции с суффиксом m, например, expm, logm и sqrtm, которые в отличие от функций exp, log, sqrt, вычисляющих значения функций поэлементно, выполняют вычисления с квадратными матрицами по правилам линейной алгебры, используя разложения математических функций в степенные ряды.

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

>> A=[2 -3 1; 1 2 -6; 5 1 1];

>> cond(A)

ans =

2.4374

>> det(A)

ans =

>> rank(A)

ans =

>> [V,D]=eig(A)

V =

-0.6085 + 0.0000i 0.1224 + 0.3873i 0.1224 - 0.3873i

0.6172 + 0.0000i 0.7355 + 0.0000i 0.7355 + 0.0000i

-0.4988 + 0.0000i 0.3184 - 0.4390i 0.3184 + 0.4390i

D =

5.8626 + 0.0000i 0.0000 + 0.0000i 0.0000 + 0.0000i

0.0000 + 0.0000i -0.4313 + 4.1075i 0.0000 + 0.0000i

0.0000 + 0.0000i 0.0000 + 0.0000i -0.4313 - 4.1075i

Пусть решается однородная СЛАУ (Ax=0). Тогда, как известно, такая СЛАУ может иметь нетривиальные (ненулевые) решения тогда и только тогда, когда определитель матрицы равен нулю. Заметим, что, если решение такой системы существует, то это некоторое пространство решений. В линейном пространстве множество решений однородной СЛАУ образует подпространство размерности n-r. Базис этого подпространства можно найти посредством применения функции null, которая допускает несколько вариантов обращения:

Z=null(A) – выполняется вычисление матрицы Z, столбцы которой представляют собой ортонормированный базис пространства решений;

Z=null(A,’r’) – выполняется вычисление матрицы Z, столбцы которой представляют собой некоторый рациональный базис пространства решений;

В обоих рассмотренных случаях справедливо равенство Az=0.

Разреженные матрицы

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

>> A=[0 1 0 0; 0 0 1 0; 0 0 0 1; -4 0 0 0]

A =

0 1 0 0

0 0 1 0

0 0 0 1

-4 0 0 0

>> S=sparse(A)

S =

(4,1) -4

(1,2) 1

(2,3) 1

(3,4) 1

>> whos S

Name Size Bytes Class Attributes

S 4x4 104 double sparse

>> whos A

Name Size Bytes Class Attributes

A 4x4 128 double

Приведенный пример демонстрирует, что функция sparse размещает в памяти компьютера только лишь индексы ненулевых элементов и их значения, при этом матрица S является переменной типа sparse array, занимает в памяти меньше места, чем исходная матрица A.

Обратное преобразование разреженной матрицы в обычную производится функцией full:

>> full(S)

ans =

0 1 0 0

0 0 1 0

0 0 0 1

-4 0 0 0

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