CC BY
77
ИНФОРМАТИКА И ВЫЧИСЛИТЕЛЬНАЯ
ТЕХНИКА
УДК 621.382.82.001
ДОСТОВЕРНОЕ И ТОЧНОЕ РЕШЕНИЕ СИСТЕМ АЛГЕБРАИЧЕСКИХ И ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙ НА ЯЗЫКЕ СИ
В.Б. Маничев1, Д.М. Жук1, В.В. Глазкова2
ХМГТУ им. Н.Э. Баумана, Москва, Российская Федерация e-mail: manichev@bmstu.ru; zhuk@bmstu.ru
2МГУ им. М.В. Ломоносова, Москва, Российская Федерация e-mail: glazkova@mlab.cs.msu.su
Приведены основные утверждения и научные результаты многолетних работ в области численного моделирования динамических систем, положенные в основу библиотеки математических программ SADEL (Sets of Algebraic and Differential Equations solvers Library) — Си-библиотеки, предназначенной для решения алгебраических и дифференциальных уравнений с максимально возможной компьютерной точностью. Показаны результаты решения тестовых и практических задач, подтверждающие эффективность разработанных программ библиотеки SADEL.
Ключевые слова: численное моделирование, динамические системы, обыкновенные дифференциальные уравнения, дифференциально-алгебраические уравнения, методы интегрирования, линейные алгебраические уравнения.
RELIABLE AND ACCURATE SOLVING OF SYSTEMS OF ALGEBRAIC AND DIFFERENTIAL EQUATIONS IN C LANGUAGE
V.B. Manichev1, D.M. Zhuk1, V.V. Glazkova2
1Bauman Moscow State Technical University, Moscow, Russian Federation e-mail: manichev@bmstu.ru; zhuk@bmstu.ru
2Lomonosov Moscow State University, Moscow, Russian Federation e-mail: glazkova@mlab.cs.msu.su
Basic statements and scientific results of many-year works in the field of numerical simulation of dynamic systems are given that underlie the mathematical program library SADEL (Sets of Algebraic and Differential Equations solvers Library). SADEL is the C language library destined for solving algebraic and differential equations with maximum possible computer precision. Results of solving the test and practical problems are given that confirm the effectiveness ofthe SADEL programs.
Keywords: simulation, dynamic systems, ordinary differential equations, differential algebraic equations, integration methods, linear algebraic equations.
Огромное значение при проектировании сложной и опасной техники (аэрокосмической техники, атомных станций и т.п.) имеют достоверность, надежность и точность компьютерных вычислений. Как показывает анализ некоторых из происшедших за последнее время аварий и катастроф, одной
из причин стали существенные погрешности компьютерных вычислений, результаты которых попали в конечную конструкторскую документацию. Главным недостатком вычислительных методов, алгоритмов и программ, используемых при проектировании и управлении и реализованных в лучших зарубежных программно-методических комплексах (ПМК) моделирования динамических систем, является возможная выдача ошибочных результатов компьютерных вычислений без предупреждения пользователей об их недостоверности. Подавляющее большинство универсальных и, особенно, специализированных программ численного моделирования разнообразных динамических систем и объектов разработано и разрабатывается с использованием языка Си (Си++). При этом пользователи ПМК моделирования динамических систем часто сами пишут программы решения систем обыкновенных дифференциальных и дифференциально-алгебраических уравнений (ОДУ-ДАУ), а также линейных алгебраических уравнений (ЛАУ), считая это простой задачей, а затем разрабатывают прекрасные интерфейсы и визуализацию неверных результатов решения сложных вычислительных задач. Библиотека SADEL (Sets of Algebraic and Differential Equations solvers Library) разработана как библиотека программ-решателей для жестких систем ОДУ-ДАУ и плохо обусловленных систем ЛАУ и предназначена для решения указанных проблем при численном моделировании динамических систем [1, 2]. В настоящей работе используется следующая терминология.
Математическая модель (mathematical model) — система ОДУ-ДАУ, отображающая функционирование моделируемой динамической системы с требуемой достоверностью и точностью.
Математическое моделирование (mathematical modeling) — процесс получения математических моделей динамических систем на основе фундаментальных законов физики, химии и др.
Компьютерная модель — алгоритмическая и программная реализация математической модели.
Численное (компьютерное) моделирование (далее для сокращения просто моделирование (simulation)) — процесс получения результатов математического моделирования динамических систем с помощью компьютерных моделей.
Жесткие системы ОДУ-ДАУ (stiff ODE-DAE systems (stiff set of ODE-DAE)) — системы ОДУ-ДАУ со степенью жесткости более 106.
Плохо обусловленные системы ЛАУ (badly conditioned (ill-conditioned) LAE systems (ill-conditioned set of LAE)) — системы ЛАУ с числом обусловленности более 106, к решению которых часто сводится решение жестких и сверхжестких систем ОДУ-ДАУ неявными методами интегрирования.
Достоверность численного моделирования (veracity of simulation) — отклонение результатов численного моделирования динамических систем от заведомо точных результатов не более чем на 20 % (точные результаты получают аналитически или экспериментально, как правило, оценивая относительную погрешность численного моделирования в процентах).
Точность численного моделирования (accuracy of simulation) — фактическая степень отклонения результатов достоверного численного моделирования от заведомо точных результатов.
Математическая точность (mathematical accuracy) решения систем ОДУ-ДАУ - заданная относительная погрешность нормы векторов решения систем ОДУ-ДАУ (параметр eps) по сравнению с заведомо точным аналитическим решением в относительных единицах (например, в ПМК MATLAB по умолчанию eps = 0,001).
Компьютерная точность (precision) решения систем ОДУ-ДАУ и ЛАУ — число полученных верных значащих цифр в каждом элементе векторов решения систем ОДУ-ДАУ и ЛАУ по сравнению с заведомо точным аналитическим (идеальным для систем ЛАУ) решением.
Решение систем ЛАУ с повышенной точностью вычислений (LAE systems solution with extra precision calculations) — применение методов получения решений с повышенной точностью вычислений (extra precision calculations), реализованных в соответствующих ПМК, а также в библиотеках математических программ и описанных в документации на соответствующие ПМК и библиотеки. Возможно также использование переменных типа Decimal на языке C# или класса java.math.BigDecimal на языке Java.
Далее приведены методы для получения решений с повышенной точностью вычислений, реализованные в наиболее известных ПМК: Maple (http://maplesoft.com — метод Software Floating Point (SFP метод)), MATLAB (http://mathworks.com — метод Variable Precision Arithmetic (VPA метод)), Mathematica (http://wolfram.com — метод Arbitrary Precision Arithmetic (APA метод));
библиотеки математических программ для получения решений с повышенной точностью вычислений: Intel для десятичных вычислений с повышенной точностью (стандарт IEEE 754-2008 Decimal Floating-Point for Intel§Architecture Processors, IEEE Standard for Float-ing-Point Arithmetic 7542008, IEEE, 2008:
http://ieeexplore.ieee.org/xpl/loginjsp?tp=&arnumber=5223319&url [3]), а также Си (Си++)-библиотеки для вычислений с повышенной точностью: http://en.wikipedia.org/wiki/Arbitrary-precision_arithmetic http://en.wikipedia.org/wiki/GNU_Multi-Precision_Library http://en.wikipedia.org/wiki/MPFR http://en.wikipedia.org/wiki/Class_Library_for_Numbers http://gmplib.org/
http://www. boost.org/doc/libs/1_53_0/libs/rational/rational.html http://speleotrove.com/decimal/dnusers.html#example7 http://www. ttmath.org/.
Основные утверждения и научные результаты, положенные в основу разработки библиотеки SADEL. Процесс численного моделирования разнообразных систем и объектов, в частности, динамических систем, в
общем случае (без использования ПМК визуального математического моделирования типа ANSYS, NX-CAE (NASTRAN), SIMULIA (ABAQUS), MATLAB-SIMULINK, Maple-MapleSim, Mathematica-SystemModeler и др.), можно представить как последовательность четырех основных этапов.
1. Получение на основе фундаментальных физических или химических законов исходной математической модели для объекта математического моделирования в форме систем уравнений, математических выражений, зависимостей и т.п.
2. Для выбранного алгоритмического языка определение численных методов и алгоритмов решения или выбор из соответствующей библиотеки готовой программы-решателя уравнений математической модели (в этом случае третий этап не нужен).
3. Разработка на выбранном алгоритмическом языке и отладка расчетной программы-решателя уравнений математической модели.
4. Ввод исходных данных для разработанной или выбранной программы-решателя уравнений математической модели, получение и обработка результатов численного моделирования.
Соответственно можно выделить четыре основных источника ошибок в процессе численного моделирования.
1. Погрешность исходной математической модели. Это очень важный источник ошибок, так как он является неустранимым. Например, второй закон Ньютона в дифференциальной форме не всегда дает корректные результаты математического моделирования, даже при отсутствии других источников ошибок математического моделирования. Анализ этих ошибок выходит за рамки настоящей статьи, поскольку относится к разным предметным областям математического моделирования (multi-physics, multi-discipline simulation).
2. Ошибки, связанные с неверным выбором численных методов и алгоритмов. Это устранимый источник ошибок, так как пользователям ПМК численного моделирования динамических систем не следует самому программировать решение вычислительных задач, а следует использовать готовые программы-решатели уравнений математических моделей из соответствующих ПМК или известных библиотек математических программ. Это устранимый источник ошибок математического моделирования.
3. Ошибки, связанные с ошибками программирования вычислений из-за ограниченной разрядной сетки компьютера. Ошибки программирования являются устранимыми и их можно исправить, поэтому в настоящей работе учитываются только неустранимые ошибки вычислений, связанные, в частности, с ошибками округления чисел при вычислениях с числами с плавающей точкой на компьютере, которые обязательно надо свести к минимуму или сообщать пользователю соответствующей программы-решателя о недопустимом значении этих ошибок.
4. Ошибки, связанные с вводом неверных исходных данных. Это устранимый источник ошибок математического моделирования.
Результатом многолетних исследований проблем численного моделирования динамических систем являются следующие утверждения.
Математические модели динамических технических систем необходимо получать на основе фундаментальных физических законов только в форме систем ДАУ, не разрешенных относительно производных, и решать эти системы без каких-либо преобразований и без получения "правой" части для производных в явном, аналитическом виде. При математическом моделировании динамических процессов в реальных технических системах и объектах состояние этих систем и объектов обычно рассматривают в пространстве дифференциальных переменных или в пространстве "переменных состояния" с получением в явном аналитическом виде "правой" части для производных этих переменных по времени, но этого недостаточно для достоверного и точного математического моделирования реальных технических систем и объектов, так как любой дифференциальной переменной состояния в технических системах обычно соответствует алгебраическая переменная. Например, в электрике и электронике дифференциальным переменным состояния — напряжениям на емкостях и токам через индуктивности — соответствуют алгебраические переменные — токи через емкости и напряжения на индуктив-ностях, в механике дифференциальным переменным состояния — скоростям тел определенной массы — соответствуют алгебраические переменные — силы инерции и другие параметры, поэтому состояние таких динамических систем и объектов следует рассматривать в пространстве дифференциально-алгебраических переменных. Основным физическим свойством алгебраических переменных является возможность идеального скачка значения этих переменных в бесконечно малый отрезок времени, в то время как для дифференциальных переменных состояния такие скачки невозможны, поэтому приведение систем ДАУ к системам ОДУ в нормальной форме Коши с преобразованием алгебраических переменных в дифференциальные переменные состояния физически неверно и полученная таким образом математическая модель динамических процессов в реальных технических системах и объектах будет физически недостоверна.
Эффективность численных методов решения систем ОДУ-ДАУ, НАУ и ЛАУ включает в себя два аспекта:
1) достоверность и точность получаемых решений данных уравнений;
2) затраты машинного времени и оперативной памяти на получение решений.
Быстрый рост производительности вычислительной техники (особенно производительности процессоров компаний Intel и NVIDIA и современных суперкомпьютеров) в последнее время выдвигает на первый план проблему получения достоверности и точности получаемых решений. При этом должна быть обеспечена гарантия требуемых достоверности и точности выдаваемых программой-решателем результатов компьютерных расчетов, математического и компьютерного моделирования, особенно инженерам-проектировщикам,
поскольку ошибки в принятии проектных решений для промышленных изделий и объектов, основанных на этих результатах, стоят очень дорого. При этом надо учесть следующее:
— большинство инженеров-проектировщиков не являются специалистами в численных методах и программах для решения различных математических уравнений, поэтому требуемая достоверность и точность должна быть обеспечена для параметров программ-решателей этих уравнений, рекомендуемых для этих программ по умолчанию;
— при базовом, начальном математическом и компьютерном моделировании динамических технических систем не требуется высокая математическая точность выдаваемых пользователю результатов, так как параметры математических моделей этих систем (и соответствующие коэффициенты систем ОДУ-ДАУ и ЛАУ) получены, как правило, экспериментально с невысокой математической точностью. Поэтому заданная математическая точность конечных результатов численного моделирования по умолчанию может быть невысокой, но она должна быть гарантировано обеспечена.
Программа-решатель систем ДАУ в библиотеке 8АББЬ рассматривает системы ДАУ в еще более расширенном координатном пространстве переменных, в котором к дифференциально-алгебраическим переменным добавлены производные дифференциальных переменных [4]:
О(РХ, X, У, г) = 0,
где X — вектор координатного базиса дифференциальных переменных состояния размерностью т, У — вектор алгебраических переменных полного координатного базиса дифференциально-алгебраических переменных размерностью к, РХ = йХ/йг — вектор производных дифференциальных переменных по времени для расширенного координатного пространства переменных размерностью т,Ь — независимая переменная (обычно время); О — вектор-функция математической модели динамической системы размерностью т + к. Заданы начальные условия Х0 = X(0) и отрезок интегрирования
г = [т о,тк ].
Новые алгоритмы для решения приведенной системы основаны на совместном решении систем НАУ Ог(Х,У, РХ,гг) = 0 и систем из т алгебраических уравнений численного интегрирования Н^ (X, РХ, Н) = 0,
сформированных на соответствующих стадиях методов интегрирования, относительно дифференциально-алгебраических переменных Х,У и производных дифференциальных переменных РХ (DABC методы Маничева-Глазковой) [2]:
Щ(РХъ,Хъ,Хп-1,РХп-1,Нп) =
в в
= Нп йгзРХз - аг3Х] - ЬгХп-1 - НпСгРХп-1 = 0;
3=1 з=1
Хп = Х3; РХп = РХ3, г = 1,... ,з, = з — ; Нп = (гп — гп-\) — п-й шаг интегрирования; йгз, (щ, Ьг,сг — параметры метода.
Параметры алгебраических уравнений численного интегрирования И (Х,РХ,Н) = 0 были получены на основе анализа пяти базовых классов задач для фундаментальных решений линейных неоднородных систем ОДУ. Было показано, что для достоверного решения всех пяти классов фундаментальных решений при больших значениях шага интегрирования (т.е. при невысокой математической точности интегрирования) численные методы решения систем ОДУ-ДАУ должны быть ЛЬ-устойчивыми, т.е. абсолютно (А-) устойчивыми строго в левой (Ь-ЬеА) полуплоскости комплексной плоскости устойчивости методов численного решения систем ОДУ-ДАУ [5]. При невысокой математической точности интегрирования шаг интегрирования может стать таким большим, что не ЛЬ-устойчивые методы интегрирования могут потерять устойчивость для некоторых фундаментальных решений. Область абсолютной устойчивости ЛЬ-устойчивых методов интегрирования показана на рис. 1.
Следует отметить, что в задачах химической кинетики и небесной механики алгебраические переменные обычно не нужны, но значения производных дифференциальных переменных требуется вычислять с такой же точностью, как и значения дифференциальных переменных.
Численное решение систем ОДУ-ДАУ во многих численных методах их решения сводится к многократному численному решению соответствующих систем нелинейных алгебраических уравнений (НАУ) и ЛАУ [5]. Высокий порядок таких систем и необходимость их многократного решения в процессе математического моделирования предъявляют жесткие требования к эффективности применяемых численных методов решения систем ОДУ-ДАУ и особенно ЛАУ.
Для решения жестких систем ОДУ-ДАУ, которым часто соответствует плохо обусловленная матрица Якоби вдоль траектории решения, рекомендуется использовать неявные методы интегрирования, которые сводятся к решению систем НАУ на каждом шаге интегрирования. Для получения решений систем НАУ с гарантированной компьютерной точностью эти системы следует решать только такими итерационными методами, которые в "точке" решения системы НАУ сводятся к решению системы ЛАУ с матрицей Якоби в качестве матрицы коэффициентов такой системы ЛАУ. Если не получать решения систем НАУ с гарантированной компьютерной точностью, то траектории решения систем ОДУ-ДАУ могут быть ошибочными. В библиотеке 8ЛББЬ используется метод Ньютона-Рафсона и решение систем НАУ и ЛАУ с гарантированной удвоенной компьютерной точностью языка Си (Си++) на каждом шаге интегрирования, т.е. с гарантированной компьютерной точностью (15 верных значащих цифр для всех элементов вектора
Рис. 1. Область абсолютной устойчивости (заштрихована) AL-устойчивых методов интегрирования
решений систем ЛАУ при вычислениях с удвоенной точностью на языке Си). Аналитически данная проблема не имеет решения, например, теоретически умножение квадратной матрицы с постоянными коэффициентами на обратную равно единичной матрице, но для плохо обусловленных матриц при вычислении на компьютере единичную матрицу мы не получим. Только многолетние вычислительные эксперименты позволили решить эту проблему в программах-решателях систем ЛАУ библиотеки SADEL на основе перечисленных методов вычислений с повышенной точностью.
При разработке математических программ на языке Си необходимо также учитывать особенности программирования арифметических операций с вещественными числами с плавающей точкой, так как именно такие операции используют упомянутые ПМК и библиотеки математических программ при решении систем ОДУ-ДАУ и ЛАУ. Институт инженеров по электротехнике и электронике (Institute of Electrical and Electronic Engineers, IEEE) разработал международные стандарты, которые описывают представление вещественных чисел с плавающей точкой:
— стандарт ANSI/IEEE 754:1985 определяет требования к реализации двоичной плавающей арифметики;
— стандарт ANSI/IEEE 754:2008 обобщает стандарт ANSI/IEEE 754:1985, допуская дополнительно, кроме двоичного, десятичное основание представлений мантиссы и экспоненты и произвольную длину машинного слова.
Данные стандарты, кроме форматов представления, описывают также основные арифметические действия, операции вычисления остатка от деления, квадратного корня, преобразования, округления и т.п. В современных платформах, таких как Intel, и в большинстве RISC-систем, аппаратно реализована плавающая арифметика, соответствующая этим стандартам. Другие стандарты на практике не используются.
Стандарты IEEE754 определяют следующие основные форматы представления вещественных чисел с плавающей точкой.
1. С одинарной точностью (соответствует типам SINGLE PRECISION (32 двоичных разряда) или REAL*4 в языке Фортран и float в языке Си).
2. С удвоенной точностью (соответствует типам DOUBLE PRECISION (64 двоичных разряда) или REAL*8 в языке Фортран и double в языке Си), этот тип представления принят по умолчанию во всех упомянутых ПМК и стандартных библиотеках математических программ.
3. С расширенной точностью (условно говоря, соответствует типам REAL*10 (80 двоичных разрядов) в языке Фортран и long double в языке Си).
Формат с учетверенной точностью QUAD PRECISION (128 двоичных разряда) или REAL*16 в языке Фортран стандарт IEEE754 отдельно не определяет, поэтому в языке Си следует создавать специальные форматы для данного представления вещественных чисел.
При реализации вычислительных методов решения жестких систем ОДУ-ДАУ и плохо обусловленных систем ЛАУ не следует использовать вычисления с одинарной точностью, причем при решении плохо обусловленных систем ЛАУ следует использовать методы вычислений с повышенной точностью, если вычисления с общепринятой удвоенной точностью не обеспечат компьютерной точности результатов. При программной реализации формул неявного интегрирования жестких систем ОДУ-ДАУ не должно быть также арифметических операций деления на шаг интегрирования, чтобы не было ограничений из-за "машинного нуля" на значение минимального шага интегрирования [6].
Результаты решения тестовых и практических задач. На основе приведенных утверждений и научных результатов была разработана библиотека SADEL — библиотека программ-решателей на языке Си для численного решения жестких систем ОДУ-ДАУ и плохо обусловленных систем ЛАУ с гарантированной достоверностью и компьютерной точностью. В результате при невысоких требованиях к математической точности конечных результатов численного моделирования получаются только корректные, достоверные и точные результаты решения данных систем по сравнению с соответствующими программами-решателями из зарубежных библиотек и пакетов математических программ: LAPACK (Linpack), MAGMA, NAG, IMSL, Intel MKL, MATLAB, Mathcad, Maple, Mathematica), что показано в работе [2]. Например, с помощью данных программ удалось решить с гарантированной компьютерной точностью ряд тестовых плохо обусловленных систем ЛАУ, которые не решаются с обычной удвоенной точностью языка Си другими известными программами-решателями систем ЛАУ из указанных библиотек и ПМК. В качестве примера рассмотрим тестовую задачу с матрицей Гильберта (Hilbert matrix) 10-го порядка в качестве матрицы коэффициентов системы ЛАУ. Все элементы матрицы коэффициентов и вектора правых частей, как суммы соответствующих строк матрицы коэффициентов, были вычислены с удвоенной точностью языка Си. Абсолютно точное решение такой задачи — единичный вектор. При решении с удвоенной точностью этой тестовой задачи с матрицей 10-го порядка с помощью программы-решателя систем ЛАУ из библиотеки MAGMA (LAPACK-Linpack) был получен следующий результат: 1.000000000464697 0.999999960244827 1.000000840984605 0.999992393461859 1.000036137385782 0.999900979283788 1.000162025520190 0.999843779478675 1.000081852296234 0.999982030574678
Курсивом выделены неверные значащие цифры в решении. Соответствующие программы-решатели из библиотеки SADEL получили единичный вектор для этой тестовой задачи с точностью в 15 верных значащих цифр для всех элементов вектора решений данной системы ЛАУ.
Сравнение программы-решателя систем ДАУ в библиотеке SADEL с зарубежными программами-решателями систем ОДУ-ДАУ было проведено
для параметров методов интегрирования программ-решателей систем ОДУ-ДАУ, рекомендуемых для этих программ-решателей по умолчанию. Все программы-решатели систем ОДУ-ДАУ сравнивались при заданной невысокой относительной точности интегрирования eps = 0,001 (в ПМК MATLAB это точность, задаваемая по умолчанию). В таблице приведены результаты сравнения только наиболее трудных для современных программ-решателей систем ОДУ-ДАУ жестких систем ОДУ с многопериодным решением. Если программа-решатель систем ОДУ-ДАУ не решит достоверно и точно хотя бы одну из приведенных в таблице тестовых задач, то эту программу-решатель не следует использовать для моделирования динамических систем, математическая модель которых представляет собой жесткую систему ОДУ-ДАУ с многопериодным решением [2].
В таблице сравнивались только методы интегрирования программ-решателей систем ОДУ-ДАУ, рекомендуемые в соответствующих ПМК и библиотеках математических программ для решения жестких систем ОДУ-ДАУ. Знак минус означает невозможность получения решения или (в большинстве случаев) качественно неверный результат без всякого предупреждения пользователя о возможных ошибках. Например, на рис. 2 приведен результат качественно неверного решения теста 4 с помощью программы-решателя систем ОДУ-ДАУ неявным методом BDF из библиотеки C-Library NAG (а), а также правильное решение этого теста с помощью методов М2 и М3 из библиотеки SADEL (б). Предупреждений о возможном недостоверном решении в программе-решателе библиотеки C-Library NAG не было.
Далее приведен пример решения практической задачи химической кинетики [7]. Данная задача является жесткой за счет большого разброса порядков в константах скоростей реакции и собственных значений матрицы Якоби, коэффициент жесткости более 106. Для решения этой задачи сначала был выбран пакет математических программ Mathematica 8.0. Для проверки достоверности и точности полученного решения была выбрана Си-библиотека программ-решателей жестких систем ОДУ-ДАУ и плохо обусловленных систем ЛАУ SADEL.
Как следует из рис.3,а, метод, выбранный в пакете Mathematica 8.0 по умолчанию выдал колебательный график решения. Неявный метод М3, реализованный в библиотеке SADEL, выдал правильный график решения — рис. 3, б. Следует отметить, что после настройки параметров методов интегрирования, реализованных в пакете Mathematica 8.0, был также получен правильный график решения, однако большинство пользователей математических программ не являются профессиональными специалистами в численных методах и программах для решения систем ОДУ-ДАУ, поэтому достоверность и требуемая точность решения систем ОДУ-ДАУ должна быть обеспечена для параметров программ-решателей этих уравнений, рекомендуемых для этих программ-решателей по умолчанию.
Заключение. Библиотека Си программ SADEL для решения жестких систем ОДУ-ДАУ и плохо обусловленных систем ЛАУ является практическим
Сравнение программ-решателей систем ОДУ-ДАУ для жестких систем ОДУ с многопериодным решением
Тесты Шрограммы-решатели систем ОДУ-ДАУ 8АБЕЬ-Си библиотека 2010. Методы М2,МЗ МАТЬАВ 2007. Метод Оёе15в Мар1е 2007. Метод ЯозепЬгоск МаИтетаИса 8 2011. Метод ВОБ МАО С-ЫЬгагу 2012. Метод е 1'>Г)К ШБЬ С-ЫЬгагу 2012. Метод ВОБ
ТЕСТ 1. Уравнения Ван дер Поля ми = 10е + + + +
ТЕСТ 2. Уравнения Ван дер Поля ми = ю9 +
ТЕСТ 3. Высокодобротный фильт Ы = 1, кг = 1, ки = 0,01 + + +
ТЕСТ 4. Высокодобротный фильтр ы = ю-104, ки = 1, Ы = 1 + +
ТЕСТ 5. Локально неустойчивая система ОДУ ми = 10е + +
ТЕСТ 6. Моделирование свечения лазера + + + +
Рис.2. Решение задачи ТЕСТ 4 при TOL = 0,001:
а — метод BDF из библиотеки C-Library NAG, б — метод М3 из Си библиотеки SADEL
результатом многолетних научно-исследовательских работ авторов в качестве математического ядра платформы математического и компьютерного моделирования разнородных (multi-physics, multi-discipline) динамических технических систем и объектов (программный комплекс ПА10 (SADEL-PA10)) [2], превосходящей подобные зарубежные программные продукты MATLAB-SIMULINK, Maple-MapleSim, C-Library NAG и C-Library IMSL, Intel MKL в части решения жестких систем ОДУ-ДАУ и плохо обусловленных систем ЛАУ со стандартной удвоенной точностью вычислений языка Си.
Разработка программ-решателей систем ОДУ-ДАУ и ЛАУ на языке Си для библиотеки SADEL позволяет также сделать следующие выводы.
1. Корректные математические модели динамических процессов в реальных технических системах и объектах необходимо получать в расширенном дифференциально-алгебраическом координатном базисе переменных на основе фундаментальных физических законов в форме неоднородных систем ДАУ, не разрешенных относительно производных, и решать эти системы без каких-либо эквивалентных математических преобразований и без получения "правой" части для производных в явном, аналитическом виде.
2. Для достоверного и точного решения систем ДАУ при невысоких требованиях к математической точности результатов следует использовать толь-
Рис.3. График зависимости дифференциальной переменной х8 системы ОДУ от времени при TOL = 10-6: а - МаЛешаНеа 8.0; б - 8ЛББЬ
ко ЛЬ-устойчивые методы интегрирования, при использовании других методов интегрирования необходимо гарантировать устойчивость численного решения при увеличении шага интегрирования в случае невысоких значений параметра ерв.
Новые научные результаты предполагается получить в направлениях разработки новых высокоэффективных неявных БЛБС методов интегрирования с полностью заполненными матрицами Б и Л в новых системах алгебраических уравнений интегрирования, разработки новых методов и алгоритмов автоматического выбора шага интегрирования для ЛЬ-устойчивых неявных методов интегрирования, а также эффективных алгоритмов и программ решения систем ЛАУ высокой и сверхвысокой размерности с повышенной точностью вычислений.
Работа выполнена при поддержке РФФИ, гранты № 12-07-00324/12 и № 12-07-31029.
ЛИТЕРАТУРА
1. Андронов А.В., ЖукД.М., Кожевников Д.Ю., Маничев В.Б. Библиотека математических программ-решателей на языке Си: SADEL // http://www.pa10.ru (дата обращения 20.10.2013).
2. Жук Д.М., Маничев В.Б., Сахаров М.К. SADEL — библиотека "сверхточных" решателей для программного комплекса ПА10 (SADEL-PA10) // Сб. науч. трудов МЭС-2012. М.: ИППМ РАН, 2012. С. 147-153.
3. IEEE Standard for Floating-Point Arithmetic 754-2008, IEEE, 2008 // http://iee-explore.ieee.org/xpl/login.jsp?tp=&arnumber=5223319&url=http%3A%2F%2Fiee-explore.ieee.org%2Fxpls%2Fabs_all.jsp%3Farnumber%3D5223319 (дата обращения 21.10.2013).
4. Норенков И.П., Трудоношин В.А., Федорук В.Г. Метод формирования математических моделей для адаптируемых программных комплексов анализа радиоэлектронных схем // Радиотехника, 1986. № 9. С. 67-72.
5. Жук Д.М., Маничев В.Б., Ильницкий А.О. Методы и алгоритмы решения дифференциально-алгебраических уравнений для моделирования систем и объектов во временной области // Информационные технологии. 2010. Ч. 1. № 7. С. 16-24. Ч. 2. № 8. С. 23-26.
6. Маничев В.Б. Новые алгоритмы для программ анализа динамики технических систем // Вестник МГТУ им. Н.Э. Баумана. Сер. Приборостроение. 1996. № 1. С. 48 -56.
7. Valerii N. Snytnikov, Tatyana I. Mishchenko, Vladimir N. Snytnikov, Sergei E. Malykhin, Vasilii I. Avdeev, Valentin N. Parmon. Autocatalytic gas-phase dehydrogenation of ethane // Research on Chemical Intermediates, March 2012. Vol. 38. Issue 3-5. Р. 1133-1147.
REFERENCES
[1] Andronov A.V., Zhuk D.M., Kozhevnikov D.Yu., Manichev V.B. Biblioteka
matematicheskikh programm-reshateley na yazyke Si: SADEL [Library of mathematical problem solvers in C: SADEL] Available at: http://www.pa10.ru (accessed 14.01.14) (in Russ.).
[2] Zhuk D.M., Manichev V.B., Sakharov M.K. SADEL — library "extremely-precise"
problem solvers for complex software PA10 (SADEL-PA10). Sb. nauch. trudov VI Vserossiyskoy nauchno-tekhn. konf. "Problemy razrabotki perspektivnykh mikro- i nanoelektronnykh sistem MES-2012 [Call. Papers of All-Russ. Sc.-Techn. Со^. "Development problems of advanced micro-and nanoelectronic systems"], Moscow, IPPM RAN Publ., 2012. pp. 147-153 (in Russ.).
[3] IEEE Standard for Floating-Point Arithmetic 754-2008, IEEE, 2008. Available at:
http://ieeexplore.ieee.org/xpl/login.jsp?tp=&arnumber=5223319&url=http%3A% 2F%2Fieeexplore.ieee.org%2Fxpls%2Fabs_all.jsp%3Farnumber%3D5223319 (accessed 28.12.2013).
[4] Norenkov I.P., Trudonoshin V.A., Fedoruk V.G. Method of mathematical models
formation for adaptable software systems of electronic circuits analysis. Radiotekhnika [Radioengineering], 1986, no. 9, pp. 67-72 (in Russ.).
[5] Zhuk D.M., Manichev V.B., Il'nitskiy A.O. Methods and algorithms for solving
differential-algebraic equations for simulation systems and objects in the time domain. Informatsionnye tekhnologii [Information technologies], 2010, no. 7, pp. 16-24 (ch. 1); no. 8, pp. 23-26 (ch. 2) (in Russ.).
[6] Manichev V.B. New algorithms for program analysis of the dynamics of technical
systems. Vestn. Mosk. Gos. Tekh. Univ. im. N.E. Baumana, Priborostr. [Herald of the Bauman Moscow State Tech. Univ., Instrum. Eng.], 1996, no. 1, pp. 48-56 (in Russ.).
[7] Snytnikov V.N., Mishchenko T.I., Snytnikov V.N., Malykhin S.E., Avdeev V.I., Parmon V.N. Autocatalytic gas-phase dehydrogenation of ethane. Research on Chem. Int., 2012, vol. 38, iss. 3-5, pp. 1133-1147.
Статья поступила в редакцию 26.11.2012
Владимир Борисович Маничев — канд. техн. наук, доцент кафедры "Системы автоматизированного проектирования" МГТУ им. Н.Э. Баумана. Автор более 80 научных работ в области математического обеспечения систем автоматизированного проектирования.
МГТУ им. Н.Э. Баумана, Российская Федерация,105005, Москва, 2-я Бауманская ул., д. 5.
V.B. Manichev — Cand. Sci. (Eng.), assoc. professor of "Systems of Computer-Aided Design" department of the Bauman Moscow State Technical University. Author of more than 80 publications in the field of mathematical support of computer-aided design systems. Bauman Moscow State Technical University, Vtoraya Baumanskaya ul. 5, Moscow, 105005 Russian Federation.
Дмитрий Михайлович Жук — канд. техн. наук, доцент кафедры "Системы автоматизированного проектирования" МГТУ им. Н.Э. Баумана. Автор более 80 научных работ в области математического обеспечения систем автоматизированного проектирования.
МГТУ им. Н.Э. Баумана, Российская Федерация,105005, Москва, 2-я Бауманская ул., д. 5.
D.M. Zhuk — Cand. Sci. (Eng.), assoc. professor of "Systems of Computer-Aided Design" department of the Bauman Moscow State Technical University. Author of more than 80 publications in the field of mathematical support of computer-aided design systems. Bauman Moscow State Technical University, Vtoraya Baumanskaya ul. 5, Moscow, 105005 Russian Federation.
Валентина Владимировна Глазкова — канд. физ.-мат. наук, ассистент кафедры "Автоматизация систем вычислительных комплексов" факультета "Вычислительная математика и кибернетика" МГУ им. М.В. Ломоносова. Автор более 10 научных работ в области вычислительной математики, программирования и систем автоматической категоризации текстовой информации.
МГУ им. М.В. Ломоносова, Российская Федерация, 119991, Москва, ГСП-1, Ленинские горы, д. 1, стр. 52.
V.V. Glazkova — Cand. Sci. (Phys.-Math.), assistant lecturer of "Automation of Computing Complex Systems" department of "Computational Mathematics and Cybernetics" faculty of the Lomonosov Moscow State University. Author of more than 10 publications in the field of computational mathematics, programming, and systems of automatic categorization of textual information.
Lomonosov Moscow State University, Leninskie Gory 1, stroenie 52, Moscow, 119991 Russian Federation.
