Чепасов В.И., Токарева М.А.
Оренбургский государственный университет
АЛГОРИТМИЧЕСКАЯ И ПРОГРАММНАЯ РЕАЛИЗАЦИЯ РЕШЕНИЯ УРАВНЕНИЯ КОЛМОГОРОВА - ЧЕПМЕНА
Рассмотрена методика получения системы конечно-разностных уравнений для уравнения Колмогорова - Чепмена.
Представлена алгоритмическая реализация методики, дано описание программной реализации алгоритма, приведены результаты решения.
Уравнение Колмогорова - Чепмена имеет вид:
А э2и(х,р + А э^и(хл)+А э3и(х,р и )_0 А1 Э\ + А2 ЭxЭt + Аз э2xэt и(хД)_0’ (1)
где в (1)
А1, А2, А3 - коэффициенты,
Л1=1/д, Л2=а/д, Л3=Ъ/(2д), Ь<0, а, Ь, д задаются.
Для получения системы конечно-разностных уравнений заменим в (1) производные конечно-разностными соотношениями:
Эи(,х)_ и( + Дt,x)-и(t,x)
Эt
Дt
Э ри(,х)ї Э
Эt I дt I Эt
U(t + Дt,x)- и(,х)
Дt
U(t + 2Дt,x)- и( + Дt,x)- (u(t + Дt,x)- и (t, х)) Дt xДt
_ и( + 2Дt,x)- 2U(t + Дt,x)+ и(,х)
_ ДхД ’
" U(t + Д1,х)- и(,х) _
(2)
Э ри(и)У Э
Эx I Эt I Эx
Дt
1
Дx xДt
(u(t + Д1,х + Дx)- и( + Д^х)-I- и(,х + Дx)+ и(х)),
(3)
Эx
Э I ЭU(t,x)
ЭХ I-ЭГ-
-[u(t + Дt,x + 2Дх)-
ДxДxДt
- и( + Д1,х + Дх)- ( + Дt,x + Дх)- U(t + Д1,х))-
- (и (, х + 2Дх)- и (,х + Дх)) + (и(,х + Дх)- и (х))] _
1
-(u(t + Дt,x + 2Дх)- 2U(t + Дt,x + Дх)+
(4)
ДxДxДt
+ и(t + Д^х)- и(t,x + 2Дх)+2И (,х + Дх)- и(t,x)) Координаты узлов сетки разбиения:
^ _(1 -l)xДt, хj _(- 1)хДх, 1 _ 1,то _ 1,п,
т - число точек разбиения по оси времени ОТ,
п - число точек разбиения по оси 0Х,
Дг - шаг разбиения по оси времени ОТ, Дх - шаг разбиения по оси ОХ.
Введем обозначение: и(^)_ ии.
С учетом этого и выражений 2, 3, 4 для производных уравнение (1) приводится к следующему конечно-разностному уравнению в узле (і,]):
Кіх иц + К 2х и^ + Кзх иад + К 4х и^ +
+ К5 хи^ + К6 хЦ+ц+2 + К7 х_ 0 , (5)
где в (5)
і - номер временного отсчета для узла сетки разбиения(отсчет с 1 = 0),
] - номер отсчета по оси ОХ для узла сетки разбиения(отсчет с х = 0),
Кі _
+
Аі
Дt xДt Дх xДt Дх хДх xДt А
-і,
К з _
Дt xДt А,
+
- 2Аі
ДtxДt ДxxДt Дх хДххД^ А 2А
К5_ - А2
Дх xДt ДххДх xДt 2А3
+
Дх xДt Дх хДх xДt’
с6 ____А_____,
6 ДххДххД^
А
Дх хДх xДt
дt - шаг разбиения по оси времени ОТ, дх - шаг разбиения по оси ОХ. Поскольку количество точек разбиения по оси ОХ п, количество точек разбиения по оси ОТ т, значения и(1, х) на границах заданы, то общее число неизвестных в системе линейных алгебраических уравнений, получающейся из (5), будет (п-2)*(т-2).
Для конкретизации построения этой системы рассмотрим следующую область разбиения:
Э
і
I т
(1=1) 5
(1=0.75) 4 . *7 *8 *9
(1=0 . 50) 3 . *4 * 5 * 6
(1=0 . 25 ) 2 . * і * 2 *3
(1=0) 1 .
1 2 3 4
X
(х=0) (х=0.25) (х=0.5) (х=0.75) (х=1) Согласно представленной области разбиения количество точек разбиения по оси ОХ п = 5, количество точек разбиения по оси ОТ т = 5.
Значения и(1, х) на границах 1 = 0, х = 0, 1 = і, х = і известны.
Тогда неизвестными значениями и(1, х) будут значения в узлах, обозначенных звездочками.
Неизвестные в этих узлах имеют номера
1, 2, 3, 4, 5, 6, 7, 8, 9.
Для нахождения неизвестных составим систему линейных алгебраических уравнений по (5) для подобласти:
0 < = х < = 0.5, 0 < = 1 < =0.5. Координаты узлов этой подобласти и номера уравнений:
(1, 1) -номер уравнения 1,
(1, 2) -номер уравнения 2,
(1, 3) -номер уравнения 3,
(2, 1) -номер уравнения 4,
(2, 2) -номер уравнения 5,
(2, 3) -номер уравнения 6,
(3, 1) -номер уравнения 7,
(3, 2) -номер уравнения 8,
(3, 3) -номер уравнения 9.
Первая координата - номер узла по оси ОТ, вторая координата - номер узла по оси ОХ.
Номер строчки матрицы системы (номер уравнения в подобласти) будет определяться:
пуг = (і-і)*(п-2)+], (6)
где в (6)
і - номер отсчета узла по оси ОТ, і = і, 2....(т-2)
] - номер отсчета узла по оси ОХ, ] = і, 2.. (п-2).
Приведем уравнение для узла с координатами 1 = 0, х = 0.
Координаты этого узла в номерах отсчетов:
1 = 1 - номер отсчета по ОТ,
] = 1 - номер отсчета по ОХ.
Номер уравнения в этом узле согласно (6): пуг = (1-1)*(п-2)+] = 1.
Тогда с учетом (5) уравнение для этого узла:
К1 х И1Д + К 2 х И31 + К 3 х И 21 + К 4 х И 22 +
+ К5 хИ12 + К6 хИ23 + К7 хИ13 _0, (7)
В уравнении (7) и1Д,и 31,и 21 ,и12 ,и13- заданные граничные значения и(1:,х).
Соответственно К1 х И1,1 + К 2 х И 31 + К 3 х и 21 + К 5 х и12 + К 7 х и13 мы переносим в правую часть уравнения 1 в качестве свободного члена.
Неизвестными в уравнении (7) будут и 22,
И 2,3.
Согласно области разбиения неизвестное И 22 будет иметь номер 1, неизвестное И 23 -номер 2.
То есть, если неизвестное имеет координаты узла (1, ]), то его порядковый номер будет определяться соотношением:
пп = (1-2)*(п-2)+]-1. (8)
С учетом (6) номер строчки (уравнения) в матрице системы линейных алгебраических уравнений для узла (1 = 1, ] = 1) будет: пуг = (1-1)*(п-2)+] = 0*3+1 = 1
Номера неизвестных в этом уравнении согласно (8):
- узел (2,2), пп = (1-2)*(п-2)+]-1 = 1,
- узел (2,3), пп = (1-2)*(п-2)+]-1 = 2.
Тогда элементы матрицы системы конечно-разностных уравнений (системы линейных алгебраических уравнений) для первого уравнения будут:
Су = К 4,
С1,2 = К6 .
Здесь первый индекс - номер строчки (уравнения), второй индекс - номер столбца (неизвестной).
Аналогично составляются уравнения для остальных восьми узлов подобласти: (1, 2), (1, 3), (2, 1), (2, 2), (2, 3), (3, 1), (3, 2), (3, 3).
При алгоритмической и программной реализации рассмотренной методики построения системы конечно-разностных уравнений вводились массивы номеров строчек, столбцов, значений граничных узлов.
Если номер строчки (по ОТ) и столбца (по ОХ) рассматриваемого узла при построении системы совпадали с номером строчки и столбца граничного узла, то соответствующий этому узлу член в уравнении переносился в правую часть уравнения (в свободный член).
Согласно выражениям (6) и (8) для определения номеров уравнений и номеров неизвестных в системе конечно-разностных уравнений минимальная область разбиения будет содержать три точки разбиения по оси ОТ и три точки разбиения по оси ОХ. Т о есть т = 3, п = 3.
В этом случае система конечно-разностных уравнений будет содержать всего одно уравнение.
Полученная система линейных алгебраических уравнений решалась методом Гаусса.
Поиск ненулевого диагонального элемента в прямом ходе метода Гаусса осуществлялся как по строчкам, так и по столбцам.
Предусмотрено решение системы, состоящей из одного уравнения.
Осуществляется проверка решения системы.
Составление матрицы системы конечноразностных уравнений и решение этой системы осуществляет программа pish.exe (исходный модуль pish.cpp).
Для функционирования программы pish.exe необходимо в текстовом редакторе создать файл wwyг следующей структуры:
- в первой строчке через пробел идут значения количества точек по горизонтали (ось ОХ) и количества точек по вертикали (ось Т). При этом значение произведения количества точек по горизонтали на количество точек по вертикали не должно превышать 9ОО. Это связано с максимальным порядком системы конечно-разностных уравнений - 9ОО.
Минимальное количество точек разбиения по оси ОТ и по оси ОХ будет 3;
- во второй строчке идут через пробел значения коэффициентов уравнения а, Ъ, д. При этом коэффициент Ъ берем со знаком минус. Значения коэффициентов должны обязательно иметь не более четырех знаков (включая сюда знак числа!) до десятичной точки.
То есть шаблон вводимых значений коэффициентов------.—...
Если шаблон не будет выдержан, то программа pish.exe работать будет, а программа virttext.exe - нет;
- в третьей строчке идут через пробел значения правых границ по координате Х>0 (отсчет от нуля! Не более пяти знаков до десятичной точки, включая в эти пять знак числа) и по координате Т>0 (отсчет от нуля! Не более пяти знаков до десятичной точки, включая в эти пять знак числа), значение на левой границе по Х (левая граница по Х = О), значение на правой границе по Х.
Сначала идет правая граница по Х, потом через пробел - правая граница по Т, далее через пробел значение на левой границе по Х (левая граница по Х = О), потом через пробел значение на правой границе по Х.
Если шаблон не будет выдержан, то программа pish.exe работать будет, а программа virttext.exe - нет;
- в четвертой строчке идут через пробел значения математического ожидания и среднего квадратического отклонения.
Далее в той же структуре (начиная с первой строчки!) в следующих строчках идут данные для других вариантов счета.
После данных по всем вариантам идет обязательная строчка с нулем.
После запуска программы pish.exe появится окно с меню, в котором будет всего один пункт run.
Прожимаем run. В появившихся окнах сообщений нажимаем ОК.
Закрываем окно программы pish.exe.
Программа pish.exe создаст два файла yrawn, grafpish.
В файле yrawn находятся результаты счета. Этот файл смотрим любым текстовым редактором.
В файле grafpish находятся данные для построения графиков.
Для вывода графиков на экран и для вывода графиков в EXCEL необходимо после программы pish.exe запустить программу virttext.exe.
После запуска программы virttext.exe появится окно с меню, в котором будут пункты TEXT и START.
Нажимаем пункт TEXT. На экране появится график.
Нажимаем пункт START. График исчезнет. Опять нажимаем пункт TEXT, потом START.
И так до тех пор, пока на экране не появится сообщение о конце просмотра всех графиков.
После этого закрываем окно программы virttext.exe.
В результате работы программы virttext.exe будет создан файл prograf.xls.
Этот файл открываем в EXCEL и строим графики.
Открыть файл prograf.xls можно следующим образом:
1. Загрузить EXCEL
2. Открыть файл prograf.xls в среде EXCEL При таком открытии файла prograf.xls в
EXCEL появится окно:
Мастер текстов шаг 1 из 3 Укажите формат данных Указываем с разделителем Нажимаем «Далее»
Появится окно Мастер текстов шаг 2 из 3 Символом разделителем является Указываем точку с запятой, знак табуляции не указываем
Нажимаем «Готово»
После такого преобразования записываем преобразованный файл и в мастере диаграмм выводим точечные диаграммы (выбираем графики).
При записи преобразованного файла появится окно, в котором надо выбрать «НЕТ» и далее записать (на все вопросы «Да»).
Если мы подвели курсорную строчку к файлу prograf.xls и нажали «Enter», то далее в среде EXCEL надо опять открыть файл prograf.xls с тем, чтобы появилось окно «Мастер текстов шаг 1 из 3». Далее работаем так, как было описано выше.
После записи преобразованного файла в мастере диаграмм выводим графики для соответствующих сечений по Т и по X.
Для этого выбираем для вывода графика два соответствующих столбика с числовыми значениями (выделяем эти столбики).
У нас первый столбик - это значения X, а второй столбик - значения Y.
После выделения двух столбиков запускаем мастер диаграмм, выбираем в этом мастере точечные диаграммы и выводим их на экран.
Например, необходимо вывести график по следующим данным в EXCEL:
10
time sec-
01 II 2 t II и a 4 0 0
0,0000 1,0000
0,0417 1,1939
0,0833 1,3828
0,1250 1,5653
0,1667 1,7399
0,2083 1,9052
0,2500 2,0599
0,2917 2,2027
0,3333 2,3323
0,3750 2,4475
1.0 b= -11.5 q=
Это данные для временного сечения под номером 2, значение t = 0.04.
То есть левый числовой столбец будет содержать значения отсчетов по оси 0X, а правый числовой столбец будет содержать значения u(x,t=0.04) по этим отсчетам.
Соответственно при построении графика по оси 0X будут значения X, а по оси 0Y - значения u(x,t=0.04).
Для построения графика выделяем в EXCEL эти два столбца и в мастере диаграмм, как было указано выше, строим график.
Рассмотрим сечение в EXCEL:
10
X sections=
0.04 a= 1.0 b= -11.5 q=
0,0000 0,2415
0,0417 1,1939
0,0833 1,1936
0,1250 1,1941
0,1667 1,1936
0,2083 1,1944
0,2500 1,1934
0,2917 1,1950
0,3333 1,1928
0,3750 1,1961
Здесь рассматривается второе сечение по 0X со значением х = 0.04.
Левый числовой столбец содержит значения отсчетов по оси 0Т, а правый числовой столбец содержит значения u(x=0.04,t) по этим временным отсчетам.
Тогда при построении графика по оси 0X будут идти временные отсчеты, а по оси 0Y -значения u(x=0.04,t) по этим временным отсчетам.
Построение идет аналогично в мастере диаграмм в EXCEL.
Так строим все графики.
В результате счета выдавались графики по границам, графики по временным сечени-
0.2
2 х=
0.2
ям, по сечениям по X, значения неизвестных с указанием координат узлов.
Рассмотрим результаты решения уравнения Колмогорова - Чепмена.
Исходные данные:
а = 1, 0,Ь = -11, 5, д = 0, 2, количество точек по X (горизонталь) = 5, шаг dx = 0, 25,
количество точек по по Т (вертикаль) = 5, шаг dt = 0, 25,
граница х = 0 - значение на границе = 1, 0, граница х = 1 - значение на границе = 2, 0, матожидание = 1, 0, среднее квадратическое отклонение = 1, 0.
Значения неизвестных в узлах: номер по t = 2, номер по х = 2, значение и = 1, 337,
номер по t = 2, номер по х = 3, значение и = 1, 622,
номер по t = 2, номер по х = 4, значение и = 1, 845,
номер по t = 3, номер по х = 2, значение и = 1, 336,
номер по t = 3, номер по х = 3, значение и = 1, 620,
номер по t = 3, номер по х = 4, значение и = 1, 842,
номер по t = 4, номер по х = 2, значение и = 1, 379,
номер по t = 4, номер по х = 3, значение и = 1, 706,
номер по t = 4, номер по х = 4, значение и = 1, 921.
Математическое ожидание и среднее квадратическое отклонение используются в программе для определения значений и(1:,х) на границах:
і = 0, і = 1.
Список использованной литературы:
1. Бендат Д. Ж., Пирсол А. Измерение и анализ случайных процессов. - М.: Мир, 1974.
2. Г. Корн, Т. Корн. Справочник по математике. - Издательство «Наука», Москва, 1973.
Статья рекомендована к публикации 11.04.08