Блочные численные операторные методы алгебраизации дифференциальных уравнений

Автор: Саух С.Е.
Институт проблем моделирования в энергетике им. Г.Е.Пухова НАНУ
http://scg.donntu.ru/2005/scg_2005.pdf

Работы академика Г.Е.Пухова в области создания новых операторных методов посвящены решению проблем, связанных с ограниченным применением классического операционного исчисления. Им предложено множество новых, более практичных и универсальных инструментов анализа систем и, самое главное, он привлек внимание своих учеников и последователей к фундаментальным проблемам, лежащим в основе разработки таких операционных методов анализа.

Так, в теории метода точек для представления непрерывных функций в области изображений Г.Е.Пухов предложил использовать изображающие векторы из конечного набора значений этих функций -«точек», а инфинитезимальным операциям над функциями-оригиналами были поставлены в соответствие матричные операторы интегрирования и дифференцирования, действующие на изображающие векторы. Поскольку матричные операторы строились на основе различных форм полиномиальных представлений непрерывных функций, то размерность и числовые значения элементов матриц определялись видом полиномиального представления и не зависели от значений функций в «точках». Такой подход позволил распространить на нелинейные системы методы исследования, характерные для традиционного операционного исчисления, что впервые было продемонстрировано в работах [1, 2] на примерах решения задач нелинейной электротехники.

В работах [1, 2] не рассматривались вопросы о полноте соотвествия основных соотношений теории метода точек основным соотношениям традиционного операционного исчисления. В частности, в теории метода точек в ряде случае не соблюдалась взаимная обратимость матричных операторов дифференцирования и интегрирования, что было показано в работе [3]. Матричные операторы оказывались взаимообратными только в случае использования для их построения тригонометрических функций. Однако такие операторы не сопрягались с начальными условиями и поэтому использовались ограниченно, например, для решения уравнений установившихся периодических режимов в нелинейных электрических цепях.

Проблема неполного соответствия в значительной степени была решена в процессе разработки теории локально-интегральных преобразований [4]. Принципиальное отличие этой теории от теории метода точек заключается в том, что при определении изображающих векторов использовались средние значения функций на некотором интервале, а не значения функций в «точках», что позволяло расширить множество преобразуемых функций и естественным образом включить в него класс обобщенных функций. В теории локально-интегральных преобразований было введено понятие определяющей матрицы линейной системы - аналога коэффициента передачи в операционном исчислении, построены взаимообратные матричные операторы интегрирования и дифференцирования, корректно учитывающие начальные условия, разработаны основы временно-частотного анализа систем.

В первоначальном виде теория локально-интегральных преобразований была основана исключительно на кусочно-постоянном представлении непрерывных функций. Поэтому для ее дальнейшего развития и создания вычислительных методов решения, например, задачи Коши для обыкновенных дифференциальных уравнений, потребовалось ввести специальные матричные операторы осреднения. Эти операторы строились независимо от значений усредняемой функции только на основе ее представления «скользящими», сдвинутыми один относительно другого, интерполяционными полиномами Ньютона. При этом на каждом последующем шаге осреднения функция была представлена новым полиномом, проходящим через узлы интерполяции, сдвинутые на один шаг вправо. Введение матричных операторов осреднения позволило построить новый класс параметрических многошаговых численных схем [5]. Однако сам способ построения матричных операторов осреднения породил проблему учета начальных условий, аналогичную проблеме «начала таблицы» многошаговых методов.

Анализируя связи между выбираемыми формами преобразования функций в область их изображений, с одной стороны, и возможностями достичь максимального упрощения трудно реализуемых операций над функциями-оригиналами при переходе в область изображений, с другой стороны, в работах [6, 7] Г.Е.Пухов

предложил использовать операционное исчисление на основе дифференциальных преобразований. Взяв за основу ряд Тейлора, как форму представления обратного преобразования изображения в функцию-оригинал, и построив формальные правила действий над коэффициентами ряда Тейлора для разных форм преобразуемых функций, он получил теоретическую основу для разработки новых аналитических, численно-аналитических и численных методов решения дифференциальных уравнений. Предложенное операционное исчисление получило широкое распространение и оказалось пригодным для анализа как линейных, так и нелинейных систем, в том числе и для некоторых одномерных по пространству систем с распределенными параметрами.

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

Методы, основанные на дифференциальных преобразованиях, используются только потому, что они позволяют эффективно решать интегро-дифференциальные уравнения. Благодаря гибкости при выборе функций, аппроксимирующих решение, построенные на основе дифференциальных преобразований вычислительные алгоритмы имеют высокую аппроксимационную точность. Вместе с тем следует отметить, что такие алгоритмы строятся на основе предварительного выполнения аналитических преобразований решаемых уравнений в форму дифференциальных изображений. Поскольку преобразованные уравнения могут коренным образом отличаться от формы исходного уравнения, то создаваемые алгоритмы не могут быть воплощены в универсальные программные коды-решатели уравнений, т.е. в коды, не требующие доопределения при изменении вида уравнений. Именно поэтому исследователи, использующие дифференциальное исчисление, разрабатывают специализированные программные средства для решения уравнений определенного вида.

Работы в области теории метода точек, локально-интегральных преобразований и преобразований Тейлора послужили базой для создания численных операторных методов на основе преобразований Ньютона [8-11]. Необходимость такой разработки вызвана нерешенностью упомянутых выше проблем объединения возмож¬ностей операторных и численных методов анализа систем, а именно:

• учета начальных условий в методе точек и «начала таблицы» в методах на основе локально-интегральных преобразований;
• взаимной обратимости операторов дифференцирования и интегрирования функций аппроксимированных полиномами;
• полиномиальной аппроксимации обобщенных функций, особенно -функции, имеющей большое значение в исследовании динамических систем;
• применимости операционного исчисления к многомерным функциям, определенным в областях сложной геометрической формы.

Использование численных операторных методов на основе преобразований Ньютона позволило комплексно решить указанные проблемы. Кроме того, эти методы в отличие от других, обладают необходимым для компьютерной реализации свойством алгоритмичности. В результате их применения для решения задачи Коши были получены новые блочные численные схемы решения систем интегро-дифференциальных уравнений с высокой аппроксимационной точностью и абсолютной устойчивостью, независящей от порядка используемых аппроксимирующих функций [8, 10, 11].

В случае целых одномерных функций преобразования Ньютона схожи с дифференциальными преобразованиями. Их принципиальное отличие состоит в использовании для представления функций ряда Ньютона, а не ряда Тейлора. Поскольку коэффициенты ряда Ньютона связаны не с производными разлагаемой в ряд функции, а с ее равноотстоящими значениями, входящими в определение разделенных разностей соответствующего порядка, изображение функции представляется вектором ее значений. Требование алгоритмичности преобразований удовлетворяет переход от ряда Ньютона к конечной сумме — интерполяционному полиному Ньютона. Естественно, при этом возникают погрешности аппроксимации преобразуемых функций, которые могут быть уменьшены до приемлемого уровня как повышением порядка интерполяционного полинома, так и дроблением интервала определения функции на локальные подинтервалы (не обязательно равной длины). Такая приближенная замена преобразуемой функции совокупностью интерполяционных полиномов Ньютона позволяет представлять функцию из ображающим вектором ее значений, подобно тому, как это делается в методе точек.

Операторы интегрирования и дифференцирования строятся на локальных подинтервалах. При построении операторов последовательно соблюдается правило соответствия порядков интерполяционных полиномов, представляющих соответствующие функции. В этом правиле учтено, что при дифференцировании полиномов порядка k получают полиномы порядка k-1. Следовательно, для своего однозначного представления изображающий вектор производной должен включать k ее значений, а изображающий вектор дифференцируемой функции должен состоять из k+1 значения. Соответственно, изображающий вектор интегральной функции должен быть образован из k+1 значений, а изображающий вектор интегрируемой функции — из k значений. Тогда связь между изображающими полиномами порядка k и k+1 может быть легко установлена путем введения матричных операторов дифференцирования и интегрирования, которые оказываются взаимообратными и позволяют корректно учитывать начальные условия с помощью векторов, сопряженных с матричными операторами [10, 11].

Введение сопряженных векторов позволило решить проблему полиномиальной аппроксимации обобщенной -функции. Поскольку соотношения, реализующие операции интегрирования и дифференцирования функций в области изображений Ньютона, совпадают по форме с аналогичными соотношениями в операционном исчислении Лапласа, выяснилось, что вектор, сопряженный с матричным оператором дифференцирования, представляет собой изображение -функции для случая ее кусочно-полиномиальной аппроксимации порядка k. Найденное аппроксимационное изображение -функции позволяло завершить разработку численных операторных методов. Так, были получены матричные аналоги коэффициентов передачи линейных систем и созданы численные алгоритмы работы с ними, что позволило предложить оригинальный численный метод приближенного обращения лапласовых изображений в область оригиналов.

Таким образом впервые удалось преодолеть трудности, ранее препятствовавшие объединению операторных и численных методов анализа систем в рамках единого подхода. Арсенал исследователя был существенно расширен не только за счет более совершенных блочных численных схем решения интегро-дифференциальных уравнений, но и за счет корректно построенных численных операторных методов анализа систем. Новые методы были исследованы на множестве тестовых задач с начальными и краевыми условиями, а также на примере решения одномерных уравнений газовой динамики, определенных на графе произвольной топологии [11]. Поскольку по своим точностным характеристикам и вычислительной устойчивости, а также благодаря возможности операторного представления моделей элементов топологического графа исследуемой системы предложенные методы превзошли известные численные методы, то они были положены в основу программного комплекса моделирования нестационарных неизотермических режимов в газопроводах.

Дальнейшая разработка численных операторных методов с использованием преобразований Ньютона была ориентирована на решение многомерных задач. Именно в этой области возможности традиционных численных методов таких, как метод сеток, проекционные методы и метод конечных элементов, оказывались значительно шире возможностей ранее разработанных операторных методов. Операторные методы использовались только для случаев простых прямоугольных форм областей определения многомерных функций.

Сложность разработки многомерного варианта преобразований Ньютона прежде всего состояла в том, что автору не удавалось найти в научной литературе готового варианта многомерного ряда Ньютона. С этой целью в работах [9, 10] были исследованы особенности построения многомерного ряда Тейлора и по аналогии с ним построен многомерный ряд Ньютона, коэффициенты которого, также как и коэффициенты одномерного ряда Ньютона, зависят от значений функции в точках ее многомерной области определения. С учетом требования алгоритмичности преобразований был осуществлен переход от ряда Ньютона к конечной его сумме — интерполяционному полиному Ньютона. В результате было обнаружено, что узлы интерполяции многомерного полинома равномерно заполняют симплексную подобласть. Это свойство многомерных интерполяционных полиномов Ньютона позволило создать многомерный вариант численных операторных методов не уступающих по своим возможностям, а в ряде случаев превосходящих традиционно применяемые численные методы.

В работах [9, 10] впервые описаны алгоритмы построения матричных операторов дифференцирования k -го порядка точности, действующие в нормированных симплексных областях произвольной мерности локальной системы координат. Эти специально разработанные алгоритмы базируются на множестве целых чисел, обладают свойством вложенности и обеспечивают получение операторов дифференцирования (N+ 1)-мерных функций из заданных операторов дифференцирования N -мерных функций. Преобразование

локальной системы координат в глобальную позволяет применять построенные операторы в симплексной области произвольных размеров при любых вариантах задания условий на ее границах.

Многомерные операторные методы на основе преобразований Ньютона имеют высокую аппроксимационную точность и просты в применении. Их использование значительно расширяет инструментальные возможности программных средств решения краевых задач для уравнений в частных производных.

1. Пухов Г.Е. Введение в теорию метода точек // Труды Таганрогского радиотехнического института. — 1954. — Вып. 1. — С. 47-77.
2. Пухов Г.Е., Борковский Б.А. К расчету электрической цепи, содержащей выпрямитель // Труды Таганрогского радиотехнического института. — 1954. — Вып. 1. —С. 78-82.
3. Пухов Г.Е. Методы анализа и синтеза квазианалоговых электронных цепей. — Киев: Наукова думка, 1967. — 568 с.
4. Береговенко Г.Я., Пухов Г.Е. Ступенчатые изображения и их применение. — Киев: Наукова думка, 1983. — 216 с.
5. Береговенко Г.Я., Саух С.Е. Параметрическая схема преобразования ступенчатых изображений и многошаговые методы численного решения дифференциальных уравнений // Электронное моделирование. — 1987. — N 5. — С. 7-11.
6. Пухов Г.Е. Преобразования Тейлора и их применение в электротехнике и электронике. — Киев: Наукова думка, 1978. — 260 с.
7. Пухов Г.Е. Дифференциальные преобразования функций и уравнений. — Киев: Наукова думка, 1980. — 420 с.
8. Саух С.Е. Локальные интегро-дифференциальные преобразования одномерных функций и уравнений. — Киев, 1988. — 53 с. (Препринт / АН УССР. Ин-т проблем моделирования в энергетике; N 142 ).
9. Саух С.Е. Многомерные локальные интегро-дифференциальные преобразования функций и дифференциальных уравнений. — Киев, 1990. — 53 с. (Препринт / АН УССР. Ин-т проблем моделирования в энергетике; 90-34).
10. Береговенко Г.Я., Пухов Г.Е., Саух С.Е. Численные операторные методы решения дифференциальных уравнений и анализа динамических систем. — Киев: Наукова думка, 1992. — 262 с.
11. Береговенко Г.Я., Гершгорин А.Е., Саух С.Е. Моделирование газотранспортных систем численными операторными методами (прямые и обратные задачи). — Киев, 1992. — 44 с. (Препринт / АН УССР. Ин-т проблем моделирования в энергетике; 92-53 ).
12. Саух С.Е. Численные операторные методы, основанные на рядах Ньютона // Электронное моделирование. — 1996. — N 4. — С. 63-69.