Системный анализ устойчивости моделирования процесса переноса разностной схемой первого порядка Текст научной статьи по специальности «Математика»

UDC 517.929

517.929

СИСТЕМНЫЙ АНАЛИЗ УСТОЙЧИВОСТИ МОДЕЛИРОВАНИЯ ПРОЦЕССА ПЕРЕНОСА РАЗНОСТНОЙ СХЕМОЙ ПЕРВОГО ПОРЯДКА

Н. В. Кудинов, О. Р. Девдариани

Донской государственный технический университет, Ростов-на-Дону, Российская Федерация

kudinov_nikita@,mail. ш zzyzx88@mail.ru

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

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

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

Ключевые слова: уравнение переноса, разностная схема, аппроксимация, устойчивость разностных схем,

неустойчивость, условие Куранта-Фридрихса-Леви, метод гармоник.

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

1) разностная аппроксимация дифференциального уравнения на сетке — написание разностной схемы,

2) решение на ЭВМ разностных уравнений, представляющих собой системы линейных алгебраических уравнений высокого порядка специального вида.

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

SYSTEM ANALYSIS OF SIMULATION STABILITY OF THE PROCESS OF TRANSFERRING A DIFFERENTIAL SCHEME OF THE FIRST ORDER

N. V. Kudinov, O. R. Devdariani

Don State Technical University, Rostov-on-Don, Russian Federation

kudinov_nikita@mail. ru zzyzx88@mail. ru

The article provides the illustrations of the profiles of the time-carrying parameter of the corresponding balance condition, approxied by different difference schemes. The analytical solution to the equation is used to estimate the magnitude of the discrepancy obtained by approximating the left-side and right-side difference schemes. The cases of the appearance of instability in the simulation of balance by difference schemes and the mapping of these modes are shown graphically. In addition to stability analysis, the article presents a model that supposedly has a more complex analysis area for stability.

Keywords: Transfer equation, difference scheme, approximation, stability of difference schemes, instability, Courant-Friedrichs-Levi condition, harmonic method.

изменять структуру разностной схемы и её параметры, оценивать критерии качества аппроксимации законов сохранения.

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

На примере задачи [2], рассматриваются процессы сушки и, одновременно, перемещения неоднородной массы. Рассматривается технология, в ходе которой сыпучий материал передвигается по конвейеру, и в процессе передвижения производится его осушение. Пусть длина конвейера — Ь, а линейная скорость его перемещения — V.. Пусть конвейер покрывается к пространственными шагами. Время Л, в течение которого лента передвигается на один пространственный шаг, составляет:

I 1

Л=г •.

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

Общий вид такого уравнения определяется законом сохранении массы:

(1)

где и = х), t — переменная параметризующая изменение времени, х — пространственная координата. Уравнение (1) называется уравнением переноса. Оно является одним из уравнений, имеющих уже доказанное решение, и позволяющее тестировать разностные схемы на свойства аппроксимации устойчивости и сходимости. Решение уравнения переноса

и(х, 0 = и0(х — ат), (2)

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

ди , ди _

+ аТ" = 0,

д1 ох

1+1 I ] ]

+ = 0

- (3)

т п 4 '

Эта конечно-разностная замена дифференциальных зависимостей называется аппроксимацией. Для компьютерного моделирования удобнее итерационные процедуры, в которых на каждом шаге итерации известно приращение параметров. Для выражения в данной итерационной форме нужно избавится от дробной части, содержащей у/+"\ путём умножения на

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

У1 Ух Я

1+1 у/ - т -у/* ■ ау{+1-у{; п

!+1 — у!* ат(у/+1-у/); 1 п '

'/+1 — у!* .. ату/+1-ту/; 1 п '

'/+1 — У!* -а(ту/+1-ту/); п '

//+1 — у!* _ -ату/+1+ату/. п ' (4)

,/+1 _ -ату/+1+ату/ + /'

у1 к

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

Параметры окружающей среды, вызывающие неравновесную динамику. В качестве модели быстропротекающего изменения параметров во времени рассматривается функция Хэвисайда [4] (5).

е={?:I > (5)

Функция (5) моделирует однократное изменение условия окружающей среды для системы, что приводит к переходу системы в другое состояние за ограниченное время, при условии системной устойчивости. Примером задачи, которая приводит к уравнению (1) могут быть процессы транспортировки в химической технологии. Движущаяся среда поступает в трубопровод со скоростью V и выходит с другой стороны, перенося в себе растворённое вещество со скоростью движения среды (рис. 1). Например, концентрация вещества ф меняется во времени, а постоянное её численное значение переносится в пространстве. Изменение режима работа химического агрегата хорошо моделируется единичной функцией Хэвисайда. Для моделируемой подсистемы она может заменить модель входного процесса, определяя граничное условие уравнения переноса (1).

¡-2 ¡-1 I ¡+11 + 2 Рис. 1. Потоковая схема моделирования распределённого объекта

Описанное математическое и алгоритмическое обеспечение реализовано в виде программной системы, допускающей имитационные эксперименты исследователя с моделью и приближающей её разностной схемой. Проведённая серия компьютерных экспериментов показала, что решением уравнения переноса при заданном входном единичном воздействии через ограниченное время является фронт волны. Изучены профили пространственного распределения, в частности, фронт волны (рис 2). Установлено, что разностная схема «деформирует» профиль фронта волны.

Рис. 2. График входного воздействия по концентрации раствора на входной (по ходу движения среды) границе распределённого объекта

Рис. 3. Графики изменения концентрации среды во множестве внутренних сечений, покрывающих весь объект

При каждом шаге вычислений процесс переноса ограничен амплитудой входного воздействия. Это ограничение наблюдается при нарушении численной устойчивости моделирования. Эксперимент показал, что разностные схемы 1 -го порядка аппроксимации могут быть как устойчивыми, так и неустойчивыми. Одна из причин такой вариации — это зависимость знака обратной связи от направления шаблона разностной схемы. Ниже показана аппроксимация модели переноса левосторонней разностной схемой и её решение относительно приращений по времени:

}+1 } } } у1 -У\ + ау/-у/-1

п

0,

,!+1

—1-— + у/.

(6)

(7)

Ух п

Правосторонняя по пространству разностная схема может проявлять неустойчивость при изменении параметра, например, при увеличении скорости переноса). На рис. 3 показаны исходные параметры и результаты модельного расчёта.

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

Уравнение переноса Уравнение Пуассона Уравнение Навье-Стокса

0 5 10 15 20 Рис. 3. Параметры модели, при которых появляется неустойчивость

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

Проведём анализ устойчивости разностной схемы, представленной формулой (4), аппроксимирующей уравнение переноса. Числом Куранта у = называется соотношение параметров разностной схемы (т — шаг по времени, h — шаг по пространству) умноженное на скорость движения среды (а). Тогда аппроксимированное уравнение записывается в следующем виде:

yL+1 = -yyLi + уу!С + yL • (8)

Общее решение ищется в виде гармоники yjl = qieikv с амплитудой, зависящей от дискретного времени t = /г, где j = 0.. п, а п — это количество шагов по времени и по пространству, т.е. частота пространственной сетки; q — амплитуда, j, к — номера узлов пространственной сетки по времени, i — мнимая единица, ф— действительная частота гармоники.

Неравенство, выражающее условие не возрастания амплитуды, имеет левую часть:

ф ф л

q = sin j у — cos j у + 1 (9)

С учётом условия устойчивости разностных схем < 1 или условием неустойчивости > 1 поэтому амплитуда гармоники должны быть распределена в отрезке [—1; 1]

— 1 < sinфу — cosфу +1 < 1. (10)

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

Свойство математической нелинейности исходной модели препятствует применению подобного анализа методом гармоник. Одна из важнейших моделей уравнение Навье-Стокса [5] содержит нелинейность типа «произведение».

Р Gf + v v) = —grad p + r\Av + pf (11)

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

Заключение. Зависимость устойчивости разностных схем от их ориентированности относительно скорости распространения возмущений в распределённых в пространстве объектах был открыт в СССР и представлен в научных трудах Самарского А.А, Годунова С.К. , Рябенького В.С. Дальнейшие исследования строятся на предположении что при компьютерном моделировании нелинейных процессов условные разностных схемы теряют численную устойчивость и это нежелательное для моделирования явление не может быть исключено изначальным выбором шаблона и параметров разностной схемы. В настоящее время тестируются разностные схемы с переменной структурой шаблона и выявляются критерии численной диагностики качества аппроксимации законов сохранения.

Описанные предварительные исследования показывают, что разделение процесса моделирования на две фазы — математическое и компьютерное моделирование осложняет анализ проблемы устойчивости имитационной модели. Компьютерная программа, которую можно использовать как исследовательскую установку для свойств тестирования разностных схем в компьютерных моделях, должна отображать все промежуточные уравнения и зависимости, включаемые в фазу математического моделирования. Такая программа должна позволят оценивать невязку при аппроксимации законов сохранения и влияние на неё изменения структуры разностной схемы.

Библиографический список

1. Самарский, А. А. Методы решения сеточных уравнений / А. А. Самарский, Е. С. Николаев. — Москва : Наука, 1978. — 592 с.

2. Моделирование систем с распределенными параметрами при перемещении материальных масс [электронный ресурс] / Мухин О. И. Курс «Моделирование систем». — Режим доступа: http://stratum.ac.ru/education/textbooks/modelir/

lection18.html (дата обращения: 24.05.2017).

3. Самарский, А. А. Теория разностных схем / А. А. Самарский. — 3-е изд., испр. — Москва : Наука., 1989. — 616 с.

4. Егерев, В. К. Методика построения графиков функций / В. К. Егерев, Б. А. Радунский, Д. А. Тальский. — 2-е изд. — Москва : Высшая школа, 1970. — 152 с.

5. Самарский, А. А. Вычислительная теплопередача / А. А. Самарский, П. Н. Вабищевич. — Москва : Едиториал УРСС, 2003. — 785 с.