Самоучитель Matlab
В принципе, похожие результаты можно получить и с помощью функции jordan(). Правда, процесс более сложный. Например, следующей командой вычисляем матрицу преобразования и соответствующую ей диагональную форму:
>> [W,D]=jordan(A) |
|
|
|
|
W = |
0.7346 |
+ 0.0000i |
-7.8362 |
+ 0.0000i |
-0.2316 - 0.0000i |
||||
-0.7009 + 0.0000i |
1.1841 |
- 0.0000i |
4.0168 |
- 0.0000i |
1.0000 |
1.0000 |
|
1.0000 |
|
D = |
0 |
|
0 |
|
-5.9000 + 0.0000i |
+ 0.0000i |
|||
0 |
8.3080 |
0 |
||
0 |
0 |
|
-0.4080 |
- 0.0000i |
На заметку
Формально элементы полученных матриц комплексные. Однако мнимая часть у них нулевая. Появление нулевых мнимых частей связано со спецификой используемого алгоритма вычисления жордановой формы и матрицы перехода.
Однако транспонированная матрица преобразования не равна обратной, в чем несложно убедиться, воспользовавшись командой
>> W*W' |
|
|
|
|
|
|
ans = |
|
|
-30.4444 |
+ 0.0000i |
-7.3333 |
- 0.0000i |
62.0000 |
|
0.0000i |
||||
-30.4444 - |
18.0278 |
- 0.0000i |
4.5000 |
+ 0.0000i |
||
-7.3333 |
+ |
0.0000i |
4.5000 |
3.0000 |
|
|
Результаты этих вычислений представлены на рис. 4.35.
Причина этой неприятности кроется в том, что хотя столбцы полученной матрицы преобразования состоят из собственных векторов исходной, преобразуемой, матрицы, они не нормированы. В этом все дело! Убедимся.
Переопределяем матрицу преобразования путем нормирования столбца матрицы на модуль этого столбца. Чтобы нормировать столбец с индексом k (возможные значения 1, 2 или 3), используем команду W(:,k)/ norm(W(:,k)). Здесь мы воспользовались тем, что ссылка на столбец матрицы с данным индексом выполняется в виде W(:,k). Вся команда переопределения матрицы W выглядит так:
>> W=[W(:,1)/norm(W(:,1)),W(:,2)/norm(W(:,2)),W(:,3)/ norm(W(:,3))]
W = |
- |
0.0000i |
0.4283 |
+ |
0.0000i |
-0.8842 |
+ |
0.0000i |
-0.1864 |
||||||||
-0.5639 |
+ |
0.0000i |
0.6904 |
- |
0.0000i |
0.4532 |
- |
0.0000i |
0.8046 |
|
|
0.5831 |
|
|
0.1128 |
|
|
196
Глава 4. Элементы матричной алгебры
Рис. 4.35. Диагонализация симметричной матрицы с помощью функции jordan()
Теперь все изменилось, и произведение матрицы W на транспонированную дает единичную матрицу:
>> W*W' |
|
|
|
|
|
|
ans = |
|
|
-0.0000 |
+ 0.0000i |
-0.0000 |
- 0.0000i |
1.0000 |
|
0.0000i |
||||
-0.0000 - |
1.0000 |
- 0.0000i |
0.0000 |
+ 0.0000i |
||
-0.0000 |
+ |
0.0000i |
0.0000 |
1.0000 |
|
|
При этом матрица W является матрицей перехода:
>> W*D*W' |
|
|
|
|
|
|
|
|
ans = |
- 0.0000i |
2.0000 |
- 0.0000i |
3.0000 |
+ 0.0000i |
|||
1.0000 |
||||||||
2.0000 |
+ |
0.0000i |
2.0000 |
+ |
0.0000i |
6.0000 |
- |
0.0000i |
3.0000 |
- |
0.0000i |
6.0000 |
+ |
0.0000i |
-1.0000 |
+ |
0.0000i |
Ее определитель равен единице:
>> det(W) ans =
1.0000 - 0.0000i
Результаты этих вычислений представлены на рис. 4.36.
197
Самоучитель Matlab
Рис. 4.36. Создание на основе матрицы преобразования матрицы перехода
В заключение отметим, что если необходимо только вычислить конечную диагональную форму, без вычисления матрицы перехода, можно воспользоваться функцией diag(), которой создается диагональная матрица. Аргументом функции передается вектор (список) диагональных значений. Так, если нужно создать диагональную форму для матрицы A, можно воспользоваться командой diag(eig(A)).
198
Глава 5 Решение уравнений и оптимизация
- Для себя я все решил. - Здесь решаю я!
К/ф "Клуб самоубийц, или Приключения титулованной особы"
В этой главе обсуждаются методы решения средствами Matlab алгебраических уравнений и систем уравнений, а также задачи оптимизации и математического (линейного и нелинейного) программирования. Сначала рассмотрим методы решения нелинейных алгебраических уравнений. Напомним, что системы линейных уравнений рассматривались в предыдущей главе в качестве примера операций с матрицами и векторами, поэтому здесь на них останавливаться не будем.
Наше изучение методов решения уравнений начнем с краткого миниобзора встроенных функций Matlab, предназначенных для этого. Затем остановимся на некоторых методах реализации поиска корней уравнений и систем уравнений средствами программирования.
Встроенные функции для решения уравнений
Курс у нас один – правильный.
В. Черномырдин
Для решения одного, отдельного уравнения можно воспользоваться встроенной функцией fzero(). У функции несколько способов вызова. Мы рассмотрим те, что наиболее востребованы с практической точки зрения. Допустим, нам необходимо решить уравнение f (x) = 0 . Функция f (x) полностью определяет решаемое уравнение. Ссылку на эту функцию необходимо передать первым аргументом функции fzero(). Ссылкой может служить указатель на функцию, имя (без расширения и в одинарных кавычках) m-файла с кодом функции или заключенное в одинарные кавычки выражение, определяющее функциональную зависимость f (x). В последнем случае независимая переменная должна называться x , и никак иначе.
Обычно для поиска решения в числовом виде (а именно об этом идет речь) указывают интервал, на котором ищется решение или, по крайней мере, начальное приближение для корня. Начальное приближение или интервал поиска (в виде списка из двух элементов – границ интервала) передается вторым аргументом функции fzero().
199
Самоучитель Matlab
На заметку
Если вторым аргументом передается интервал поиска корня, то на границах интервала функция уравнения должна принимать значения разных знаков, иначе функцией fzero() результат вычислен не будет.
Возвращает функция fzero() много полезной информации. Но самое главное – она возвращает корень уравнения. Рассмотрим пример, представленный в документе на рис. 5.1.
Рис. 5.1. Решение уравнения с помощью функции fzero() на интервале
Командой x=fzero('2*sin(x)^2+sin(x)-1',[0 pi/2]) вычисляется корень уравнения 2 sin(x)2 + sin(x) −1 = 0 на интервале значений переменной x от 0 до π 2 (на этом интервале точное решение x = π 6 ). Если воспользоваться командой x=fzero('2*sin(x)^2+sin(x)-1', [pi/2 pi]), получим решение этого же уравнения на интервале от π 2 до π (решение x = 5π
6 ).
На рис. 5.2 приведен фрагмент документа, в котором вместо интервала поиска решения указывается начальная точка для поиска решения.
Рис. 5.2. При вычислении корня указывается начальное значение для поиска решения
200