Содержание к диссертации
Введение
Глава 1. Уравнения первого порядка 21
1.1. Точные формулы для решения линейных дифференциальных уравнений первого порядка 21
1.2. Оценки приближения 29
1.3. Оценки спектра матрицы Л 37
1.4. Способы приближения 40
1.4.1. Приближение многочленами 40
1.4.2. Приближение рациональными функциями 45
1.5. Случай необратимой матрицы D 51
Глава 2. Уравнения второго порядка 54
2.1. Факторизация уравнения второго порядка с помощью матричного корня его пучка 54
2.2. Факторизация уравнения второго порядка с помощью матричного корня обратного пучка 60
2.3. Связь между матрицами Z и V 62
2.4. Сведение дифференциального уравнения второго порядка к двум уравнениям первого порядка . 64
Глава 3. Аналитические функции от квадратичного пучка 69
3.1. Точные формулы для решения линейных дифференциальных уравнений второго порядка 71
3.2. Функции от факторизованного матричного пучка 77
3.2.1. Представление функций от пучка с помощью жордановой формы 77
3.2.2. Произведение функций от квадратичного пучка . 85
3.2.3. Связь функций от пучка с уравнением Сильвестра . 86
3.3. Функции от нефакторизованного матричного пучка . 90
Литература 95
Приложение 105
- Способы приближения
- Факторизация уравнения второго порядка с помощью матричного корня обратного пучка
- Сведение дифференциального уравнения второго порядка к двум уравнениям первого порядка
- Функции от факторизованного матричного пучка
Введение к работе
Настоящая диссертация посвящена построению и анализу приближенных операторных методов решения систем линейных дифференциальных уравнений первого и второго порядка вида Dx + Вх = /(*) (1) Nx + Dx + Вх = f(t) (2) на основе аналитических функций от матриц D, В и N. Характерная особенность рассматриваемых уравнений заключается в том, что они являются неразрешенными относительно старшей производной и все матричные коэффициенты D, В и N являются самосопряженными и неотрицательно определенными. В диссертации основное внимание уделяется случаю, когда матрицы D, В и N имеют большой порядок.
Такие уравнения возникают, например, в теории линейных электрических цепей [2, 17, 33, 41, 42, 79]. Одним из наиболее сложных примеров таких цепей являются дискретные модели длинных связанных линий, т.е. нескольких параллельных проводов, расположенных на близком расстоянии друг от друга. Расчет больших линейных электрический цепей, в частности, связанных линий, актуален в связи с задачами проектирования микросхем, см. например, [12, 22, 79].
Современное развитие электроники движется в сторону создания все более миниатюрных и высокоскоростных микросхем. Такая схема представляет собой сложную систему, состоящую из миллионов полупроводниковых устройств и соединяющих их проводов. С целью уменьшения размера микросхемы, провода стараются располагать как можно ближе друг к другу, а с целью увеличения быстродействия — проходящие по ним импульсы стараются сделать как можно более короткими и тем самым с большим содержанием высоких частот. В результате взаимное влияние сигналов, проходящих по соседним проводам, оказывается существенным и его приходится учитывать.
С целью учета влияния соединительных проводов их рассматривают как часть принципиальной электрической схемы. При этом длинные линии обычно заменяют дискретными моделями, содержащими достаточно большое количество секций. Использование ЯС-моделей приводит к линейному дифференциальному уравнению первого порядка (1), а і?С-моделей — к линейному дифференциальному уравнению второго порядка (2).
Для анализа работы микросхемы приходится решать систему, состоящую из огромного числа уравнений, как линейных, так и нелинейных. С целью упрощения этой системы в инженерных расчетах пользуются следующей ее особенностью. Оказывается, нелинейные уравнения можно отделить от линейных. Поэтому линейные уравнения обычно решают независимо, а затем результат подставляют в оставшиеся нелинейные. Однако размерность системы уравнений, описывающих линейные части схемы, все равно остается большой.
Свойства матричных коэффициентов возникающих линейных уравнений зависят от способа составления уравнений. Существует несколько вариантов составления уравнений, что соответствует различным модификациям уравнений Кирхгофа [2, 17, 33, 41, 42, 79]. Для нас важно, что уравнения можно выписать таким образом, чтобы матричные коэффициенты D, В и N в (1) и (2) были самосопряженными и неотрицательно определенными.
Существует обширная литература, в которой линейные дифференциальные уравнения (в основном, вида (1)), неразрешенные относительно старшей производной, изучаются со спектральной точки зрения, см. например, работы А.Г. Баскакова, СТ. Крейнаи К.И. Чернышова [3, 4, 26], А.Г. Руткаса [34, 37, 38], Г.А. Свиридюка [40] и литературу в них. Эти работы имеют некоторые идейные точки соприкосновения с настоящей диссертацией. Основное отличие настоящей работы от перечисленных заключается в том, что в ней дифференциальные уравнения вида (1) рассматриваются с точки зрения возможности эффективного построения приближенного решения.
Абстрактная спектральная теория дифференциальных уравнений второго порядка вида (2) и более высокого порядка разрабатывалась в работах М.В. Келдыша [23], Ф.Р. Гантмахера [10], М.Г. Крейна и Г.К. Лангера [25], А.С. Маркуса [30], И.М. Глазмана и Ю.И. Люби-ча [11] и многих других. Важную роль в этой теории играют так называемые операторные пучки (в случае уравнения (2) это функция А і-)- X2N + XD + В). Некоторые из этих результатов, а также общая идеология используются в настоящей диссертации в главах 2 и 3. В диссертации уравнения вида (2) также рассматриваются с несколько иной точки зрения, чем в перечисленных работах, а именно, обсуждается возможность получения с помощью абстрактной теории эффективных методов приближенного решения.
Напомним еще раз, что диссертация посвящена приближенным операторным методам решения уравнений (1) и (2). Классический подход в теории дифференциальных уравнений [32] рекомендует для решения линейного дифференциального уравнения первого порядка вида (1) в случае обратимой матрицы D заменить исходное уравнение (1) на уравнение ~x{t) = Ax(t)+ 0^/(1), (3) где А = — D~XB, и с помощью жордановой формы матрицы А найти фундаментальное решение eAt. К сожалению, на практике процесс вычисления жордановой формы оказывается неустойчивым [14, 21]: даже незначительные возмущения могут сильно испортить результат вычисления. Известны различные модификации этого подхода [14, 21, 66], основанные на других спектральных разложениях А — ТАТ~1У например, на представлении матрицы А в форме Шура. В результате получаются более устойчивые процедуры. Однако в случае матриц высокого порядка вычисление функции от матрицы типа eAt спектральными методами является довольно трудоемким. Поэтому в общем случае уравнения (1) и (3) приходится решать приближенно.
Следует отметить, что даже в случае обратимой матрицы D при замене уравнения (1) на уравнение (3) и последующем его решении традиционными методами теории дифференциальных уравнений или численными методами упускается серьезная деталь, а именно, теряется свойство самосопряженности коэффициентов: для самосопряженных матриц D и В произведение A = —D-1B уже самосопряженным не является. В то же время естественно ожидать, что такая важная особенность исходного уравнения, как самосопряженность и неотрицательная определенность матричных коэффициентов, может быть использована при построении решения. В диссертации предлагается приближенный метод, не требующий вычисления собственных значений и собственных векторов и позволяющий учесть указанные свойства матриц D и В. Метод основан на приближении аналитической функции от матрицы А (матричной экспоненты eAt или некоторой специальной функции dt[A)) многочленом или рациональной функцией rt(A).
В случае необратимой матрицы D ситуация усложняется. Классический метод решения уравнений первого порядка с необратимым старшим коэффициентом основан на возможности следующего преобразования, восходящего к Вейерштрассу [80] и Кронекеру [70] (см. также [49, с. 497] и [10, с. 315]): для любой пары матриц D и В существуют невырожденные квадратные матрицы Ти*?, такие что где матрица D\ нильпотентна. Прозрачное объяснение существования такого преобразования, использующее контурные интегралы от резольвенты, можно найти в упомянутых выше работах [3, 4, 26, 40]. Современный численный алгоритм построения подобного преобразования берет начало от работы Стюарта [78], см. также его изложение и обсуждение в [49, с. 511], [14, 7.7] и [18, 4.5]. К сожалению, этот алгоритм в случае матриц большого порядка является довольно трудоемким и поэтому не всегда удобен. Отметим, что предлагаемый в диссертации метод применим и в случае необратимой матрицы D, см. 1.5.
В литературе по методам вычислений [6, 13, 48, 49] для уравнений вида (1) и (3) обычно рекомендуют использовать различные разностные схемы, наиболее популярными из которых являются методы Рунге-Кут-ты. Схемы Рунге-Кутты образуют целое семейство [48, 49], включающее в себя огромное количество численных методов приближенного решения дифференциальных уравнений. Это семейство распадается на два больших класса — явные методы и неявные. Можно показать [49], что применение явных методов Рунге-Кутты для решения однородного уравнения jtx{t) = Ax{t) (и тем самым для нахождения импульсной характеристики уравнения (3)) на первом шаге эквивалентно вычислению некоторого многочлена от матрицы Л, а затем на последующих шагах возведению этого многочлена в степень (в качестве такого многочлена чаще всего используются многочлены Тейлора для экспоненты). Аналогичным образом неявные методы Рунге-Кутты сводятся к вычислению рациональной функции от матрицы А и возведению ее в степень (для большинства неявных методов такими рациональными функциями являются аппроксимации Паде [5] для экспоненты). Поэтому использование методов Рунге-Кутты для матриц высокого порядка связано со значительными вычислительными затратами: чтобы найти, например, решение в точке, отстоящей от начальной на 100 шагов, по сути дела, приходится матрицу высокого порядка возводить в сотую степень.
В диссертации, приближенное решение также строится с помощью некоторой функции от матрицы Л, но выбор такой функции оптимизируется. Поэтому подход, предлагаемый в диссертации, позволяет избежать вычисления высоких степеней и ограничиться вычислением рациональных функций сравнительно небольшого порядка (со степенью знаменателя не выше 15). Кроме того, если необходимо найти решение на отрезке, не содержащем начальную точку, метод, предлагаемый в диссертации, не требует вычисления решения в точках, лежащих между начальной точкой и левой границей интересующего нас отрезка, тогда как в методах Рунге-Кутты расчет всегда необходимо выполнять последовательно от начальной до конечной точки. Это преимущество становится особенно ощутимым, если интересующий нас отрезок расположен достаточно далеко от начальной точки.
В приложениях, связанных с расчетом микросхем, для решения уравнений типа (1) традиционно используются проекционные методы Крылова Арнольди-Ланцоша [54]-[56], [58]-[65], [67]-[69], [71]-[77], [81], которые позволяют понизить порядок системы. На практике эти методы обычно неплохо работают, но присутствие некоторого элемента случайности в выборе приближения иногда приводит к сбоям (например, приближенная модель может оказаться неустойчивой [71]). Кроме того, недостатком этих методов является невозможность оценить точность результата. При таком способе расчета неясно, как именно необходимо построить приближение, какой размерности должна быть новая система, чтобы получить решение, удовлетворяющее заранее выбранной точности. Метод приближенного решения уравнения (1), излагаемый в настоящей диссертации, по объему вычислений сопоставим с методами Крылова-Арнольди-Ланцоша, но в принципе позволяет в явном виде получить оценку точности вычислений.
Некоторые уравнения вида (1) возникают как результат применения метода прямых [6, 7, 27] к некоторому параболическому уравнению, тем самым методы решения уравнения (1) можно интерпретировать как схемы приближенного решения уравнений в частных производных. В терминах нашего основного приложения это означает, что для расчета длинных линий можно использовать не только дискретные, но и непрерывные модели. Во втором случае длинные линии описываются системой телеграфных уравнений ТЭ1 пт dU л 9 с самосопряженными неотрицательно определенными матричными коэффициентами, которая при L — 0 имеет параболический тип. Среди численных методов для уравнений в частных производных к обсуждаемому в работе подходу наиболее близки разностные схемы Кранка-Николсон [7, 57] для параболических уравнений, см., например, работы В.В. Смагина [43, 44, 45]. Отметим, что разностные схемы Кранка-Николсон аналогичны неявным методам Рунге-Кутты, поскольку эти методы также основаны на вычислении рациональной функции от некоторого оператора, возникающего после дискретизации по переменным s и t, и последующем возведении ее в степень. Поэтому сделанные выше замечания о методах Рунге-Кутты относятся и к методам Кранка-Николсон.
Для решения уравнений второго порядка вида (2) в классической теории дифференциальных уравнений обычно рекомендуют [32] заменить исходное уравнение (2) на уравнение первого порядка удвоенной размерности. В приложениях с помощью специальных способов записи уравнений Кирхгофа (метод переменных состояния [79], модифицированный метод узловых напряжений [17, 33]) вместо уравнения второго порядка выписывают уравнение первого порядка с матрицами того же порядка. И в том, и в другом случае матричные коэффициенты уравнения первого порядка, как правило, уже не являются неотрицательно определенными. После этого к полученным уравнениям первого порядка применяют те или иные из упомянутых выше методов решения. Отметим также, что известны [6] разностные схемы типа схем Рунге—Кутты, применимые непосредственно к уравнениям второго порядка и не требующие предварительного сведения этих уравнений к системам первого порядка, однако в существенном (возведение в большую степень рациональной функции от некоторой матрицы) эти методы очень близки к методам Рунге-Кутты для уравнений первого порядка. В диссертации предлагается альтернативный подход для приближенного решения уравнения (2), учитывающий свойства самосопряженности и неотрицательной определенности матричных коэффициентов N, D и В, основанный на факторизации квадратичных операторных пучков [10, 11, 25, 30].
Первая глава посвящена построению приближенного решения линейного дифференциального уравнения первого порядка вида (1) с самосопряженными неотрицательно определенными матричными коэффициентами D и В.
Содержание 1.1 носит вспомогательный характер. В нем определяется изучаемая начальная задача на оси: до некоторого момента времени система находилась в нулевом положении равновесия x(t) = 0, t < to, затем поступил входной импульс в виде свободного члена /, который вызвал некоторый ответный отклик x(t), t > to- С математической точки зрения такая постановка задачи сводится к нахождению решения ж, удовлетворяющего дифференциальному уравнению на всей оси в предположении, что f(t) — 0 левее некоторой точки to> в классе функций х, обладающих тем лее свойством. Кроме того, для приложений обычно интерес представляет решение не самого уравнения (1), а связанной с ним линейной динамической системы Dx(t) + Bx(t) ==bu(t), y(t) = {x(t),l), (4) в которой правая часть имеет вид f(t) — bu{t) (для электрической цепи это соответствует тому, что схема содержит лишь один входной источник напряжения или тока), а вместо решения х достаточно искать лишь одну его координату, имеющую смысл выхода и порождающую функцию у. Здесь b, I — фиксированные векторы, а скалярная функция и задает форму входного сигнала.
Приводятся определения и свойства импульсных характеристик системы (4). Матричной гилпульсной характеристикой называют матричное решение ti-b-H(t) уравнения H{t) = AH{t) + D-l5(t), teR, где S — дельта-функция Дирака, равное нулю на (—оо,0). (Скалярной) импульсной характеристикой системы (4) называют решение h системы (4), соответствующее входу и — 6 и равное нулю на (—оо,0). Для импульсных характеристик справедливы представления eMD~x при t > О, Я (*) ={
О при t < О,
А(*) = Лй(«) = (Я(*)М).
В предложении 3 напоминается, что решение, соответствующее произвольной функции и, может быть получено как свертка импульсной характеристики h с входной функцией и. В предложении 4 показано, что точное решение системы (4) представимо в виде y(t) = {#t(A)D-4j), где fit — специальная аналитическая функция, явно выражающаяся через функцию и.
В 1.2 излагается основная идея метода приближенного решения уравнения первого порядка (1), состоящая в приближении множителя ем в импульсной характеристике аналитической функцией Tt{A) от матрицы А. Обсуждаются способы оценки точности такого приближения. Основными результатами являются теоремы 8 и 12, в которых точность приближения импульсной характеристики h{i) оценивается в терминах максимума разности г*(А) — ех\ где А пробегает спектр <т{А) матрицы А. При этом в теореме 8 оценивается абсолютная ошибка, а в теореме 12 — относительная ошибка приближения.
Теорема 8. Пусть на спектре матрицы А — —D~lB приближение A i-> rt(X) отличается от функции А і-)- ем на є(і) > 0, т.е. \rt(X) - eAf| < є(і) для всех А сг{А).
Тогда для приближенной импульсной характеристики hb(t) = (rt(A)D-4,b) справедлива оценка \kb{t) - hlb{t)\ < s(t)^/hu(+Q)hbb(+Q).
Теорема 12. Пусть на спектре матрицы А = —D В приближение А і-» Г((А) действительно для действительных А и отличается от функции А н-> еЛі на є(і)ел*, т.е. |г((А) - eAl| < e{t)eM для всех А Є о-(А).
Тогда для приближенной "диагональной" импульсной характеристики hu(t) - (r^D-Hj) справедлива оценка \ha(t) - hu(t)\ < e(t)htt{t).
А для приближенной импульсной характеристики hib(t) имеем \kb(t) - ha(t)\ < e(t)(|A»(t)| + >Ди(*)Ы*))-
Справедливы аналогичные утверждения (теоремы 9 и 13) об оценке точности приближенного решения y(t) с помощью разности Г((А) — $((А).
Предполагается, что точно спектр о-(А) не известен. В связи с этим в 1.3 обсуждаются методы его оценки, а именно, способы нахождения отрезка [—/3, —а], содержащего о-(А).
В 1.4 обсуждаются способы выбора приближающей функции Г(. Простейшим способом построения приближения является использование в качестве rt многочлена rt(A) — co(t)l 4- c\(t)A + * + с^^А*1 с коэффициентами Cfc(t), зависящими от t. В пункте 1.4.1 обсуждаются три основных способа построения многочлена rt — использование многочленов Тейлора, интерполяционных многочленов и многочленов наилучшего приближения. Последние дают лучшие результаты, но требуют более сложных вычислений. Интересно также отметить, что выбор многочлена rt зависит от типа оценки точности, которую необходимо получить, иначе говоря, минимизация абсолютной или относительной ошибки требует использования разных многочленов rt. Отметим, что при использовании многочленов хорошего результата можно добиться лишь на некотором отрезке вида [0,to]- Величина to зависит от количества степеней Ак, которые мы готовы вычислить.
В пункте 1.4.2 вместо многочлена г* используются аппроксимации Паде и рациональные функции наилучшего приближения. Оказывается, использование рациональных функций не требует оценки спектра А. Кроме того, такой подход позволяет получать высокую точность далее при небольшом порядке рациональной функции. К сожалению, прямое вычисление рациональной функции г* требует в каждой точке t заново обращать матрицу, стоящую в знаменателе. Частично эту проблему удается решить, если "замораживать" знаменатель и строить приближение на отрезке. Такой отрезок никогда не содержит ноль. Поэтому для построения приближения в окрестности нуля приходится пользоваться многочленами или аппроксимациями Паде без "замораживания". Аналогично случаю многочленов, наилучшее рациональное приближение обеспечивает значительно лучшее качество приближения, но является более трудоемким.
При изложении предлагаемого в главе 1 метода в 1.2—1.4 предполагается, что матрица D обратима. В 1.5 обсуждается случай, когда D вырождена. Здесь показано, что этот случай полностью сводится к рассмотренному.
Во второй главе рассматриваются уравнения вида (2), обладающие специальным свойством: одна из матриц N или В мала, а матрица D обратима. Ограничение на малость матричных коэффициентов обычно выполняется для RLC-цеией с малыми индуктивностями, а обратимость матрицы D соответствует тому, что в цепи достаточно много сопротивлений. В этих предположениях предлагается эффективный метод решения уравнения (2), основанный на его сведении к двум уравнениям первого порядка. Как и в первой главе, для уравнения (2) также рассматривается начальная задача на оси, а правая часть предполагается имеющей вид /() = bu(t), ср. с (4).
Первые три параграфа содержат результаты, необходимые для доказательства основных теорем 38 и 39, приведенных в четвертом параграфе. В 2.1 с помощью теории операторных пучков [11, 30] показано, что дифференциальный оператор L = ЛГ-^ + >4 + В, соответствую- щий уравнению (2), можно представить в виде произведения 1~М + MZ + lW^l — Z\ дифференциальных операторов первого порядка. В 2.2 аналогичным образом для оператора L строится факторизация — ( ^{AV +l) + jllf^V — 1 к В 2.3 устанавливается взаимосвязь между этими представлениями. Наконец, в 2.4 дифференциальное уравнение второго порядка заменяется на два уравнения первого порядка х — Zx = / и Vx ~ х = ft к которым применимы идеи построения приближенного решения из главы 1.
Теорема 38. Решение начальной задачи на оси для уравнения Nx + Dx + Bx=: bu{t) можно представить в виде X — Х\ + Х2, где х\ и Х2 — решения уравнений І1 - Zxi = {NZ + ZN + D)~4u Vi2 -х2 = (BV + D)'lN{NZ + ZN + D)~lbu.
Приводятся формулы (предложения 25 и 29) для эффективного вычисления матриц Z hV.
В приложениях часто возникают уравнения вида (2) с правой частью вида f(t) ~ bu(t). В частности, уравнениями такого вида описываются RLC-цепя. Следующая теорема является аналогом теоремы 38 для таких уравнений.
Теорема 39. Решение начальной задачи на оси для уравнения Nx + Dx + Bx^= bu(t) можно представить в виде х = Xi + Х2, где xi и Х2 — решения уравнений її - Zxi = -{NZ + D)~1B(BV + VB + D)~4u VX2 - X2 = -(BV + 1/B + Z))-1^.
Построение решения для уравнения второго порядка тесно связано с аналитическими функциями от матричного пучка. В связи с этим, в третьей главе изучаются аналитические функции от квадратичных матричных пучков специального вида. А именно, обсуждаются функции от факторизованного матричного пучка F{AB) = ^-J /(А)(А1 - Л)-\Х1 - ВГЧХ (5) и функции F{N, D, В) = і / /(A) {X2N + XD + В)~г dA, (6) тп J г где матричный пучок не факторизован, но матрица N является обратимой. Функции вида (5) и (6) можно рассматривать как обобщение функций от одной матрицы f(A) = ^-^f(x)(xi-ArUx.
В 3.1 показывается, что матричная импульсная характеристика уравнения й-в) (г-*)*<*> = *"« <7> при t > О имеет вид (предложение 42) H(t) = Exp(t, Л, #) = -^7 f ext(Xl - Д)-Х(А1 - В)-1 dA, а решение уравнения (7) представимо в виде (предложение 44) s(i) = МЛ В)Ь = (^т [ 0*(А)(А1 - Л)_1(А1 - #) * <гл)&.
Таким образом, функции H(t) — Ехр(, .4, В) и x(t) = <&t(A,B)b являются примерами аналитических функций от матричного пучка А >—>-(А1-В)(А1-Л).
Основное содержание главы составляет 3.2, в котором для функций вида (5) удается построить достаточно развитую теорию.
В пункте 3.2.1 приводятся явные формулы для нахождения функций от пучка, основанные на жордановой форме матриц Л и В. Пусть Т — матрица перехода к жордановой форме матрицы Л, а Н — матрица перехода к жордановой форме матрицы В. Обозначим через S матрицу Т 1Н. Разобьем ее на прямоугольные блоки S = {S*j} размера uxqjt где Гі — порядок жорданова блока Лі матрицы Л, a qj — порядок жорданова блока Bj матрицы В. Нам понадобится также вспомогательная блочная матрица L = {l
Здесь С+р — биномиальные коэффициенты, а о» и / — собственные значения матриц Лий, отвечающие жордановым блокам Лі u Bj соответственно. Размер блока Ьц равен г» х qj.
Теорема 46. Пусть спектры матриц Л и В не пересекаются. Тогда для F{Ab В) справедливо представление F(A,B) = TQH~\ где Q = Q{f,A,B) = {Qij} — блочная матрица с блоками Qij = f{Ai)Lij - UjfiBj), а блоки Lij определены формулой (8).
Эту теорему можно интерпретировать как аналог представления функций от одной матрицы [53] с помощью ее жордановой формы.
Альтернативная формула для F(A> В) предлагается в теореме 47, где используются функции от матриц «Аи^а также специальная матрица Р, которая зависит только от Л и В и не зависит от /. Эту формулу удобно использовать, если необходимо вычислить несколько функций / от одного и того же пучка.
Теорема 47. Пусть спектры матриц ЛиВ не пересекаются. Тогда функцию F(A,B) можно представить в виде F(A,B) = f(A)P-Pf(B),
Р = TLH~\ (9) а матрица L определена формулой (8).
Приводятся также аналоги этих теорем для диагонализуемых матриц А и В7 а также для случая, когда спектры матриц Аи В пересекаются.
Для функций от одной матрицы хорошо известна теорема о функциональном исчислении [36, 53], утверждающая, в частности, что сумма функций / + g переходит в сумму матриц f(A) + д(А), а произведению функций / - д соответствует произведение матриц f(A) д(А). В случае квадратичного пучка эта теорема перестает быть справедливой. Сумме функций / + д> очевидно, по-прежнему соответствует сумма матриц F(A} В) + С?(Д,), но произведению функций / д уже не соответствует произведение матриц F(A, В) G(A, В). В теореме 52 (пункт 3.2.2) приводится представление для произведения функций от квадратичного матричного пучка.
Теорема 52. Для двух аналитических функций fug имеет место тождество ^1 J ДАМА)(А1 - А)-\Х1 - By1 dX = F(Ay B)g{B) + f{A)G(A> В).
В пункте 3.2.3 предлагается метод нахождения F(A,В), основанный на решении непрерывного уравнения Сильвестра. Он позволяет обойти некоторые вычислительные трудности, связанные с использованием жордановой формы. Дело в том, что для матриц, далеких от нормальных, плохая обусловленность матрицы перехода к жордановой форме — распространенное явление. В этом случае даже незначительные погрешности округления могут сильно испортить результат вычисления по явным формулам. Предлагаются два способа преодоления этой проблемы.
Первый из них основан на представлении Р(А, Б) = f{A)P — Pf{B) из теоремы 47. При этом для нахождения матрицы Р вместо явной формулы (9) рекомендуется решать уравнение ЛР — РВ = 1. Это уравнение представляет собой непрерывное уравнение Сильвестра, для решения которого известны численно устойчивые ортогональные методы [21].
Теорема 53. Матрица Р является решением непрерывного уравнения Сильвестра
АР-РВ = 1.
Как и в пункте 3.2.1, этот способ рекомендуется для вычисления нескольких функций от одного и того же пучка.
Второй способ описывается в теореме 55. Оказывается, функция F является решением уравнения Сильвестра AF(A, В) — F(A, В)В ~ f(A) — f{B)7 которое можно упростить (следствие 56), воспользовавшись разложением Шура [14, 21] матриц Лий.
Теорема 55- Функция F(A,B) удовлетворяет непрерывному уравнению Сильвестра AF(A, В) - F(A, В)В = f(A) - f(B).
В 3.3 обсуждаются матричные функции вида (6) в предположении, что матрица N обратима. Примером функций такого вида являются матричная импульсная характеристика уравнения (2) (предложение 42)
Я () = Ехр(, N,D,B) = ^-( eXt(X2N + XD + Б)"1 dX,
27гг уг а также решение уравнения (2), соответствующее правой части вида bu(t) (предложение 44), x(t) = tit{N, D, B)b = Г-ї-т f tf*(A)(A2iV + XD + B)_1 dx\ 6.
Приводятся явные формулы для вычисления многочленов и рациональных функций от матричного пучка A t-> X2N + AD + В, которые можно использовать для построения приближенной импульсной характеристики и приближенного решения уравнения (2). Отметим, что результаты
3.3 могут быть применены в ситуациях, когда дифференциальный оператор L = N^+D^+B факторизовать не удается, и поэтому ни методы из главы 2, ни идеи 3.2 применить нельзя.
В приложении приводятся результаты численных экспериментов, иллюстрирующие применение методов настоящей диссертации к расчету некоторых линейных электрических цепей.
Основные результаты работы докладывались на: международной научной конференции "Нелинейный анализ и функционально-дифференциальные уравнения" — Воронеж, 2000; Воронежской зимней математической школе "Современные методы теории функций и смежные проблемы" — Воронеж, 2001; межвузовской научно-технической конференции "Новые технологии в научных исследованиях, проектировании, управлении, производстве" — Воронеж:, 2001; Воронежской зимней математической школе — Воронеж, 2002; международной научной конференции "Современные проблемы функционального анализа и дифференциальных уравнений" — Воронеж, 2003; Воронежской зимней математической школе "Современные методы теории функций и смежные проблемы" — Воронеж:, 2005; 12-ой Всероссийской межвузовской научно-технической конференции студентов и аспирантов "Микроэлектроника и информатика - 2005" — Москва, 2005; семинаре ВТ. Курбатова (2000-2005, ЛГТУ), семинаре Б.Н. Садовского (2004, ВГУ).
Основные результаты работы опубликованы в [82]-[88]. В совместной с ВТ. Курбатовым работе [85] автору диссертации принадлежит метод построения приближенного решения дифференциального уравнения второго порядка, а ВТ. Курбатову — утверждения о соотношении между корнями квадратичного операторного пучка. В совместной с ВТ. Курбатовым статье [86] автору диссертации принадлежит метод построения приближенного решения дифференциального уравнения первого порядка и оценки его точности, а ВТ. Курбатову — результаты о связи операторных методов с методами Арнольди-Ланцоша. Только результаты, принадлежащие лично автору диссертации, включены в диссертационную работу.
Способы приближения
В этом и следующем пунктах этого параграфа обсуждается, как можно выбирать функцию г . В этом пункте в качестве г используется многочлен с коэффициентами с (), в общем случае зависящими от t. Рассматриваются три варианта выбора приближающего многочлена г$: многочлен Тейлора, интерполяционный многочлен с узлами интерполяции в нулях многочлена Чебышева первого рода и многочлены наилучшего приближения. С вычислительной точки зрения многочлен является самой удобной функцией: в этом случае для вычисления rt(A) = co(f)l + ci(t)A + + Cn(t)An необходимо лишь найти степени А, А2, ..., Ап матрицы А и подставить их в многочлен. Замечание 7. Если нас интересует приближенное решение (1.12) или приближенная импульсная характеристика (1.14), не обязательно вычислять степени Ак целиком. Достаточно ограничиться итерациями AkD b или 1 Ак матрицы А. Действительно, имеем, например, hib(t) = (г((А)1? 16, /) = EIUctW M) = J =QCk{t)(D lbt(Akyi). Это наблюдение хорошо известно в теории проекционных методов Крылова [14]. Пусть задана функция / действительного или комплексного переменного. Ниже под / обычно подразумевается еА или #t(A). является использование многочлена Тейлора (в нуле): Если, например, /(Л, ) еА , то р«(А,) — Х) =о її А И В качестве приближения для функции Н{) = eAtD l можно брать функцию Н{) = Другим способом построения приближающего многочлена является использование интерполяционных многочленов. Пусть задана п+1 точка AQ, AI, ..., АЛ, лежащая в области определения функции /. Напомним, что интерполяционным многочленом в форме Лагранжа (функции /, построенным по точкам Ао, Ai, ..., Ап) называют (см., например, [51]) многочлен р(А) = ао + ai А + аг А2 Н \- ап Ап степени п, совпадающий в точках Ао, Ai,..., Aft с функцией /: При этом точки Ао, Aj, ..., А„ называют узлами интерполяции. Рассмотрим фундаментальные интерполяционные многочлены Нетрудно видеть, что Jjt(Aj) = 0 при к ф j, и /fc(Aj) = 1 при к = j, A;, j = 0,1,..., п. Очевидно, Вернемся к нашей основной задаче приближения аналитической функции / многочленом rt на отрезке [—/?, — х]. Одним из сравнительно хороших универсальных методов построения многочлена г является следующий. В качестве многочлена rt можно взять интерполяционный многочлен Лагранжа, а в качестве узлов интерполяции — нули многочленов Чебышева первого рода, приведенные к отрезку [—Д —а]. Напомним (см. [1, с. 99], [46, с. 75]), что многочлены Чебышевапервого рода Тп на отрезке [—1,1] задаются формулой а их нулями являются точки Следовательно, для отрезка [—/3, — а] имеем следующие узлы интерполяции:
Использование в качестве узлов интерполяции именно нулей многочленов Чебышева целесообразно по следующим причинам. Во-первых, известно [1, с. 190], [20], что узлы интерполяции являются асимптотически оптимальными (в смысле оптимального поведения константы Лебега при п -+ со; константой Лебега называют L -норму проектора / i- р), если они совпадают с нулями многочленов Чебышева первого рода Тп. Во-вторых, в силу аналитичности функции / во всей комплексной плоскости рассматриваемые многочлены равномерно сходятся (см., например, [15, с. 94] и [47]) к / на любом компактном подмножестве комплексной плоскости. Наконец, в качестве многочлена rt, приближающего функцию /, молено взять многочлен наилучшего приближения. Говорят (см., например, [1, с. 141], [19, с. 18]), что непрерывные функции ірії [а,Ъ] — R, і = 0,1,..-п, образуют систему Чебышевау если при любых постоянных с ЄМ, і = 0,1,... п, среди которых хотя бы одна отлична от нуля, обобщенный многочлен р(Х) = =0 (А) имеет на [а, 6] не более п нулей. Говорят, что обобщенный многочлен р является многочленом наилучшего приближения для данной функции / на [а, 6], если на нем достигается инфимум величины зирдфм /(А) — р(Х) \ по всем многочленам р данной степени п. Предложение 15 ([19, с. 27]). Для любой системы Чебышева на отрезке и любой непрерывной функции f существует единственный {обобщенный) многочлен р ее наилучшего приближения. Следующий критерий обычно называют теоремой Чебышева об аль-тернансе. Предложение 16 ([1, с. 144], [15, с. 274], [19, с. 28]). Для того чтобы {обобщенный) многочлен р был многочленом наилучшего приближениям для f на [а, Ь], необходимо и достаточно, чтобы нашлось не менее п + 2 точек, в которых разность f — р принимает свое максимальное по модулю значение с чередующимися знаками. Для нас важны следующие два конкретных следствия этого утверждения. В первом описывается критерий оптимальности оценок (1.16) и (1.17), а во втором — оценок (1.23), (1.24), (1.26) и (1.27). Следствие 17. Пусть t фиксировано. (a) Для того чтобы многочлен А н- г (А) степени п был многочленом наилучшего приближения для функции А и-»- еЛ на отрезке [а,Ь], необходимо и достаточно, чтобы на отрезке [а, Ъ] нашлось не менее п + 2 точек, в которых разность А \-у р{Х, і) — еЛ( принимает свое максимальное по модулю значение с чередующимися знаками. (b) Для того чтобы многочлен А »- г((А) степени п был многочленом наилучшего приближения для функции A i-J- $t(A) на отрезке [а, Ь\, необходимо и достаточно, чтобы на отрезке [а, Ъ] нашлось не менее тг + 2 точек, в которых разность A i + rt{\) — #t(A) принимает свое максимальное по модулю значение с чередующимися знаками.
Факторизация уравнения второго порядка с помощью матричного корня обратного пучка
Рассмотрим теперь квадратное уравнение и его корни Предложение 29. Пусть выполнено условие (2.5). Тогда матричное уравнение имеет единственное решение V, удовлетворяющее условию VJ vi. Более того, V может получено как предел последовательности Доказательство. Доказательство практически дословно повторяет до казательство предложения 25. П Предложение 30. Пусть V произвольное решение уравнения (2.11). Тогда дифференциальное уравнение (2.3) эквивалентно уравнению Доказательство. Поскольку М + V + AV2 — 0, имеем Предложение 31. Пусть где т Є (иі г) — произвольное число, р = 1/т. Тогда матрицы G и F являются DG l-самосопряженными, причем G положительно определена. Более того, матрица V = FG l есть решение уравнения (2.11), и ее спектр лежит в {(і: \it\ v\ = 1/ }. Доказательство. Очевидно, это частный случай предложения 27. Отметим лишь, что интегралы связаны заменой переменных ц = 1/А. Следствие 32. Матрица V = FG 1 является самосопряженной и неположительно определенной относительно DG l-скалярного произведения, т.е. при всех р,ф ЄІ" Спектр матрицы V содержится в [—v\,0] и существует DG l-opmo-нормированный базис, состоящий из собственных векторов матрицы V. Доказательство. Доказательство полностью аналогично доказатель Ниже нам понадобится вспомогательное уравнение Его роль аналогична роли уравнения (2.8). Предложение 33. Пусть Z — решение уравнения (2.6) со спектром, лежащим в круге {А: А Zi}t а V — решение уравнения (2.11) со спектром, лежащим в {/І: р г і}. {Напомним, что \[v\ = z% и 1/zx = г?2-) Тогда матрицы MZ + 1 и AV + 1 обратимы и еде запись V означает матрицу, сопряженную к V относительно D-скалярного произведения. Доказательство. Доказательство принадлежит Курбатову В.Г., поэто му здесь оно не приводится, см. по этому поводу статью [85]. Следствие 34. Решение Z уравнения (2.6) со спектром, лежащим в круге {А: [А] z\}, является единственным; и решение V уравнения (2.11) со спектром, лежащим в круге { : // vi}, также является единственным. Доказательство. Действительно, представления Z ,D — A{AV -\-I)-1 и V yD = —M(MZ+ I)-1 справедливы для любой пары корней Z и V. Поэтому, если зафиксировать один из корней, то второй определяется единственным образом. П Предложение 35. ДЛЯ любой матрицы Т имеем T iD = D lrT D.
Доказательство. Действительно, Доказательство. Эквивалентность трех представлений для Q и 13 следует из предложений 33 и 35, а также определений М = D lN и А = D B. Формула (2.17) вытекает из преобразований Формула (2.18) доказывается аналогично- Предложение 37. Имеют место формулы ПриСледствие 32. Матрица V = FG 1 является самосопряженной и неположительно определенной относительно DG l-скалярного произведения, т.е. при всех р,ф ЄІ" Спектр матрицы V содержится в [—v\,0] и существует DG l-opmo-нормированный базис, состоящий из собственных векторов матрицы V. Доказательство. Доказательство полностью аналогично доказатель Ниже нам понадобится вспомогательное уравнение Его роль аналогична роли уравнения (2.8). Предложение 33. Пусть Z — решение уравнения (2.6) со спектром, лежащим в круге {А: А Zi}t а V — решение уравнения (2.11) со спектром, лежащим в {/І: р г і}. {Напомним, что \[v\ = z% и 1/zx = г?2-) Тогда матрицы MZ + 1 и AV + 1 обратимы и еде запись V означает матрицу, сопряженную к V относительно D-скалярного произведения. Доказательство. Доказательство принадлежит Курбатову В.Г., поэто му здесь оно не приводится, см. по этому поводу статью [85]. Следствие 34. Решение Z уравнения (2.6) со спектром, лежащим в круге {А: [А] z\}, является единственным; и решение V уравнения (2.11) со спектром, лежащим в круге { : // vi}, также является единственным. Доказательство. Действительно, представления Z ,D — A{AV -\-I)-1 и V yD = —M(MZ+ I)-1 справедливы для любой пары корней Z и Vведем, наконец, формулы, выражающие решение уравнения (2.1) через специальные решения уравнений (2.8) и (2.14). Теорема 38. Предположим, что обратная к матрице (2.15) существует. Тогда решение начальной задачи на оси для уравнения (2.1) можно представить в виде
Сведение дифференциального уравнения второго порядка к двум уравнениям первого порядка
Доказательство. Очевидно, это частный случай предложения 27. Отметим лишь, что интегралы связаны заменой переменных ц = 1/А. Следствие 32. Матрица V = FG 1 является самосопряженной и неположительно определенной относительно DG l-скалярного произведения, т.е. при всех р,ф ЄІ" Спектр матрицы V содержится в [—v\,0] и существует DG l-opmo-нормированный базис, состоящий из собственных векторов матрицы V. Доказательство. Доказательство полностью аналогично доказатель Ниже нам понадобится вспомогательное уравнение Его роль аналогична роли уравнения (2.8). Предложение 33. Пусть Z — решение уравнения (2.6) со спектром, лежащим в круге {А: А Zi}t а V — решение уравнения (2.11) со спектром, лежащим в {/І: р г і}. {Напомним, что \[v\ = z% и 1/zx = г?2-) Тогда матрицы MZ + 1 и AV + 1 обратимы и еде запись V означает матрицу, сопряженную к V относительно D-скалярного произведения. Доказательство. Доказательство принадлежит Курбатову В.Г., поэто му здесь оно не приводится, см. по этому поводу статью [85]. Следствие 34. Решение Z уравнения (2.6) со спектром, лежащим в круге {А: [А] z\}, является единственным; и решение V уравнения (2.11) со спектром, лежащим в круге { : // vi}, также является единственным. Доказательство. Действительно, представления Z ,D — A{AV -\-I)-1 и V yD = —M(MZ+ I)-1 справедливы для любой пары корней Z и V. Поэтому, если зафиксировать один из корней, то второй определяется единственным образом. П Предложение 35. ДЛЯ любой матрицы Т имеем T iD = D lrT D. Доказательство. Действительно, Доказательство. Эквивалентность трех представлений для Q и 13 следует из предложений 33 и 35, а также определений М = D lN и А = D B. Формула (2.17) вытекает из преобразований Формула (2.18) доказывается аналогично- Предложение 37. Имеют место формулы Приведем, наконец, формулы, выражающие решение уравнения (2.1) через специальные решения уравнений (2.8) и (2.14). Теорема 38. Предположим, что обратная к матрице (2.15) существует. Тогда решение начальной задачи на оси для уравнения (2.1) можно представить в виде х = хг + х2, где х\ и Х2 — решения уравнений х[ - Zxl = a D-Hu = {NZ + ZN + D) 4u и Vx 2 x2 = {AV + l)-xMSl-lD-lbu = (BV D)-lN(NZ + ZN + D) 1bu. Доказательство. Подставим x = x\ +x% в левую часть уравнения (2.3), равносильного (2.1), и воспользуемся представлениями (2.7) и (2.12): Итак, видно, что х удовлетворяет уравнению (2.3). П Напомним, что уравнение (2.1) имеет вариант (2.2). Приведем для его решения явную вычислительную формулу. Теорема 39. Предположим, что обратная к матрице (2.16) существует. Тогда решение начальной задачи на оси для уравнения (2.2) можно представить в виде где х\ и Х2 — решения уравнений Доказательство. Подставим x = xi + X2 в левую часть уравнения равносильного (2.2), и воспользуемся представлениями (2.7) и (2.12): Мх" Тем самым, x удовлетворяет уравнению (2.21). Итак, уравнения (2.1) и (2.2) свелись к уравнениям первого порядка со специальными векторами Ь\, Ъ% и исходной функцией и.
Обсудим способы решения этих уравнений. Как обычно, будем предполагать, что функция г непрерывна и равна нулю левее некоторой точки t$. Кроме того, будем считать, что нас интересуют не сами решения уравнений (2.22) и (2.23), а функция вида (ср. (1.2)) Тем самым, вместо решений уравнений (2.22) и (2.23) нас будут интересовать функции В силу предложения 4 (при D = 1, В = — Z) функцию yi(t) = {xx(t),l) можно искать с помощью формулы Если порядок матрицы Z является небольшим, то для построения решения по этой формуле можно воспользоваться предложением б (напомним, что в силу следствия 28 спектр матрицы Z содержится в [—zi,0] и является простым). В случае, когда порядок матрицы Z велик, можно воспользоваться методом, описанным в 1.2- 1.4. Чтобы получить приближенное решение уравнения (2.22), выберем функцию rt, приближающую функцию i?t на спектре матрицы Z, а затем в качестве приближенного решения возьмем функцию Поскольку в силу следствия 28 спектр матрицы Z содержится в [—z\y 0], в качестве rf разумно брать многочлен или рациональную функцию, приближающую fit на [—2i,0]. Обсудим теперь способ решения уравнения (2.23). Рассмотрим вспомогательную функцию Заметим, что ядро P\(s) = ХеХз интеграла обладает следующими свойствами: ip\(s) 0 при s [0,-foo) и А О, А Ж; /0+ AeAi ds = —1 при А О, А Є R; для любого є 0 и любого Предложение 41. Решение начальной задачи на оси для уравнения (2.23) можно представить в виде Доказательство. В силу следствия 32 спектр матрицы V содержится в (—оо,0] и существует базис из собственных векторов матрицы V. В базисе из собственных векторов матрицы V уравнение (2.23) распадается на независимые скалярные уравнения где /xjt — собственные значения матрицы V, а4 . Ь2 — координаты векторов Х2 и 6г соответственно. Если щ. 0, то имеем решение (ср. с предложением 4)
Функции от факторизованного матричного пучка
В этом пункте приводятся явные формулы для нахождения функций от пучка. Основными результатами являются теоремы 46 и 47. Теорему 46 можно интерпретировать как аналог представления функций от одной матрицы с помощью ее жордановой формы [53]. В теореме 47 предлагается альтернативная формула для F(A,В), в которой используются функции от матриц А и В, а также специальная матрица Р, которая зависит только от Л и В и не зависит от /. Эту теорему удобно использовать, если необходимо вычислить несколько функций F от одного и того же пучка. Рассмотрим матричную функцию вида где контур Г охватывает спектр квадратичного матричного пучка А м-(А1 — #)(Л1 — А). Здесь А и В — квадратные матрицы, а / — функция, аналитическая в окрестности множества, границей которого является контур Г. Обсудим сначала случай, когда спектры матриц А и В не пересекаются. Будем предполагать, что матрицы Л и В могут иметь кратные собственные значения. Жордановы разложения [14, 53] матриц ЛиВ имеют вид — матрицы соответствующего размера. Очевидно, Nk — 0 при к г, где г — порядок матрицы N. Под № будем подразумевать единичную матрицу 1. Нетрудно видеть, что столбцы матриц Т и Н являются собственными и присоединенными векторами матриц А и В. Отметим, что в ситуации общего положения присоединенных векторов не бывает, т.е. все блоки имеют размер 1 X 1, а матрицы N являются нулевыми матрицами размера 1x1. Обозначим через S матрицу Т гН. Разобьем ее на прямоугольные блоки S = {Sij} размера r$ х qj, где г» — порядок жорданова блока АІ матрицы A, a qj — порядок жорданова блока Bj матрицы В. Нам понадобится также вспомогательная блочная матрица L = { }» Здесь С+ — биномиальные коэффициенты, а и j3j — собственные значения матриц Лий, отвечающие жордановым блокам АІ и Bj соответственно. Размер блоков L(j равен г І X qj. В частном случае, когда матрицы диагонализуемы (и, соответственно, все матрицы N = Nixi = 0), ситуация упрощается, и матрица L состоит из скалярных элементов (не блоков) где Sij — скалярные элементы матрицы S. Предложение 45 ([14]). Пусть матрица Л представляет собой один оюорданов блок порядка г, на главной диагонали которого стоят числа а, тогда Доказательство. Следует из формулы Неймана [36]. Напомним, что для функции от одной матрицы справедлива [14, 53] формула (3.19) позволяющая легко находить /(-4), если известна жорданова форма матрицы А. Следующий результат можно рассматривать как аналог формулы (3.19) для случая матричного пучка. Теорема 46. Пусть спектры матриц А и В не пересекаются. Тогда для F(A, В) справедливо представление Ниже в теореме 47 приведена еще одна формула для функции от фак-торизованного матричного пучка. В этом представлении для F(A, В) используются функции от матриц А и В, а также некоторая матрица Р, которая зависит только от А и В и не зависит от /.
Теорема 47. Пусть спектры матриц Л и В не пересекаются. Тогда а матрица L определяется формулой (3.17). Доказательство. Перепишем правую часть равенства (3.22) в виде Рассмотрим отдельно случай, когда матрицы А и В молено привести к диагональному виду. Это условие заведомо выполнено в ситуации общего положения, для самосопрялсенных матриц, а также для пучка, все собственные значения которого являются простыми. В этом случае формулы (3.20) и (3.22) упрощаются. Теорема 48. Пусть матрицы Л и В диагонализуемы и не имеют общих собственных значений. Тогда P = TLH lf матрица L определяется формулой (3.18), a Sij — скалярные элементы (не блоки) матрицы S = Т гН. Доказательство. Следует из теорем 46 и 47. Предположим теперь, что спектры матриц Ли В пересекаются. Аналоги формул (3.20) и (3.22) для случая одного общего собственного значения для диагонализуемых матриц приведены в теореме 49. Они очевидным образом обобщаются на случай большего числа совпадений. Если хотя бы одна из матриц А или В имеет нетривиальные жордановы блоки, формулы становятся еще более громоздкими. Этот случай обсуждается в теореме 50. Предельная ситуация, когда матрицы Л is. В совпадают, рассматривается в теореме 51. Теорема 49. Пусть матрицы Л и В диагонализуемы и имеют одно общее собственное значение, т.е щ0 = /?J0 для некоторых г0 и jo Тогда где матрица Q l получается из матрицы Q (см. теорему 4$) заменой элемента с индексом IQJQ на число и? = Wi0j0 = Sioj0f (ai0); матрица Qbofa] получается из нулевой матрицы заменой элемента с индексом Wo на ш7 РІ 0 0 — Ті№іІН 1, а матрица ІА 0 0 получается из матрицы L (см. формулу (3.18)) в результате обнуления элемента с индексом го Jo
