Раздел: Документация
0 ... 15 16 17 18 19 20 21 ... 78 Определение. Числом обусловленности простого собственного значения X называется секанс утла между строчными и столбцовыми собственными векторами, соответствующими X. Подробнее см. [19, гл. 2]. Это число есть также норма спектрального проектора на собственное пространство, соответствующее X, и является коэффициентом при старшем члене при разложении в ряд возмущения собственного значения X. Число обусловленности собственного значения X, связанного с жордановой клеткой, равно -°°. Все собственные значения всех симметричных матриц полностью обусловлены (т.е. число обусловленности равно 1). Изменение собственного значения ограничено спектральной нормой изменения матрицы. Пакет программ EISPACK обеспечивает почти предельно возможную точность при вычислениях малых матриц. Действительно, неопределенность в плохо обусловленных задачах может быть очень велика. Тем не менее точность вычисления хорошо обусловленных собственных значений матрицы не ухудшается из-за наличия плохо обусловленных значений в той же матрице. 6.4. ПОЛОЖЕНИЕ ДЕЛ ПРИ РЕШЕНИИ ПРОБЛЕМ, ПРИВОДЯЩИХ К БОЛЬШИМ МАТРИЦАМ Картина резко меняется, когда матрица (на входе) и результаты (на выходе) не могут быть размещены в запоминающем устройстве с произвольной выборкой. Такие задачи, как уже упоминалось во введении, весьма актуальны для пользователей. К настоящему времени разработано несколько хороших методов, продолжается активный поиск новых методов, но пока в нашем распоряжении нет ничего подобного пакету программ EISPACK. Некоторые хорошие методы реализованы в пакетах программ Университета шт. Индиана, предназначенных для специальных целей (например, NAST-RAN, ADINA для структурного анализа) и анализа процессов обмена в квантовой химии [3J. В общем случае пользователю следует обращаться к каталогу [18J, начатому в 1982 г. Необходимо подчеркнуть два аспекта проблемы. 1.Следует использовать любую структуру разреженности матрицы для достижения успеха в решении задачи. 2.Эффективность метода в значительной степени зависит от самой вычислительной системы (от наличия виртуальной памяти или обмена информацией между запоминающими устройствами первого и второго уровней). Информация, которую требуется получить при исследовании больших матриц, менее разнообразна, чем при исследовании малых матриц. В основном преобладают случаи анализа симметричных матриц. Автору неизвестны случаи необходимости определения всех собственных значений и всех собственных векторов, правда, известна одна прикладная задача, при решении которой используется весь спектр матрицы. В основном требуется только несколько собственных значений (их может быть 3, а может быть и 50) в левой части спектра (вблизи 0) и соответствующих собственных векторов. В 1978 г. одна техническая компания потратила 12000 дол. на машинное 112 время для определения 30 собственных значений и собственных векторов при решении задачи с п = 12000. Эта сумма не включает стоимость разработки программы. . Есть прекрасный способ использования разреженности матрицы, позволяющий не разрабатывать различные методы для каждого практически важного распределения в ней нулевых элементов. Пользователь составляет программу, в которой каждый и-мерный вектор v на входе преобразуется в «-мерный выходной вектор и по правилу u=Av. Затем рассчитываются собственные значения и собственные векторы при введении в программу соответствующих для данной задачи векторов v. При данном методе доступ к линейному оператору А и его преобразование невозможны. В этом случае пользователь должен привлечь к вычислениям все, что известно относительно оператора А, чтобы сделать вычисление u( = Av) наиболее эффективным. На практике оператор А не будет простей матрицей. В типичном случае оператор А содержит диагональную матрицу D, структурированную матрицу Н и параметр сдвига о, так что A = D(H - aD2)-1D. Выходной вектор и определяется в три этапа: 1)формируется вектор w = Dv; 2)решается матричное уравнение (Н - aD2 ) х = w относительно х; 3)формируется вектор u=Dx. На этане 2 для решения матричного уравнения можно использовать какой-либо итерационный метод, если другие методы непригодны. Выполнение этапа 2 для больших матриц иногда требует использования дополнительного ЗУ. Это простое, но существенное преобразование уравнения (Н - XD2 ) г -= 0 к стандартному виду (А - ХГ) у = 0 всегда используется в структурном анализе. Не удивительно, что векторные компьютеры (в частности, типа Сгау-1 и Cyber 205) очень удобны для решения задач, в которых векторы длинные и заполненные, а матрицы - разреженные. Надежный путь решения матричного уравнения (Н - ХМ) г - 0 состоит в использовании алгоритма QZ в виде серийного пакета прикладных программ. Однако согласно этой процедуре определяются все собственные значения независимо от того, нужны они или нет, и разрушается симметрия. Кроме того, она требует размещения в быстродействующем ЗУ двух матриц размера пХп. Поэтому эта процедура используется только для малых матриц, как будет описано в следующем разделе. На этом заканчиваем обзор. В разд. 6.5 опишем в общем виде основные методы, а сейчас дадим несколько общих советов "охотникам" за собственными значениями. Если ваши матрицы малы, то (в США) пользуйтесь пакетом программ EISPACK и руководством к нему [5, 16]. Для большинства задач требуемое число арифметических операций лежит в пределах от А3 до 10А3. Если вы имеете дело с большими матрицами, то следует обратиться к каталогу программного обеспечения [18] для разреженных матриц. 6.5. МЕТОДЫ АНАЛИЗА МАЛЫХ МАТРИЦ Все методы, реализованные в пакете программ EISPACK, основаны на двухстороннем преобразовании исходной матрицы А. Такие методы целе- 113 сообразно применять при наличии ЗУ большой емкости с произвольной выборкой, которым снабжен последовательный центральный процессор ЭВМ. Более того, как для симметричных (что не удивительно), так и для несимметричных (что удивительно) матриц преобразования представляют собой исключительно ортогональные операции конгруенции вида где Р - ортогональное преобразование (т.е. PPf = PrP=I). Остановимся на этом последнем утверждении более подробно. Из линейной алгебры известно, что для упрощения задачи можно осуществить преобразование подобия В-э-FBF-1, так как собственные значения матрицы В при таких преобразованиях сохраняются. Однако в том случае, когда элементы матрицы В являются неопределенными и присутствуют погрешности округления, необходимо учитывать характер преобразования F. Оказывается, что любая неопределенность в элементах матрицы В увеличивается примерно в cond(F) : Ff!-HF"11 раз, где cond(F) - число обусловленности матрицы F. Поскольку 1 <cond(F) < 00, то при вычислениях на ЭВМ следует ограничиваться преобразованиями с небольшими значениями величины cond(F). Самая надежная стратегия, конечно, состоит в том, чтобы обеспечить cond(F) = 1. Достоинство ортогональной матрицы Р состоит в том, что: 1)P~, = Pf (т.е. обращение матрицы является тривиальным); 2)cond (Р) = 1; действительные симметричные матрицы всегда могут быть приведены к диагональному виду соответствующим ортогональным преобразованием подобия, но лучшая форма, к которой можно привести в общем случае матрицы, отличные от нормальных, - это треугольная. И последнее замечание относительно терминологии. Конгруенция - это преобразование вида B-*FBF?; в том случае, когда преобразование F ортогонально, конгруенция является преобразованием подобия. Методы Якоби [12, с. 202] состоят в выполнении последовательности простых преобразований Р (называемых плоским вращением), в результате каждого из которых становится равной нулю (или сильно уменьшается) одна пара внедиагональных элементов матрицы, например пара, расположенная на позиции (?,/) и (/,/). Симметричная матрица постепенно приводится к диагональному виду. Конечно, элементы матрицы, которые становятся нулевыми на одном этапе преобразования, теряются на последующих этапах. Тем не менее процесс сходится достаточно быстро при правильной стратегии выбора последовательности пар элементов (/, i). Для матриц, в которых доминирует главная диагональ, методы Якоби являются идеальными, но для обычных матриц их быстродействие в два-три раза меньше [12, гл. 9], чем обеспечивает наиболее подходящий метод, называемый QR-методом Хаусхолдера. Последнее справедливо для обычных последовательных ЭВМ. Метод Хаусхолдера реализуется в два этапа [12, с. 202]. 1.С помощью последовательности (п - 1) операций ортогональной конгруенции симметричная матрица А приводится к трехдиагональному виду Т, а несимметричная матрица В приводится к матрице Хассенберга 114 Н=  2. С помощью другой последовательности матриц Р уменьшается число ненулевых внедиагональных элементов матриц Т и Н без "повреждения" нулевых элементов. Последний в недиагональный элемент уменьшается очень быстро. Когда он становится пренебрежимо малым, согласно программе просто игнорируются последние строка и столбец матрицы (такая процедура называется автоматическим понижением порядка матриц). К концу данного этапа матрица Т становится диагональной. Даже в случае несимметричной исходной матрицы ее собственные значения лежат на диагонали результирующей верхней треугольной матрицы. Для неспециалиста важным является то, что для получения матрицы Т требуется (2/3) и3 операции умножения с плавающей запятой, а для получения матрицы Н - (5/3)и3 операций умножения [20]. Для выполнения второго этапа требуется приблизительно 10и2 операций умножения в случае симметричной исходной матрицы и примерно (5/2)и3 операций умножения в случае несимметричной матрицы. Не было случая, чтобы при выполнении программы возникали какие бы-то ни было затруднения. Методы устойчивы и обеспечивают сходимость, а программы неоднократно испытывались. Необходимо обратить внимание на один технический момент. Когда действительная матрица имеет комплексные собственные значения, матрица Н будет действительной, все операции вычислений на втором этапе будут действительными, но не все внедиагональные элементы матрицы будут стремиться к нулю. На самом деле матрица приводится к блочно-треугольному виду с диагональными блоками размера 2X2 для каждой комплексно-сопряженной пары собственных значений. Это очень эффективный способ сохранения арифметических операций действительными, но получения при этом комплексных величин. В том случае, когда необходимо получить собственные векторы (при симметричной исходной матрице), произведения всех Р матриц могут быть накоплены по мере их образования для получения таким образом ортогональной системы собственных векторов, как бы близки ни были при этом собственные значения. Стоимость накопления этих произведений превышает стоимость остальных вычислений. Иногда требуется найти всего несколько собственных значений; в этом случае можно использовать другие методы. Если при симметричной исходной матрице необходимо точно определить более одной четвертой части спектра, то для увеличения быстродействия лучше прибегнуть к QR-методу Хаусхолдера для определения сразу всех собственных значений, чем находить нужные собственные значения матрицы по одному [16, с. 80]. При несимметричной исходной матрице ситуация еще сложнее. Если нужно найти только одно собственное значение, то в соответствии с принципом действия 1 IS программ EISPACK необходимо определить их все. Эта ситуация встречается в известной задаче о показателе устойчивости: лучший из известных способов определения одного из собственных значений с наибольшей действительной частью малой матрицы состоит в вычислении всех собственных значений; а затем в их исследовании. При других методах тратится много усилий на то, чтобы убедиться, является ли полученное собственное значение искомым, из-за чего они не могут соперничать с упомянутым выше. (Показатель устойчивости - это действительная часть этих собственных значений). Все это весьма удивительно! 6.6. ОБРАТНАЯ ИТЕРАЦИЯ Существует одна техническая задача, которую часто неправильно интерпретируют. Представьте, что пользователь ищет все собственные значения в заданном интервале симметричной трехдиагональной матрицы размера 100X100. Предположим, что их всего три. Эти значения можно найти с помощью простого метода, называемого методом деления пополам. Проблема состоит в вычислении соответствующих собственных векторов. Наиболее предпочтительный метод заключается в использовании одного шага обратной итерации, а именно решении матричного уравнения (Т - X) s = b относительно s, где X - найденное собственное значение, a b - вектор, выбранный так, чтобы он зависел от X и не был ортогонален искомому собственному вектору. Полагаем, что значение X вычислено с достижимой по практике точностью, так что матрица Т -X является практически сингулярной. Тем не менее единственным компьютерным решением матричного уравнения (Т - X) х = 0 будет нулевое, так что необходимо использовать некоторое очень малое, но ненулевое значение правой части Ь. Программа TINVIT из пакета EISPACK, которая реализует этот метод, дает хорошую аппроксимацию собственного вектора, соответствующего X, Это можно проверить апостериори, так как величина ll(T-X)s b INIlis Н известна и она должна иметь один порядок малости с б • Т, где е - единица округляемого разряда. В чем трудности? Предположим, что два собственных значения Xj и Х2 почти совпадают. В этих случаях вектор b изменяется очень мало при замене Xi на Х2 и, следовательно, вычисленный собственный вектор s2, хотя и не очень близок к вектору sj, к сожалению, далеко не ортогонален вектору Si. Конечно, как показывает результат их вычитания, векторы Sj и s2 прекрасно аппроксимируют собственные векторы матрицы Т. Плохо то, что векторы Sj и s2 лежат в плоскости, которая является собственным пространством с достижимой на практике точностью. Поэтому в действительности любая линейная комбинация векторов sx и s2 будет такой же хорошей аппроксимацией, как и вектор s2 116 Отсюда возникает миф о существенной трудности вычисления ортогональных собственных векторов для очень близких собственных значений матрицы, в то время как на самом деле трудности связаны только с методом обратной итерации. Недостаток заключается в том, что вектор b является очень гладкой функцией X. Раз причина выяснена, ее нетрудно устранить. Действительно, в последней версии программы TINVIT учтены труд-юности получения почти ортогональных последовательностей собственных векторов. Однако, чтобы включить процедуры оценок и ортогонализации, в жертву пришлось принести элегантность и простоту первоначального метода. Как упоминалось в разд. 6.5, накопление произведений всех QR-преоб-разований дает полное множество практически ортогональных собственных векторов. 6.7. О МЕТОДАХ РЕШЕНИЯ ПРОБЛЕМЫ СОБСТВЕННЫХ ЗНАЧЕНИЙ БОЛЬШИХ СИММЕТРИЧНЫХ МАТРИЦ Преобразования подобия разрушают любую структуру разреженности матрицы А и, как известно, для больших матриц в явном виде не используются. Кажется, все специалисты по численному анализу отдают предпочтение алгоритму Ланцоша. Однако большинство инженеров пользуется своими любимыми вариантами метода одновременной итерации, который они называют итерированием подпространства. Здесь будут более или менее подробно рассмотрены оба зти метода. Вычислители-химики используют либо метод итерирования подпространства, либо специальные методы, которые не являются общими, но применимы для решения их специфических задач [3]. Пока нет ничего сравнимого с пакетом программ EISPACK и не ясно, появится ли вообще, если учесть сильное влияние на эффективность вычислений конкретной компьютерной системы. Нужно подчеркнуть, что как метод Ланцоша (LAN), так и метод итерирования подпространства (SI) относятся к обобщенным методам и допускают существенные вариации при их использовании. В этом они противоположны QR-алгоритму, который оставляет мало возможностей для вариаций. Проблема выбора сдвига при реализации QR-алгоритма решена, в то время как в методе SI в нем есть что-то от искусства. Многие (но не все) большие матрицы являются разреженными (т.е. имеют менее 20 ненулевых элементов в одной строке). Такие матрицы часто можно хранить в быстродействующем (сверхоперативном) ЗУ вычислительной системы при условии, что используются компактные структуры данных для записи ненулевых элементов. Тем не менее важные вычисления, в которых каждое произведение Аи требует обмена данными между первичным (сверхоперативным) и вторичным (оперативным) запоминающими устройствами, находятся на пределе возможностей современных компьютерных систем. Специалисты в области квантовой химии имеют дело с простыми симметричными матрицами с иЮ4 и 107, которые хранятся в файлах или на лентах. Химикам помогает то, что грубая аппроксимация нескольких 117 0 ... 15 16 17 18 19 20 21 ... 78
|
