Метод моделирования динамики технических систем на основе формальных схем

МЕТОД МОДЕЛИРОВАНИЯ ДИНАМИКИ ТЕХНИЧЕСКИХ СИСТЕМ НА ОСНОВЕ ФОРМАЛЬНЫХ СХЕМ

, ,

МГТУ имени , [e-mail]

Введение

Системы моделирования динамических процессов в технических объектах прошли достаточно большой путь развития (работы в этой области начались еще в 60-х годах прошлого века) и использовались при решении ряда сложных технических задач. Однако в настоящее время подобные системы эффективно и широко используются только в ECAD (электронных CAD) при проектировании изделий электроники, микроэлектроники, радиоэлектроники и вычислительной техники. При проектировании же в других областях, и, в первую очередь, в MCAD (машиностроительных CAD) системах, а также при моделировании разнородных технических систем (электро-механических, электро-тепловых и т. п.) они используются эпизодически и не вносят существенного вклада в решение проектных задач. Можно выделить два подхода к моделированию динамики технических объектов.

Первый подход предполагает разработку специализированных программ моделирования для конкретных предметных областей: PSPICE в электронике; ADAMS, DADS, EULER в механике и гидравлике; МВТУ (моделирование в техничеких устройствах) для моделирования систем управления, в первую очередь, ядерными реакторами и т. п.

Второй подход предполагает разработку универсальных программ моделирования технических объектов разной физической природы на основе метода аналогий переменных и уравнений математических моделей таких объектов: программы ПА4, ПA6, ПA7, ПA8, ПA9, разработанные в МГТУ им. . Поскольку в настоящее время встают задачи комплексного моделирования сложных технических систем, а сложные технические системы, в свою очередь, состоят из физически разнородных технических объектов, то актуальность второго подхода значительно возрастает. В этих программах использовался метод формирования математических моделей динамических процессов в технических системах, при котором заключительная модель для разработки рабочих подпрограмм представляла собой систему нелинейных алгебраических уравнений (НАУ). Однако разработанные на основе этих методов перечисленные выше программы имеют ряд принципиальных недостатков. В первую очередь это ограниченность используемых в этих программах базисов переменных моделирования, что не позволяет практически реализовать новые неявные А-устойчивые методы интегрирования дифференциальных уравнений с порядком точности выше второго и ограничивают возможности разработчиков математических моделей элементов моделируемых объектов. В данной работе описывается обобщенный подход к моделированию динамики разнородных технических систем на основе метода аналогий и обобщенных формальных схем (ФС), а также описывается обобщенный базис переменных моделирования таких схем, при этом заключительная математическая модель для разработки рабочих подпрограмм представляет собой систему нелинейных дифференциально-алгебраических уравнений общего вида.

Обобщенные формальные схемы

При математическом моделировании состояние любого технического объекта (ТО) в конкретный момент времени можно описать множеством некоторых физических величин, имеющих как постоянные значения, так и переменные значения в течение заданного времени моделирования. Обозначим множество переменных модели ТО вектором Vто(t). Конечной математической моделью динамических процессов в ТО в общем случае будет некоторая система дифференциально-алгебраических уравнений (ДАУ), замкнутая относительно этих переменных:

F(Vто, dVто/dt, t)=0 (1)

Вместе с тем любой ТО можно представить как множество элементов, связанных между собой множеством связей, поэтому его можно отобразить схемной моделью в виде некоторой обобщенной формальной схемы (ФС), состоящей из набора элементов, связанных между собой в некоторых узлах этой схемы. Такая схемная математическая модель технического объекта должна отображаться аналогичным множеством переменных Vфс(t) -> Vто(t), а конечной математической моделью такой ФС должна быть система ДАУ, аналогичная системе (1). Для моделирования динамики разнородных технических систем разработатывается программная система моделирования обобщенных ФС (ПА10), с помощью которой можно будет моделировать разнообразные технические объекты, описываемые системой ДАУ вида (1).

Элементный состав обобщенной ФС и типы связей между элементами должны давать возможность моделирования как любых дифференциально-алгебраических уравнений вида (1), так и разнообразных элементов технических объектов любой физической природы (электрических, магнитных, механических, гидравлических, тепловых, оптических, акустических и др.) и разнообразных связей между такими элементами. В работе [1] предлагается строить ФС на основе 6 типов двухполюсных элементов и универсального многополюсного элемента.

Состояние каждого узла связи в ФС в некоторый момент времени t описывается двумя типами безразмерных переменных – рис. 1.

Рис. 1. Узел формальной схемы и переменные, описывающие его состояние

1) Переменная типа потока, втекающего в узел связи или вытекающего из узла связи Imk(t). Эти переменные удовлетворяют уравнению непрерывности потоков для каждого узла ФС, которое будет базовым для получения модели ФС на основе моделей отдельных элементов ФС:

Im1+Im2+…+Imk+… = 0

2) Переменная типа потенциала Pm(t), которая определяется по отношению к некоторому базовому узлу, потенциал которого принимается равным нулю. Разновидностью переменной этого типа будет переменная типа напряжения (разности потенциалов) Umn(t)=Pm(t)-Pn(t), которая определяется как разность потенциалов между узлами m и n.

Простейшими базовыми элементами ФС будут 6 широко известных типов двухполюсников.

1. Элемент типа зависимый источник потока.

2. Элемент типа емкость

3. Элемент типа сопротивление.

4. Элемент типа проводимость.

5. Элемент типа зависимый источник напряжения.

6. Элемент типа индуктивность.

Динамику любой технической системы, описываемой обыкновенными дифференциальными уравнениями, можно промоделировать с помощью определенных выше двухполюсников, однако для моделирования систем большой размерности вводится еще один базовый элемент – функциональный многополюсник (ФМ), с помощью которого при моделировании систем больщой размерности можно применять метод декомпозиции такой системы на части.

Функциональный многополюсник

Функциональным многополюсником (ФМ) назовем элемент ФС, который представляет собой схему соединения двухполюсников и внутренних ФМ, в которой выделены внешние полюса Pi для связи с другими элементами ФС – рис. 2.

Рис. 2. Функциональный многополюсник ФС и переменные, описывающие его полюса

Каждый i-ый полюс описывается двумя переменными:

1. Поток iPi, втекающий в i–ый полюс, который назовем полюсным потоком.

2. Потенциал i-го полюса uPi, который назовем полюсным потенциалом.

Обозначим через N общее количество внешних полюсов ФМ. Предположим, что ФМ состоит только из n двухполюсников, которые описываются вектором внутренних переменных ФМ – Vвнфм, размерностью 2n+k+m, где k – число внутренних узлов ФМ, а m – число элементов типа C и L внутри ФМ. Математической моделью такого ФМ будут внутренние уравнения для внутренних двухполюсников ФМ, размерностью 2n+k, а также N полюсных уравнений для суммы потоков, втекающих в каждый i-ый внешний полюс внутри ФМ. Таким образом, общее количество переменных ФМ будет на N больше, чем общее число уравнений, составляющих математическую модель ФМ. При включении ФМ в общую схему ФС недостающие N линейно-независимых уравнений будут получены для суммы потоков, втекающих в каждый i-ый внешний полюс снаружи ФМ.

В качестве математического ядра разрабатываемых программ моделирования динамики технических систем на основе ФС предполагается использовать стандартную программу решения систем ДАУ общего вида – DMAN [1], которая была разработана, не привязываясь к выбранному базису переменных моделирования.

Программа DMAN для интегрирования систем ДАУ

Программа DMAN предназначена для интегрирования систем ДАУ общего вида, не разрешенных относительно производных:

F(XP, X,Y,t) = 0 (2)

где X - вектор дифференцируемых переменных системы (2) размерностью m; XP - вектор производных дифференцируемых переменных по времени размерностью m, т. е. XP = dX/dt; Y - вектор алгебраических переменных системы (2) размерностью k; t - независимая переменная (например, время); F - вектор-функция системы (2) размерностью n, где n = (m + k). Заданы начальные условия для вектора дифференцируемых переменных: X0 = X(t0), где t0 - начальный момент времени интегрирования. Начальные значения остальных переменных, т. е. XP0 и Y0 рассчитываются в программе DMAN перед началом интегрирования автоматически.

В процессе интегрирования система (2) на каждом шаге интегрирования дополняется системой:

G(XP, X,h) = 0, (3)

где h - шаг интегрирования; G - вектор-функция размерностью m, которая соответствует выбранному методу интегрирования.

В программе DMAN реализованы 3 метода интегрирования систем ДАУ:

М1 - А-устойчивый неявный метод Эйлера первого порядка точности;

М2 - А-устойчивый неявный метод второго порядка точности;

М3 - новый А-устойчивый неявный метод четвертого порядка точности;

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

Во всех методах на каждом шаге интегрирования решается система нелинейных алгебраических уравнений относительно переменных X, XP, Y методом Ньютона, для которого необходимо вычислять матрицу Якоби для системы (2) и системы (3). Матрица Якоби для системы (3) вычисляется в программе DMAN автоматически. Матрица Якоби для системы (2) должна вычисляться в специальной подпрограмме пользователя, наряду с вычислением вектор-функции F.

Начальным приближением для метода Ньютона являются значения переменных в начале шага интегрирования, поэтому, если на некотором шаге нет сходимости итераций, то текущий шаг h уменьшается вплоть до минимального, после чего выдается сообщение об ошибке (обычно из-за ошибок в формулах для определения элементов матрицы Якоби в подпрограмме пользователя или из-за расходимости самой системы ДАУ).

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

Примеры получения моделей простых электронных схем в новом обобщенном базисе

Рассмотрим примеры получения моделей ФС для реальных технических объектов - двух электрических схем, одной линейной и одной нелинейной. Сначала рассмотрим пример получения модели ФС для линейной электрической схемы, приведенной на рис. 3, в новом обобщенном базисе.

Рис. 3. Пример линейной электронной схемы

Вектор базисных переменных для этой ФС будет иметь размерность 12 и включает в себя следующие переменные:

Vфс = (uE, iE, uR, iR, uC, iC, uL, iL, ф1,ф2,xpC, xpL),

здесь и далее xpC – производные от напряжений на соответствующих емкостях – duC/dt, а xpL – производные от токов на соответсвующих индуктивностях – diL/dt.

Система ДАУ для этой ФС будет иметь размерность n=10, m=2 и включает в себя следующие уравнения:

uE=f(t) или uE-f(t)=0

uE=ф1 или uE-ф1=0

uR=R iR или uR-R iR=0

uR=ф1-ф2 или uR-ф1+ф2=0

iC=C xpC или iC-C xpC=0

uC=ф2 или uC-ф2=0

uL=L xpL или uL=L xpL=0

uL=ф2 или uL-ф2=0

iE-iR=0

iR-iC-iL=0.

Отметим, что в базисе узловых потенциалов размерность системы ДАУ будет равна 2 (в пять раз меньше!), но при этом возможности моделирования отдельных элементов ФС будут ограничены (например, нельзя моделировать идеальные источники напряжения).

Для этой задачи вектор X=(uC, iL), вектор XP=(xpC, xpL), а вектор Y= (uE, iE, uR, iR, iC, uL, ф1,ф2). Вектор-функция F включает в себя 10 вышеприведенных уравнений. Полученная система ДАУ была решена с помощью программы DMAN. Полученное решение сравнивалось с аналитическим, при этом были получены численно более устойчивые решения при изменении коэффициентов этой системы ДАУ, по сравнению с системой ДАУ, полученной в базисе узловых потенциалов.

Рассмотрим пример получения модели ФМ для нелинейной модели диода.

На рис. 4 показаны эквивалентная схема диода из двухполюсных элементов и его представление как ФМ.

Рис. 4. Диод как функциональный многополюсник

Вектор переменных для ФМ диода будет иметь размерность 9 (из них первые 5 переменных являются внутренними) и включает в себя следующие переменные:

Vфмдиод = (uCд, iCд, xpCд, uIд, iIд, uк, iк, uа, iа)

Система ДАУ для ФМ диода (математическая модель диода в обобщенном базисе) будет иметь размерность 6 и включает в себя следующие уравнения (первые 4 уравнения – для 2 внутренних двухполюсников ФМ диода):

iCд-Cд0 xpСд=0

uCд-uа+uк=0

iIд-(iт exp(uIд/фт)+1)=0

uIд-uа+uк=0

iа-iCд-iIд=0

iк+iCд+iIд=0

В этих уравнениях iт, Cд0, и фт – это параметры данной модели диода, Сд – значение емкости диода.

Такие уравнения будем включать в библиотеку моделей ПA10 и использовать при получении моделей сложных ФС, содержащих диоды.

Рассмотрим пример получения в обобщенном базисе с использованием ФМ диода модели ФС для нелинейной электрической схемы, приведенной на рис. 5.

Рис. 5. Нелинейная электронная схема с диодом

Вектор переменных Vфс включает в себя первую очередь переменные диода, для которых узловые потенциалы заменяются на потенциалы узлов подключения диода:

Vфмдиод = (uCд, iCд, xpCд, uIд, iIд, ф3, iк, ф2, iа)

Затем в вектор переменных ФС добавляются остальные переменные, в результате вектор ФС будет иметь размерность 19 и включает в себя следующие переменные:

Vфс=( uCд, iCд, uIд, iIд, ф3, iк, ф2, iа uE, iE, uR1, iR1, uC, iC, uR2, iR2, ф1, xpCд, xpC).

Система ДАУ для этой ФС будет иметь размерность 17 и включает в себя следующие уравнения:

Сначала в эту систему добавляются 6 уравнений диода, в которых узловые потенциалы заменяются на узлы подключения диода:

iCд-Cд xpСд=0

uCд-ф2+ф3=0

iIд-(iт e(uIд/фт)+1)=0

uIд-ф2+ф3=0

iа-iCд-iIд=0

iк+iCд+iIд=0

Затем добавляются 2 уравнения для суммы токов в узлах подключения диода:

iR1-iа=0

iк+iC1+iR2=0

Затем добавляются 8 уравнений для 4 двухполюсников:

uE-f(t)=0

uE-ф1=0

uR1-R1 iR1=0

uR1-ф1+ф2=0

iC-C xpC=0

uC1-ф3=0

uR2-R2 iR2=0

uR2-ф3=0

И в завершении формирования замкнутой системы ДАУ добавляется 1 уравнение для суммы токов в узле ф1:

iE-iR1=0

Для этой задачи вектор X=(uCд, uC), вектор XP=(xpCд, xpC), а вектор Y= (iCд, uIд, iIд, ф3, iк, ф2, iа uE, iE, uR1, iR1, iC, uR2, iR2, ф1). Вектор-функция F включает в себя 17 вышеприведенных уравнений. Полученная система ДАУ была решена с помощью программы решения систем ДАУ общего вида DMAN. Полученное решение сравнивалось с решением, полученным для метода узловых потенциалов в программе ПA7 [2], при этом были получены численно более устойчивые решения при изменении коэффициентов этой системы ДАУ, по сравнению с системой ДАУ, полученной в базисе узловых потенциалов.

Заключение, выводы

В результате проведенных исследований выполнен анализ и обоснование выбора нового обобщенного набора базовых элементов ФС и нового обобщенного базиса переменных ФС, в качестве основы для разработки новой, мобильной, открытой программы моделирования динамики технических систем - ПА10. Программа ПA10 будет базовым математическим ядром системы моделирования динамических процессов нового поколения и будет иметь следующие принципиальные особенности:

- использование новых неявных А-устойчивых методов интегрирования систем ДАУ с порядком точности выше второго;

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

Дальнейшие исследования будут включать в себя:

- разработка структуры всей системы моделирования динамики технических объектов и организации межмодульных связей в этой системе;

- разработка новых методов решения систем линейных алгебраических уравнений с учетом очень сильной разреженности этих систем для обобщенного базиса при сохранении численной устойчивости вычислений;

- разработка новой методики тестирования программы, основанной на всестороннем тестировании математического ядра на основе методов тестирования математических программ типа Mathlab, Mathematica и т. п.;

- использование универсальных СУБД для хранения библиотек моделей элементов моделируемых технических объектов;

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

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

ЛИТЕРАТУРА

1.  , , Папсуев моделирования динамики технических объектов в программе ПА10 // Информационные технологии№ 8. - С. 6–14.

2.  Маничев алгоритмы для программ анализа динамики технических систем // Вестник МГТУ, сер. ПриборостроениеВып. 1. - С. 48-56.