8.2.3. Вырожденные и плохо обусловленные системы
Вернемся вновь к СЛАУ Aх=b
с квадратной матрицей А размера MхN
, которая, в отличие от рассмотренного выше "хорошего" случая (см. разд. 8. Г), требует специального подхода. Обратим внимание на два похожих типа СЛАУ:
- вырожденная система (с нулевым определителем |А|=0 );
- плохо обусловленная система (определитель А не равен нулю, но число обусловленности очень велико).
Несмотря на то, что эти типы систем уравнений существенно отличаются друг от друга (для первого решение отсутствует, а для второго существует и единственно), с практической точки зрения вычислителя, между ними много общего.
Вырожденные СЛАУ
Вырожденная система - это система, описываемая матрицей с нулевым определителем
|A|=0
(сингулярной матрицей). Поскольку некоторые уравнения, входящие в такую систему, представляются линейной комбинацией других уравнений, то фактически сама система является недоопределенной. Несложно сообразить, что, в зависимости от конкретного вида вектора правой части ь, существует либо бесконечное множество решений, либо не существует ни одного. Первый из вариантов сводится к построению нормального псевдорешения (т. е. выбора из бесконечного множества решений такого, которое наиболее близко к определенному, например, нулевому, вектору). Данный случай был подробно рассмотрен в разд. 8.2.2 (см. листинги 8.11-8.13).
Рис. 8.7. Графическое представление несовместной системы двух уравнений с сингулярной матрицей
Рассмотрим второй случай, когда СЛАУ Aх=b с сингулярной квадратной матрицей А не имеет ни одного решения. Пример такой задачи (для системы двух уравнений) проиллюстрирован рис. 8.7, в верхней части которого вводятся матрица А и вектор b , а также предпринимается (неудачная, поскольку матрица А сингулярная) попытка решить систему при помощи функции isolve . График, занимающий основную часть рисунка, показывает, что два уравнения, задающие систему, определяют на плоскости (x0,xi) две параллельные прямые. Прямые не пересекаются ни в одной точке координатной плоскости, и решения системы, соответственно, не существует.
ПРИМЕЧАНИЕ
Во-первых, заметим, что СЛАУ, заданная несингулярной квадратной матрицей размера 2x2, определяет на плоскости пару пересекающихся прямых (см. рис. 8.9 ниже). Во-вторых, стоит сказать, что, если бы система была совместной, то геометрическим изображением уравнений были бы две совпадающие прямые, описывающие бесконечное число решений.
Рис. 8.8. График сечений функции невязки f (х) = |Ах-b|
Несложно догадаться, что в рассматриваемом сингулярном случае псевдорешений системы, минимизирующих невязку |Ax-b| , будет бесконечно много, и лежать они будут на третьей прямой, параллельной двум показанным на рис. 8.7 и расположенной в середине между ними. Это иллюстрирует рис. 8.8, на котором представлено несколько сечений функции f(x)= = | Аx-b | , которые говорят о наличии семейства минимумов одинаковой глубины. Если попытаться использовать для их отыскания встроенную функцию Minimize , ее численный метод будет всегда отыскивать какую-либо одну точку упомянутой прямой (в зависимости от начальных условий). Поэтому для определения единственного решения следует выбрать из всего множества псевдорешений то, которое обладает наименьшей нормой. Можно пытаться оформить данную многомерную задачу минимизации в Mathcad при помощи комбинаций встроенных функций Minimize , однако более эффективным способом будет использование регуляризации (см. ниже) или ортогональных матричных разложений (см. разд. 8.3).
Плохо обусловленные системы
Плохо обусловленная система - это система, у которой определитель А не равен нулю, но число обусловленности |А -1 | |А| очень велико. Несмотря на то, что плохо обусловленные системы имеют единственное решение, на практике искать это решение чаще всего не имеет смысла. Рассмотрим свойства плохо обусловленных СЛАУ на двух конкретных примерах (листинг 8.14).
Листинг 8.14
. Решение двух близких плохо обусловленных
СЛАУ
Каждая строка листинга 8.14 содержит решение двух очень близких плохо обусловленных СЛАУ (с одинаковой правой частью ь и мало отличающимися матрицами А). Несмотря на близость, точные решения этих систем оказываются очень далекими друг от друга. Надо заметить, что для системы двух уравнений точное решение получить легко, однако при решении СЛАУ большой размерности незначительные ошибки округления, неминуемо накапливаемые при расчетах (в том числе и "точным" алгоритмом Гаусса), приводят к огромным погрешностям результата. Возникает вопрос: имеет ли смысл искать численное решение, если заранее известно, что, в силу неустойчивости самой задачи, оно может оказаться совершенно неправильным?
Еще одно соображение, которое вынуждает искать специальные методы решения плохо обусловленных СЛАУ (даже приведенной в качестве примера в листинге 8.14 системы двух уравнений), связано с их физической интерпретацией как результатов эксперимента. Если изначально известно, что входные данные получены с некоторой погрешностью, то решать плохо обусловленные системы не имеет вовсе никакого смысла, поскольку малые ошибки модели (матрицы А и вектора ь) приводят к большим ошибкам решения. Задачи, обладающие такими свойствами, называются некорректными.
Чтобы лучше понять причину некорректности, полезно сравнить графическую интерпретацию хорошо (рис. 8.9) и плохо (рис. 8.10) обусловленной системы двух уравнений. Решение системы визуализируется точкой пересечения двух прямых линий, изображающих каждое из уравнений.
Рис. 8.9. График хорошо обусловленной системы двух уравнений
Рис. 8.10. График плохо обусловленной системы двух уравнений
Из рис. 8.10 видно, что прямые, соответствующие плохо обусловленной СЛАУ, располагаются в непосредственной близости друг от друга (почти параллельны). В связи с этим малые ошибки в расположении каждой из прямых могут привести к значительным погрешностям локализации точки их пересечения (решения СЛАУ) в противоположность случаю хорошо обусловленной системы, когда малые погрешности в наклоне прямых мало повлияют на место точки их пересечения (рис. 8.9).
ПРИМЕЧАНИЕ
Плохая обусловленность матрицы типична и при реконструкции экспериментальных данных, задаваемых переопределенными (несовместными) СЛАУ (например, в задачах томографии). Именно такой случай приведен в качестве примера в следующем разделе (см. листинг 8.16 ниже).
Метод регуляризации
Для решения некорректных задач, в частности, вырожденных и плохо обусловленных СЛАУ, разработан очень эффективный прием, называемый регуляризацией. В его основе лежит учет дополнительной априорной информации о структуре решения (векторе априорной оценки хо), которая очень часто присутствует в практических случаях. В связи с тем, что регуляризация была подробно рассмотрена в разд. 6.3.3, напомним лишь, что задачу решения СЛАУ Аx=b можно заменить задачей отыскания минимума функционала Тихонова:
Ω (х,λ ) = |Ах-b| 2 +λ |х-х0| 2 . (8.3)
Здесь Я, - малый положительный параметр регуляризации, который необходимо подобрать некоторым оптимальным способом. Можно показать, что задачу минимизации функционала Тихонова можно, в свою очередь, свести к решению другой СЛАУ:
(А T А+λ I)-х=А T В+λ х0, (8.4)
Которая при λ ->0 переходит в исходную плохо обусловленную систему, а при больших x , будучи хорошо обусловленной, имеет решение х 0 . Очевидно, оптимальным будет некоторое промежуточное значение А , устанавливающее определенный компромисс между приемлемой обусловленностью и близостью к исходной задаче. Отметим, что регуляризационный подход сводит некорректную задачу к условно-корректной (по Тихонову) задаче отыскания решения системы (8.4), которое, в силу линейности задачи, является единственным и устойчивым.
Приведем, без излишних комментариев, регуляризованное решение вырожденной системы, которая была представлена на рис. 8.8. Листинг 8.15 демонстрирует отыскание решения задачи (8.4), а полученная зависимость невязки и самого решения от параметра регуляризации Я показана на рис. 8.11 и 8.12 соответственно. Важно подчеркнуть, что решения исходной системы и, следовательно, системы (8.4) при λ =0 не существует.
Листинг 8.15
.Регуляризация вырожденной СЛАУ
Заключительным шагом регуляризации является выбор оптимального λ . Имеется, по крайней мере, два соображения, исходя из которых, можно выбрать параметр регуляризации, если опираться на зависимость от него невязки. В рассматриваемом примере применим критерий определения λ , опирающийся -на подбор нормы невязки, равной априорной оценке погрешностей задания входных данных: матрицы А и вектора b , т. е. величине |δA | + |5λ| . Например, можно выбрать норму невязки и, соответственно, параметр λ и решение х(λ ), которые отмечены на рис. 8.11 и 8.12 пунктирами.
ПРИМЕЧАНИЕ 1
Другим вариантом выбора λ , не требующим никаких априорных соображений относительно погрешностей модели, является так называемый квазиоптимальный метод, рассмотренный в разд. 6.3.3.
ПРИМЕЧАНИЕ 2
Полезно убедиться в том, что формула (8.4) в случае линейной задачи дает тот же результат, что и общая формула (8.3). Для этого достаточно изменить в листинге 8.15 последнюю строку, выражающую решение СЛАУ (8.4), на код, реализующий минимизацию численным методом, как это показано в листинге 8.16. Расчеты (которые требуют значительно большего компьютерного времени) дают тот же самый результат, что был приведен на рис. 8.11 и 8.12.
ПРИМЕЧАНИЕ 3
Попробуйте в расчетах листингов 8.15 и 8.16 взять иную, например, более реалистичную, априорную оценку решения (вместо использованного в них нулевого вектора х0) и посмотреть, как это повлияет на результат.
Рис. 8.11. Зависимость невязки регупяризованного решения вырожденной СЛАУ от параметра А. (продолжение листинга 8.15)
ПРИМЕЧАНИЕ 4
Любопытно также применить вместо формулы (8.3) в качестве функционала Тихонова другую зависимость: Ω (х, λ ) = |Ах-b|+ λ |х-х0 | . Соответствующий пример расчетов вы найдете на компакт-диске.
Рис. 8.12. Регуляризованное решение в зависимости от λ (продолжение листинга 8.15)
Листинг 8.16.
Регуляризация СЛАУ при помощи алгоритма минимизации (продолжение листинга 8.15)
УДК 519.61:621.3
В.П. ВОЛОБОЕВ*, В.П. КЛИМЕНКО*
ОБ ОДНОМ ПОДХОДЕ К РЕШЕНИЮ ПЛОХО ОБУСЛОВЛЕННОЙ СИСТЕМЫ ЛИНЕЙНЫХ АЛГЕБРАИЧЕСКИХ УРАВНЕНИЙ, ОПИСЫВАЮЩИХ ФИЗИЧЕСКИЙ ОБЪЕКТ
Институт проблем математических машин и систем НАН Украины, Киев, Украина
Анотація. Показано, що вірогідність результатів моделювання фізичних об"єктів, дискретна модель яких описується системою лінійних алгебраїчних рівнянь (СЛАР), залежить не від поганої обумовленості матриці, а від некоректного вибору змінних СЛАР на етапі складання рівнянь методом вузлових потенціалів або його аналогів, а сам метод є окремий випадок методу коректної постановки завдання. Запропоновано методику перевірки на коректність СЛАР, складеної методом вузлових потенціалів, що має невироджену й симетричну матрицю, і якщо необхідно перетворення її до коректного виду.
Ключові слова: система, моделювання, некоректне завдання, погана обумовленість, система лінійних алгебраїчних рівнянь, метод вузлових потенціалів, метод коректної постановки завдання, перевірка на коректність.
Аннотация. Показано, что достоверность результатов моделирования физических объектов, дискретная модель которых описывается системой линейных алгебраических уравнений (СЛАУ), зависит не от плохой обусловленности матрицы, а от некорректного выбора переменных СЛАУ на этапе составления уравнений методом узловых потенциалов или его аналогами, а сам метод есть частный случай метода корректной постановки задачи. Предложена методика проверки на корректность СЛАУ, составленной методом узловых потенциалов, имеющей невырожденную и симметричную матрицу, и если необходимо преобразование её к корректному виду.
Ключевые слова: система, моделирование, некорректная задача, плохая обусловленность, система линейных алгебраических уравнений, метод узловых потенциалов, метод корректной постановки задачи, проверка на корректность.
Abstract. The paper shows that reliability of results of simulation of the physical objects, which discrete model is described by a system of the linear algebraic equations (SLAE) depends not on poor-conditioned matrix but on an incorrect choice of variable SLAE at generation of equations stage by a node potential method or its analogues, and the method is a special case of a method of correct statement of a problem. It was suggested the check-out method on a correctness of SLAE, made by a node potential method, having nonsingular and a symmetric matrix and if it is necessary its transformation to a correct form.
Keywords: system, simulation, incorrect problem, poor-conditioned, system of the linear algebraic equations, node potential method, method of correct statement of a problem, check-out on a correctness.
1. Введение
Многие задачи моделирования физических (технических) объектов сводятся к решению систем линейных алгебраических уравнений (СЛАУ). Поскольку все вычисления при решении таких систем выполняются с конечным числом значащих цифр, то точность может значительно теряться из-за ошибок округления. Плохо обусловленной (неустойчивой) системой или в более общей формулировке - некорректно поставленной задачей принято считать ту задачу, которая при фиксированном уровне ошибок входных данных и точности вычислений не гарантирует в решении никакой точности. В качестве априорной наихудшей оценки возможных ошибок решения СЛАУ используется число обусловленности . Как следует из литературы, разработка методов решения некорректно поставленных задач рассматривается как чисто математическая задача, в которой не учитываются особенности физических (технических) объектов, несмотря на то, что численное решение многих задач математической физики и математического моделирования сложных физических процес-
© Волобоев В.П., Клименко В.П., 2014
сов и технических систем является неиссякаемым источником задач линейной алгебры. Для перечисленного класса задач при разработке методов решения не рассматривается этап составления СЛАУ, на котором тем или иным способом можно учесть особенности конкретной задачи. В том, что этот этап необходимо учитывать, подтверждают результаты следующих работ.
Прежде всего, следует отметить работу , где приведены примеры матриц, для которых потеря точности при решении СЛАУ невелика, а величина числа обусловленности огромна, то есть показано, что общепринятый критерий априорной оценки точности решения СЛАУ по числу обусловленности есть необходимый, но недостаточный. Совершенно новый подход к решению некорректно поставленной задачи предложен в работах . Он заключается в том, что с целью повышения точности решения СЛАУ даже при большой величине числа обусловленности на этапе описания дискретной модели физического объекта предлагается корректно составлять СЛАУ. Это означает не только то, что такие матрицы существуют, как об этом сообщалось в работе , но и то, что предложен метод корректного составления матрицы СЛАУ, описывающей дискретную модель объекта. Метод составления матрицы СЛАУ рассмотрен применительно к задачам моделирования поведения электрических цепей , энергосистем , стержневых систем механики и эллиптических уравнений матфизики .
Суть данного метода заключается в том, что, в отличие от существующих методов, при формировании СЛАУ целенаправленным выбором переменных учитываются параметры дискретной модели физического объекта. Следует заметить, что метод применим только к тем объектам, топология дискретной модели которых представлена графом.
Этому требованию удовлетворяет расчетная модель электрической цепи и энергосистемы. Для многих задач математического моделирования сложных физических процессов, технических систем и математической физики представление топологии дискретной модели в виде графа не применяется. В работах показано, что вышеприведенное ограничение снимается за счет представления топологии элементов расчетных схем дискретной модели физического объекта в виде графа. Там же приведена методика представления топологии элементов в виде графов.
В данной работе будет предложен метод корректировки некорректно поставленной задачи для случая, когда топология дискретной модели не представлена в виде графа. При разработке метода учитывается тот факт, что общепринятый метод описания дискретных моделей задач математической физики и сложных физических процессов и технических систем (метод узловых потенциалов ) есть частный случай метода корректного составления матрицы СЛАУ.
2. Связь между точностью решения СЛАУ, описывающей дискретную модель объекта, и методом составления уравнений
Академик Воеводин В.В. в работе показал, что наибольшая точность результатов решения СЛАУ методом Гаусса достигается в случае применения метода с выбором главного элемента. На основе этой идеи было опубликовано огромное число работ. Однако решение практических задач показало, что точность решения СЛАУ, особенно в случае плохо обусловленных матриц, значительно теряется из-за ошибок округления, то есть для повышения точности результатов на этапе решения недостаточно одного применения метода Гаусса с выбором главных элементов.
Дальнейшим развитием этой идеи есть метод, предложенный в работе , где предлагается на этапе составления описания дискретной модели объекта формировать диагональные элементы матрицы как главные. Для этого при составлении описания используется дополнительная информация, а именно, параметры дискретной модели. Эффективность данного подхода, а именно, зависимость точности решения СЛАУ, описывающей дискрет-
ISSN 1028-9763. Математичні машини і системи, 2014, № 4
ную модель объекта, от метода составления уравнений будет продемонстрирована на модельном примере. Ниже будет рассмотрено составление описания модельного примера методом, описанным в , с выбором главного элемента и без выбора и его решение.
В качестве модельного примера выбрана электрическая цепь, приведенная на рис. 1.
Рис. 1. Электрическая цепь
Известно , что обусловленность СЛАУ, описывающей электрическую цепь, зависит от диапазона разброса величин проводимостей (сопротивлений) компонент цепи. Выбранный диапазон изменения проводимостей компонент электрической цепи, равный 15 порядкам, обеспечивает плохую обусловленность СЛАУ и тем самым, как принято считать, некорректность задачи. На примере расчета потенциала узла 2 (напряжение на компоненте G2) будет анализироваться зависимость достоверности результатов расчета от способа формирования диагонального элемента при составлении описания электрической цепи.
Ниже приведены основные положения, необходимые для решения модельного примера методом корректной постановки задачи . Построение математической модели электрической цепи данным методом базируется на основной системе уравнений электрической цепи, куда входят компонентные уравнения и уравнения, составленные на основе законов Кирхгофа. Для модельного примера компонентное уравнение имеет вид
где U i - напряжение, падающее на компоненте, I - ток, протекающий через компоненту, Gt - проводимость компоненты.
Для описания графа электрической цепи и, соответственно, уравнений на основе законов Кирхгофа применяются топологические матрицы контуров и сечений. Граф цепи совпадает с электрической цепью. Составление топологических матриц контуров и сечений включает выбор дерева графа цепи и составление контуров для выбранного дерева. Дерево графа электрической цепи выбирается таким образом, чтобы все источники напряжения включались в дерево, а все источники тока в хорды. Элементы в векторах напряжений U и токов I компонент цепи группируются на входящие в дерево (индекс Д), то есть ветви и хорды (индекс X), таким образом:
Контуры образуются присоединением хорд к дереву графа схемы. В этом случае
топологическая матрица контуров имеет вид
где 1 - единичная подматрица хорд, t
Обозначает транспонирование матрицы, а топологическая матрица сечений - вид |1 -F , где 1 - единичная подматрица ветвей. Как следует из , диагональные члены матрицы
ISSN 1028-9763. Математичні машини і системи, 2014, № 4
будут главными в том случае, когда проводимости компонент дерева в контурах имеют максимальную проводимость. Учитывая вид топологических матриц, уравнения цепи, составленные на основе законов Кирхгофа, в матричном виде можно записать следующим образом:
их =-ґид, (3)
Переменные составляемой системы уравнений выбираются из напряжений и/или токов компонент в результате анализа основной системы уравнений. В случае выбора в качестве переменных напряжений компонент, входящих в ветви дерева, компонентные уравнения (1) и уравнения (3), (4) можно преобразовать к следующему виду:
Gд U д - F(Gx (- FUд)) = 0.
Ниже будет приведено составление уравнений для модельного примера. Вначале составляется описание электрической цепи таким образом, чтобы диагональные члены матрицы были главными. Этому требованию удовлетворяет набор компонент E1, G6, G3, G2, входящих в дерево (на рис. 1 ветви дерева выделены жирной линией). Выбранному дереву соответствуют следующие векторы напряжений, токов компонент:
и топологические матрицы
Уравнение (5), с учетом (6), (7) и компонентных уравнений после преобразований, имеет следующий вид:
- (G4 + G5) (G4 + G5) G1 + G2 + G4 + G5
СЛАУ (8) есть плохо обусловленная, так как собственные числа матрицы \= 1,5857864376253, Я2 = 5,0E +14+j5,0E +14, А, = 5,0E +14 - j5,0E +14. Для того, чтобы определить, как зависит точность результатов решения системы от выбора варианта составления уравнений, расчет потенциала Uq узла 2 будет выполнен в общем виде:
ISSN 1028-9763. Математичні машини і системи, 2014, № 4
(g1+g2 +g4 +g5)-
Из анализа вычислительного процесса (9-11) следует, что, несмотря на большой диапазон изменения величин проводимостей (15 порядков), нет жестких требований к конечной точности представления чисел как при составлении уравнений, так и при их решении. Для получения достоверного результата достаточно выполнить вычислительный процесс составления и решения СЛАУ с точностью представления чисел в две значащие цифры.
Следует заметить, что в СЛАУ (8) диагональный элемент второй строки (столбца) матрицы G+G4+G5I значительно больше (на 15 порядков) суммы остальных членов
строки (столбца) | G4 + 2G51. Это означает, что, приняв UG = 0, можно упростить СЛАУ
(8), сохранив достоверность результатов. В эпоху ручного счета этому приему соответствовало объединение узла 2 с 3 (рис. 1).
Во втором случае (без выбора диагонального элемента главным) достаточно выбрать в дерево компоненты Ex, G6, G4, G2 (на рис. 1 ветви дерева помечены прерывистой
линией). Падения напряжений на этих компонентах соответствуют узловым потенциалам 1, 4, 3, 2, отсчитываемым от нулевого узла. Это означает, что при таком выборе компонент в дерево метод корректного составления матрицы СЛАУ совпадает с методом узловых потенциалов. Выбранному дереву и хордам соответствуют следующие векторы напряжений, токов компонент:
U Д = UG UG G4 , Ux = G1 UG3 UG G Д G ig G4 , Ix = G1 IG3 IG
UG G2 G5 ig G2 G5
и топологических матриц
Уравнение (5), с учетом (12), (13) и компонентных уравнений, примет следующий
ISSN 1028-9763. Математичні машини і системи, 2014, № 4
G5 + G6 -G5 0 UG G6 0
G5 G3 + G4 + G5 -G3 Uo. = 0
0 - G3 G1 + G2 + G3 Uo2 G1E1
Система уравнений (14) есть плохая обусловленная, так как имеет следующие собственные числа матрицы: 1 = 1.0,1 =1015 +у1015,1 =1015-/1015 . Как в первом варианте примера, будет выполнен расчет потенциала UG узла 2 в общем виде:
(G + G + G)------------
V 3 4 У (G + G)
+ (G1 + G2 + G3)
3 4 5" (G5 + G6)
Из анализа вычислительного процесса решения системы уравнений (15-17) следует, что достоверность результатов зависит как при составлении, так и решении уравнений от конечной точности представления чисел. Так, если вычислительный процесс решения системы (15-17) выполняется с точностью меньше 15 значащих цифр, то результат будет
1015 +1015 ~ o,
а в случае, когда с точность больше 15 значащих цифр, будет
1030 + 2*1015 +1030 + %+ 3/1015)
Из сравнения матриц (8) и (14), а также вычислительных процессов решения систем уравнений, вытекают следующие выводы.
Метод узловых потенциалов есть частный случай метода, предложенного в , а именно, в методе узловых потенциалов всегда в дерево выбраны ребра графа, связывающие базовый узел с остальными.
Диагональные элементы матрицы по модулю больше остальных элементов, как в строках, так и столбцах, независимо от того, матрица составлена с или без выбора максимальных диагональных. Различие заключается только в том, насколько диагональные элементы больше недиагональных. Это означает, что решение такого типа СЛАУ методом Г аусса с выбором главного элемента не повышает точности результатов для данного класса задач.
ISSN 1028-9763. Математичні машини і системи, 2014, № 4
Конечное число значащих цифр, используемых при решении методом Г аусса, существенным образом зависит от того, матрица составлена с или без выбора максимальных диагональных элементов. Отличие одного варианта задачи от другого заключается только в том, что на этапе составления уравнений в одном случае компонента с максимальной проводимостью выбирается в дерево и тем самым напряжение этой компоненты выступает в качестве переменной составляемой СЛАУ. Проводимость этой компоненты участвует только при формировании диагонального элемента матрицы. В другом случае эта компонента попадает в хорды. Как следует из уравнения (3), напряжение компоненты определяется через напряжения компонент дерева. Из уравнения (4) следует, что проводимость компоненты участвует в формировании элементов строк и столбцов и тем самым проводимость хорды определяет величину этих элементов матрицы.
3. Преобразование матрицы СЛАУ, составленной методом узловых потенциалов, к виду, соответствующему корректной постановке
При численном решении задач математической физики и математического моделирования сложных физических процессов и технических систем для составления СЛАУ, описывающих дискретные модели этих задач, в основном применяется метод узловых потенциалов или его аналоги. Отличительной чертой этого метода есть то, что в качестве переменных СЛАУ применяются потенциалы расчетной схемы дискретной модели, отсчитываемые от базового узла к остальным узлам, простой алгоритм составления уравнений, слабо заполненная матрица СЛАУ. Платой за такую эффективность может быть некорректность поставленной задачи. Учитывая, что метод узловых потенциалов есть всего лишь один из вариантов метода корректной постановки задачи, то некорректно поставленную задачу можно откорректировать, применив преобразование матрицы. Ниже будет рассмотрен алгоритм преобразования некорректно составленной методом узловых потенциалов задачи.
Из всего многообразия физических объектов будут рассмотрены только те объекты, линейная дискретная модель которых описывается СЛАУ с невырожденной и симметричной матрицей.
3.1. Алгоритм преобразования матрицы
При разработке алгоритма преобразования матрицы используется тот факт, что j -ый недиагональный элемент i -ой строки матрицы входит в матрицу со знаком минус и содержит параметр дискретной модели, описывающий связь между i -ым и j -ми узлами дискретной модели. Диагональный элемент входит в матрицу с положительным знаком, содержит сумму недиагональных элементов и параметр дискретной модели, описывающий связь между i -ым узлом и базовым. Обычно при нумерации узлов дискретной модели базисный узел принято считать нулевым.
Как следует из исследования, выполненного выше, некорректность поставленной задачи на уровне составленной СЛАУ возникает только в том случае, если хотя бы один из недиагональных элементов строки будет значительно больше параметра дискретной модели, входящего только в диагональный элемент. Ниже приводится методика проверки на корректность составленной СЛАУ.
Пусть СЛАУ имеет вид
где x - вектор узловых потенциалов (узловых воздействий), у - вектор внешних потоков, A - матрица вида
ISSN 1028-9763. Математичні машини і системи, 2014, № 4
а11 а1і a1j a1n
аі1 а,і aj ain , (21)
aJ1 an1 аі aJJ ann
где n - размер матрицы. Элементы матрицы удовлетворяют следующим требованиям:
аи > 0, а. < 0, а. = а]г,1 < i < n, 1 < j < n при j Ф і. (22)
Ниже будет рассмотрена проверка на корректность і -ой строки матрицы и, если необходимо, её корректировка.
Прежде всего определяется параметр дискретной модели ait, входящий только в диагональный элемент і -ой строки матрицы,
Считается, что і -ая строка матрицы корректно составлена, если параметр ait удовлетворяет условию
1 < j < n, при j Ф і.
При невыполнении условия (24) выполняется корректировка і -ой строки. Вначале из недиагональных элементов выбирается наибольший. Пусть это будет j -ый элемент і -ой строки. Нетрудно убедиться, что в силу специфики составления матрицы (условие (22)) параметр дискретной модели, который участвует в формировании элементов о. и а.^ і -ой и j -ой строк, входит как составная часть в элементы aii и а. . Суть корректировки і -ой строки заключается в преобразовании і -ой и j -ой строк матрицы таким образом, чтобы величина элемента а. входила только в элемент aii. Нетрудно убедиться, что, представив переменную xi в виде
X = xj + xj (25)
и выполнив следующее преобразование элементов j -ого столбца матрицы СЛАУ
о = аі. + ai, 1 < 1 < n , (26)
получим новый j -ый столбец матрицы, в котором преобразованные элементы а. и а. не содержат параметр дискретной модели, формировавший элементы а. и а. .
На следующем шаге выполняется преобразование j -ой строки по формуле
aji = а.і + аіі, 1 < l < n . (27)
Элементы a і преобразованной j -строки уже не содержат параметр дискретной модели, соответствующий элементу а і .
ISSN 1028-9763. Математичні машини і системи, 2014, № 4
Проверка на корректность матрицы СЛАУ и корректировка некорректных строк выполняются для всей матрицы. В данной работе рассмотрен только подход к построению алгоритма преобразования матрицы к корректному виду. Вопросы, связанные с разработкой эффективно работающего алгоритма преобразования матрицы к корректному виду, в этой работе не рассматриваются. Ниже будет приведен пример преобразования матрицы СЛАУ (14), составленной методом узловых потенциалов.
3.2. Демонстрационный пример
Прежде всего следует отметить, что матрица (14) симметричная и невырожденная. Коэффициенты матрицы удовлетворяют условию (22). Узловые потенциалы соответствуют падению напряжений на компонентах
U4 = UG^, U3 = UG,U2 =UG
Учитывая (28), СЛАУ (14) можно представить следующим образом:
G5 + G6 - G5 0 U 4 0
G5 G3 + G4 + G5 - G3 U3 = 0
0 - G3 G + G2 + G3 U2 GA
Проверка на корректность матрицы включает следующие операции.
Определение по формуле (23) параметра дискретной модели ait, входящего только
в диагональный элемент. Для первой строки матрицы это будет G6, для второй строки G4 и для третьей - (Gl + G2) .
Проверка строк матрицы на корректность выполняется в соответствии с формулой (24). В результате этой проверки оказывается, что вторая строка не удовлетворяет требованию корректности, так как (G4 = 1) ^ (G3 = 1015) . Параметр G3 входит также в третью строку матрицы, поэтому в соответствии с формулой (25) выбирается представление переменной U3 в виде
U3 = U2 + U23 , (30)
В результате преобразования элементов 3 -го столбца, в соответствии с формулой (26), получаем матрицу (29) следующего вида:
G5 + G6 - G5 - G5
G5 g3 + g4 + g5 g4+g5
а после преобразования третьей строки, в соответствии с формулой (27), матрица (31) будет иметь вид
(G5 + G6) - G5 - g5 U 4 0
G5 (G3 + G4 + G) (G4 + G5) U 23 = 0 . (32)
G5 (G4 + g5) (G + G2 + G4+g5) U2 G E
СЛАУ (32) удовлетворяет требованию корректности, поэтому корректировка считается завершенной. Переменные СЛАУ (32) соответствуют переменным СЛАУ (8), то есть в
ISSN 1028-9763. Математичні машини і системи, 2014, № 4
результате преобразования в дерево были выбраны те же компоненты, что и в методе корректной постановки задачи. Из сравнения СЛАУ (8) и (32) следует, что недиагональные элементы матрицы (32) второго столбца и второй строки отличаются по знаку от матрицы (8). Это есть результат того, что при преобразовании матрицы (14) было выбрано направление тока компоненты G3, противоположное направлению, выбранному при составлении СЛАУ (8). Выполнив замену переменной U23 на U23 = -U23 и поменяв во втором уравнении знаки элементов на противоположные, получим матрицу (8).
4. Заключение
Моделирование стало неотъемлемой частью интеллектуальной деятельности человечества, а достоверность результатов моделирования - основным критерием оценки результатов моделирования. Для обеспечения достоверности результатов требуются новые подходы к разработке методов и алгоритмов описания сложных объектов и их решения.
В отличие от существующего подхода к разработке методов решения некорректных задач, в данной работе предлагается некорректно поставленную задачу (плохо обусловленную) приводить к корректному виду. Показано, что не плохая обусловленность матрицы затрудняет получение достоверных результатов при решении СЛАУ, описывающих дискретные модели физических объектов, а некорректный выбор переменных СЛАУ на этапе составления уравнений, а метод узловых потенциалов и его аналоги, которые применяются для составления СЛАУ, описывающих дискретную модель, есть частный случай метода корректной постановки задачи. Предложена методика проверки на корректность составленной методом узловых потенциалов СЛАУ для случая, когда матрица СЛАУ есть невырожденная и симметричная. Рассмотрен алгоритм преобразования матрицы к корректному виду.
СПИСОК ЛИТЕРАТУРЫ
1. Калиткин Н.Н. Количественный критерий обусловленности систем линейных алгебраических уравнений / Н.Н. Калиткин, Л.Ф. Юхно, Л.В. Кузьмина // Математическое моделирование. - 2011. Т. 23, № 2. - С. 3 - 26.
2. Волобоев В.П. Об одном подходе к моделированию сложных систем / В.П. Волобоев, В.П. Клименко // Математичні машини і системи. - 2008. - № 4. - С. 111 - 122.
3. Волобоев В.П. Об одном подходе к моделированию энергосистем / В.П. Волобоев, В.П. Клименко // Математичні машини і системи. - 2009. - № 4. - С. 106 - 118.
4. Волобоев В.П. Механика стержневых систем и теория графов / В.П. Волобоев, В.П. Клименко // Математичні машини і системи. - 2012. - № 2. - С. 81 - 96.
5. Волобоев В.П. Метод конечных элементов и теория графов / В.П. Волобоев, В.П. Клименко // Математичні машини і системи. - 2013. - № 4. - С. 114 - 126.
6. Пухов Г.Е. Избранные вопросы теории математических машин / Пухов Г.Е. - Киев: Изд-во Академии наук УССР, 1964. - 264 с.
7. Сешу С. Линейные графы и электрические цепи / С. Сешу, М.Б. Рид. - М.: Высшая школа, 1971. - 448 с.
8. Зенкевич О. Конечные элементы и аппроксимация / О. Зенкевич, К. Морган. - М.: Мир, 1986. -318 с.
9. Воеводин В.В. Вычислительные основы линейной алгебры / Воеводин В.В. - М.: Наука, 1977. -304 с.
10. Теоретические основы электротехники: учебник для ВУЗов / К.С. Демирчян, Л.Р. Нейман, Н.В. Коровкин, В.Л. Чечурин. - . - Питер, 2003. - Т. 2. - 572 с.
Искомый вектор
Если , то система (1) называется плохо обусловленной. В этом случае погрешности коэффициентов матрицы и правых частей или погрешности округления при расчетах могут сильно исказить решение.
При решении многих задач правая часть системы (1) и коэффициенты матрицы А известны приближенно. При этом вместо точной системы (1) имеем некоторую другую систему
такую, что
Полагаем, что величины и d известны.
Так как вместо системы (1) имеем систему (2), то можем найти лишь приближенное решение системы (1). Метод построения приближенного решения системы (1) должен быть устойчивым к малым изменениям исходных данных.
Псевдорешением системы (1) называется вектор , минимизирующий невязку на всем пространстве .
Пусть х 1 –некоторый фиксированный вектор из , определяемый обычно постановкой задачи.
Нормальным относительно вектора х 1 решением системы (1) называется псевдорешение х 0 с минимальной нормой , то есть
где F-совокупность всех псевдорешений системы (1).
Причем
где ¾компоненты вектора х.
Для любой системы вида (1) нормальное решение существует и единственно. Задача нахождения нормального решения плохо обусловленной системы (1) является некорректно поставленной.
Для нахождения приближенного нормального решения системы (1) воспользуемся методом регуляризации.
Согласно указанному методу построим сглаживающий функционал вида
и найдем вектор , минимизирующий на этот функционал. Причем параметр регуляризации a однозначно определен из условия
где .
Вырожденные и плохо обусловленные системы могут быть неразличимы в рамках заданной точности. Но если имеется информация о разрешимости системы (1), то вместо условия (5) следует использовать следующее условие:
Компоненты вектора являются решениями системы линейных алгебраических уравнений, которая получается из условия минимума функционала (4)
и имеет вид
где Е-единичная матрица,
¾эрмитово сопряженная матрица.
На практике для выбора вектора нужны дополнительные соображения. Если их нет, то полагают =0.
Для =0 систему (7) запишем в виде
где
Найденный вектор будет являться приближенным нормальным решением системы (1).
Остановимся на выборе параметра a. Если a=0, то система (7) переходит в плохо обусловленную систему. Если a велико, то система (7) будет хорошо обусловлена, но регуляризованное решение не будет близким к искомому решению системы (1). Поэтому слишком большое или слишком малое a не пригодны.
Обычно на практике проводят расчеты с рядом значений параметра a. Например,
Для каждого значения a находят элемент , минимизирующий функционал (4). В качестве искомого значения параметра регуляризации берется такое число a, для которого с требуемой точностью выполняется равенство (5) или (6).
III. ЗАДАНИЕ
1. Построить систему линейных алгебраических уравнений, состоящую из трех уравнений с тремя неизвестными, с определителем, величина которого имеет порядок 10 -6 .
2. Построить вторую систему, аналогичную первой, но имеющую другие свободные члены, отличающиеся от свободных членов первой системы на величину 0,00006.
3. Решить построенные системы методом регуляризации (полагая =0 и d=10 -4) и каким-либо другим методом (например, методом Гаусса).
4. Сравнить полученные результаты и сделать выводы о применимости использованных методов.
IV. ОФОРМЛЕНИЕ ОТЧЕТА
В отчете должны быть представлены:
1. Название работы.
2. Постановка задачи.
3. Описание алгоритма (метода) решения.
4. Текст программы с описанием.
5. Результаты работы программы.
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Тихонов А.Н., Арсенин В.Я. Методы решения некорректных задач. - М.: Наука, 1979. 286 с.
2. Бахвалов Н.С., Жидков Н.П., Кобельков Г.М. Численные методы. - М.: БИНОМ. Лаборатория Знаний, 2007 636с.
Лабораторная работа № 23
Вернемся вновь к СЛАУ Aх=b с квадратной матрицей А размера MхN , которая, в отличие от рассмотренного выше "хорошего" случая (см. разд. 8.Г), требует специального подхода. Обратим внимание на два похожих типа СЛАУ:
- вырожденная система (с нулевым определителем |А|=0 );
- плохо обусловленная система (определитель А не равен нулю, но число обусловленности очень велико).
Несмотря на то, что эти типы систем уравнений существенно отличаются друг от друга (для первого решение отсутствует, а для второго существует и единственно), с практической точки зрения вычислителя, между ними много общего.
Вырожденные СЛАУ
Вырожденная система - это система, описываемая матрицей с нулевым определителем |A|=0 (сингулярной матрицей). Поскольку некоторые уравнения, входящие в такую систему, представляются линейной комбинацией других уравнений, то фактически сама система является недоопределенной. Несложно сообразить, что, в зависимости от конкретного вида вектора правой части ь, существует либо бесконечное множество решений, либо не существует ни одного. Первый из вариантов сводится к построению нормального псевдорешения (т. е. выбора из бесконечного множества решений такого, которое наиболее близко к определенному, например, нулевому, вектору). Данный случай был подробно рассмотрен в разд. 8.2.2 (см. листинги 8.11-8.13).
Рис. 8.7 . Графическое представление несовместной системы двух уравнений с сингулярной матрицей
Рассмотрим второй случай, когда СЛАУ Aх=b с сингулярной квадратной матрицей А не имеет ни одного решения. Пример такой задачи (для системы двух уравнений) проиллюстрирован рис. 8.7, в верхней части которого вводятся матрица А и вектор b , а также предпринимается (неудачная, поскольку матрица А сингулярная) попытка решить систему при помощи функции isolve . График, занимающий основную часть рисунка, показывает, что два уравнения, задающие систему, определяют на плоскости (x0,x1) две параллельные прямые. Прямые не пересекаются ни в одной точке координатной плоскости, и решения системы, соответственно, не существует.
Примечание
Во-первых, заметим, что СЛАУ, заданная несингулярной квадратной матрицей размера 2x2, определяет на плоскости пару пересекающихся прямых (см. рис. 8.9 ниже). Во-вторых, стоит сказать, что, если бы система была совместной, то геометрическим изображением уравнений были бы две совпадающие прямые, описывающие бесконечное число решений
.
Рис. 8.8
. График сечений функции невязки f (х) = |Ах-b|
Несложно догадаться, что в рассматриваемом сингулярном случае псевдорешений системы, минимизирующих невязку |Ax-b| , будет бесконечно много, и лежать они будут на третьей прямой, параллельной двум показанным на рис. 8.7 и расположенной в середине между ними. Это иллюстрирует рис. 8.8, на котором представлено несколько сечений функции f(x)= | Аx-b | , которые говорят о наличии семейства минимумов одинаковой глубины. Если попытаться использовать для их отыскания встроенную функцию Minimize , ее численный метод будет всегда отыскивать какую-либо одну точку упомянутой прямой (в зависимости от начальных условий). Поэтому для определения единственного решения следует выбрать из всего множества псевдорешений то, которое обладает наименьшей нормой. Можно пытаться оформить данную многомерную задачу минимизации в Mathcad при помощи комбинаций встроенных функций Minimize , однако более эффективным способом будет использование регуляризации (см. ниже) или ортогональных матричных разложений (см. разд. 8.3).
Две на первый взгляд похожие системы линейных уравнений могут обладать различной чувствительностью к погрешностям задания входных данных. Это свойство связано с понятием обусловленности системы уравнений .
Числом обусловленности линейного оператора A , действующего в нормированном пространстве а также числом обусловленности системы линейных уравнений Ax = у назовем величину
Таким образом, появляется связь числа обусловленности с выбором нормы.
Предположим, что матрица и правая часть системы заданы неточно. При этом погрешность матрицы составляет dA , а правой части - dу . Можно показать, что для погрешности dx имеет место следующая оценка ():
В частности, если dA = 0, то
При этом решение уравнения Ax = у не при всех у одинаково чувствительно к возмущению dу правой части.
Свойства числа обусловленности линейного оператора:
1.
причем максимум и минимум берутся для всех таких x , что Как следствие,
3
где и - соответственно минимальное и максимальное по модулю собственные значения матрицы A . Равенство достигается для самосопряженных матриц в случае использования евклидовой нормы в пространстве
4.
Матрицы с большим числом обусловленности (ориентировочно ) называются плохо обусловленными матрицами . При численном решении систем с плохо обусловленными матрицами возможно сильное накопление погрешностей, что следует из оценки для погрешности dx . Исследуем вопрос о погрешности решения, вызванной ошибками округления в ЭВМ при вычислении правой части. Пусть t - двоичная разрядность чисел в ЭВМ. Каждая компонента вектора правой части округляется с относительной погрешностью Следовательно,
Таким образом, погрешность решения, вызванная погрешностями округления, может быть недопустимо большой в случае плохо обусловленных систем.
Итак – принципиально остаются две проблемы –
1 .не обеспечивается обоснованная сходимость алгоритма к единственной (в случае модельного примера- истинной) структуре и
2 . Не разрешено противоречие о неадекватности моделей шаговой регрессии на новых точках, не участвовавших при оценке параметров модели. Возможно ли, если не обеспечить такую адекватность при других способах синтеза моделей, то хотя бы найти путь к решению такой задачи (возможно и адекватность определить другим способом)
Для АШР даже в случае применения для МНК оценки процедуры Грамма-Шмидта не разрешается вопрос о единственности модели – просто оценки параметров становятся наиболее точными и несмещенными
Т.о. гарантированное нахождение всего множества подходящих решений в реальных задачах (при - количество линейных входных аргументов и степени ПП p >3) получим только после полного
перебора всех подструктур полной структуры как в методе всех регрессий (у Дрейпера и Смита). Тогда мы найдем всете модели , в которых все аргументы входят с уровнем значимости не менее чем заданный. Со всеми выше описанными проблемами – а какая же из них, из этого множества та, которая действительно наша.
Можно еще добавить камень в огород АШР о неиспользуемой возможности вариации уровнем значимости для учета уровня шума в данных
Именно эту проблему предлагает решать МГУА с помощью введения понятия внешних критериев .
Необходимое примечание .
при все типы АШР МВИ, МГУА, другие целесообразные подходы, дают практически одинаково эффективные (или неэффективные) решения. Кривые критериев одинаково асимптотически стремятся к некоторому ненулевому уровню, при подходе к которому и определяется единственная модель.
Каждый из них это делает по своему, и определить адекватность метода по сходимости к нужной модели можно только построив соответствующий вашей задаче модельный пример.
Однако наиболее распространенный случай – это когда число точек невелико , тогда мы здесь решаем не переопределенную задачу (здесь точного решения нет и мы ищем среди плохих решений наилучшие), а близкую к определенной - вернее даже когда неизвестно, задача переопределена – определена или недоопределена. То есть включается сюда и совсем, казалось бы некорректная задача.
И наиболее эффективный подход к решению структурно-параметрического синтеза при данных условиях демонстрирует МГУА
Как видим нарушение уже первого условия порождает необходимомость разрешения проблемы множественности моделей не прибегая к процедуре полного перебора –надо предложить какой-то принцип, позволяющий найти путь к истинной или квазиистинной модели без полного перебора претендентов моделей.
Следующая проблема не менее реальна и еще более запутывает задачу поиска структуры. – проблема шума в данных – мы помним что при это нарушаются свойства проекционности аппарата МНК – нарушаются свойства оценок, но проблема в том что на зашумленных данных найти истинную структуру вообще может бытььпроблематично – если неизвестны х-ки шума и точки их приложения алгоритм будет тупо подстраиваться под шум.
Основная проблема – проблема необоснованности выбора структуры модели классическими АШГ многократно обостряется в связи с тем что порог используемый критерием Фишера в виде уровня значимости
на самом деле регулирует не только риск ошибки
– его выбор должен учитывать уровень шума и точки его приложения.
Ведь увеличение уровня шума например на выходе неизбено требует загрубить модель (не подстраивать ее под шум) а значит изменить уровень увеличить значимости с тем чтобы более жестко фильтровать апгументы в модель и ее упрощать.
Гораздо сложнее учесть шумы на входе тем более если они проходят нелинейное преобразование модели.
Однако методологически механизмов учета данных коррекций в агоритмах нет что делает выбор структур в условиях шума необоснованным.