Болтенков Степан Анатольевич, Сибирский федеральный университет (СФУ)

Введение

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

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

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

Одним из таких подходов может стать объектно-ориентированный подход (ООП) с концепцией слоев [1]. Концепция слоев (layers) – одна из моделей, используемых разработчиками программного обеспечения для разделения сложных систем на более простые части. Описывая систему в терминах архитектурных слоев, удобно воспринимать составляющие ее подсистемы в виде «слоеного пирога».

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

В данной работе будет рассмотрена модель одной из химических реакций и на ее примере представлены все пункты применения ООП. Будет проведен сравнительный анализ возможных и известных подходов для решения поставленной задачи.

1. Постановка задачи

Для формирования четкого представления о предложенном методе необходимо подробно рассмотреть предметную область, а именно, параметрический анализ структуры Тьюринга [2]. В общем случае под термином структура Тьюринга понимают систему дифференциальных уравнений определенного вида. Для реакции двух веществ с одномерной диффузией система уравнений будет иметь следующий вид [3]:

 (1)

где

 (2)

 (3)

Начальные данные:

 .

Краевые условия для отрезка (0,l) с непроницаемыми стенками

 .

Дополнительное условие, следующее из предметной области, состоит в том, что концентрация вещества не может превышать 1 или быть отрицательной.

 (4)

В качестве примера рассмотрим реакцию вида:

1. Z↔X1

2. X1+2Z→3Z

3. Z↔X2

которую описывает система дифференциальных уравнений [4]

  (6)

где

  (7)

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

Разделим параметрический анализ на три основных этапа.

1. Нахождение стационарных точек системы

2. Исследование устойчивости стационарных точек

3. Бифуркационный анализ.

Нахождение стационарных точек ( далее - ст.с) заключается в поиске решений системы уравнений:

  (8)

Вторым шагом исследования системы (1) является определение характера особых точек и построение параметрических кривых.

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

  (9)

где - это параметры из (2), (3).

Последним этапом параметрического анализа является построение бифуркационных кривых: кривой кратности стационарных состояний LΔ: Δ=0 и кривой нейтральности Lσ,: σ=0.

Опишем процедуру построения этих кривых.

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

Для этого запишем линеаризованную относительно отклонений систему:

  ,

  (10)

Сформируем матрицу Якоби с элементами

   (11)

где

  

   (12)

Будем искать решение в виде:

 ,  (13)

при котором характеристическое уравнение примет вид:

  (14)

  где

  (15)

 (16)

Значение x определено ст.с. Устойчивость ст.с. определяется собственными числами матрицы Якоби. Для исследования устойчивости достаточно исследовать знак σ и Δ.

Выделим из два параметра p1 и р2 и построим линии LΔ и Lσ в плоскости этих параметров. Граница области множественности LΔ определяется, как решение системы уравнений:

H(x,p1,p2) = 0 

Δ(x,p1,p2)=0  

Таким образом, кривая краткости стационарных состояний LΔ в плоскости параметров (p1,p2) выписана в параметрическом виде:

P2=ξ2(x)

P1= ξ1(x, ξ2(x,)), 

2. Объектная схема программного обеспечения

Общую структуру программного обеспечения можно отразить в виде схемы (Рис №1).

Рис. №1. Объектная схема приложения

Рассмотрим последовательно все блоки: Модель, Численные алгоритмы, Интерфейс

2.1 Модель

Данный блок представляет собой класс, методами которого являются функции (2,3) и два дополнительных метода (15,16). Параметрами класса выступают коэффициенты уравнений (2,3), представленных в виде динамических массивов, а также их предельные значения, полученные из предметной области. В результате, универсальным интерфейсом класса будет являться набор методов и параметров, которые инкапсулируют все необходимые данные. Приведем пример в синтаксисе языка Pascal.

TFtp = function(X : TPoint) : Real of Object;

TModel = class (TObject)

vF : TFtp;

vG : TFtp;

vSigmaF : TFtp;

vDeltaF : TFtp;

constructor Create;

procedure AddK(pK,pKmin,pKmax : Real; Ks : ShortString);

function GetK(pI : integer) : TK;

procedure SetK(pI : integer; pK : TK);

function KCount : integer;

function GetKParam(pI : integer; pX : Real) : Real;

procedure GetKFSigma(pI,pJ : integer; pX: Real; var pK,pS : Real);

procedure GetKFDelta(pI,pJ : integer; pX: Real; var pK,pS : Real);

Из выше приведенного примера можно видеть, что в интерфейсе класса TModel отражены все необходимые методы для решения поставленной задачи (Таблица №1).

Метод класса TModel	Пояснение
GetKParam	Функция предназначенная для проведение параметрического анализа и построения зависимостей по формуле (9)
GetKFSigma	Процедура предназначена для получения параметрических зависимостей кривой нейтральности (Lσ) кратности относительно двух параметров
GetKFDelta	Процедура предназначена для получения параметрических зависимостей кривой кратности (LΔ) относительно двух параметров.

Таблица №1. Основные методы класса TModel

Остальные методы и параметры класса TModel предназначены для задания начальных условий и взаимодействия с внешними классами.

2.2 Численные методы

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

2.2.1 Нахождение стационарных точек системы

В этом случае достаточно применение метода последовательных разбиений на отрезки. Условие (4) позволяет взять достаточно малый шаг h от которого и зависит точность решения.

2.2.2 Исследование устойчивости стационарных точек

Задачу по исследованию стационарных точек необходимо разбить на три стадии. На первой стадии выясняется условие вхождения исследуемого параметра в систему. Для большей наглядности рассмотрим данный этап на примере построения параметра из (7). Благодаря тому, что параметр  входит лишь в одно уравнение системы (7), мы можем получить значения x1 x2 из второго уравнения методом дихотомии, с высокой степенью точности. Подставляя, таким образом, x1 x2 в функцию F и повторно применяя метод дихотомии, но уже для  мы получим искомую зависимость.

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

2.2.3 Бифуркационный анализ

Для бифуркационного анализа рассмотрим методику построения кривой нейтральности, используя систему уравнений (7). Для примера возьмем зависимость:

. (17)

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

2.3 Интерфейс

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

 I. Класс TGraph - предназначен для отображения графиков получаемых при построении параметрических зависимостей (Рис №2), а так же фазовых портретов при решении системы дифур и графиков динамики изменения концентрации вещества при численном счете (1). Данный класс во внешнем интерфейсе имеет лишь один метод AddXY(X,Y : Real), предназначенный для внесения очередной точки-параметра на плоскость. Масштабирование, фокусировку, выделение области и смещение класс TGraph осуществляет самостоятельно.

Рис. №2. Пример работы класса TGraph

II. Класс TPlant – предназначен для отображения степени концентрации вещества на пластине при решении (1) в двумерном случае. Данный класс, как и предыдущий содержит в себе только один метод предназначенный для внесения точки на плоскость и применяется уже не при исследовании, а расчете модели. На Рис №3 показан пример такого отображения.

Рис. №3. Пример работы класса TPlant

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

Рис. №4. Пример работы класса TCube

3. Сравнительный анализ предложенного метода и альтернативных решений

Альтернативных методов решения поставленной задачи может быть два, а именно: аналитические решения и решения с помощью прикладных продуктов, например Matlab.

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

К примеру, функция кривой нейтральности  для предложенного примера (6), выглядит следующим образом:

где

.

Поскольку такие зависимости нужно вывести не только для каждого параметра, но и для каждой бифуркационной кривой, то процесс расчета занимает значительную часть времени. При любом изменении уравнений (7) всю работу предстоит делать заново. Помимо «неудобства» использования аналитического метода, существует ряд параметров которые вообще невозможно получить аналитическим способом. В рассматриваемом примере (6) такой зависимостью является  для кривой кратности.

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

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

Заключение

На сегодня ООП является одним из основных методов в написании приложений, представленная работа показывает возможности которые открываются перед исследователем при его применении в задачах математического моделирования процессов химической кинетики. Механизм ООП реализован в программе «Параметрический анализ структуры Тьюринга», где он применяется не только для расчетов рассмотренных выше уравнений, но и для получения различных решений таких как: фазовые портреты и временные зависимости. Используя созданный механизм, специалист исследующий свою задачу, уже не тратит время на выведения нужных формул или подборке необходимого численного алгоритма для решения того или иного уравнения, а сосредотачивается непосредственно на рассматриваемом объекте.

Список литературы

1. Фаулер М. Архитектура корпоративных программных приложений. –М.: Издательский дом «Вильямс», 2006.-544с.

2. Пушкарёва Т.П. Параметрический анализ простейшей модели автоколебаний в реакции ассоциации: Препринт /Т.П. Пушкарёва, В.И. Быков // ВЦ СО АН СССР N13. – Красноярск, 1985.

3. Жаботинский А.М. Концентрационные колебания – М.:Наука,1974. 176 с.

4. Быков В.И. Моделирование критических явлений в химической кинетике /В.И. Быков. – М.: Наука,1988.

5. Быков В.И. Механизмы нелинейной диффузии и качественные свойства уравнений «кинетика + диффузия» /В.И. Быков, Г.П. Саркисян // Математические методы в химии: Сборник тезисов V Всесоюзной конференции. – Грозный, 1985. – С. 10-11.

6. Киселев Н.В Критические явления в системе «Автокаталитический осциллятор + диффузия» /Н.В. Кисилев // Вестник Красноярского государственного университета. – Красноярск, 2004. – С. 249-256.

7. Болтенков С.А. Бифуркационный анализ структуры Тьюринга / С.А. Болтенков // Молодежь и современные информационные технологии: Сборник трудов IV Всероссийской научно-практической конференции студентов, аспирантов и молодых ученых. – Томск, 2006. – С. 87-89.