Решение систем дифференциальных уравнений методом Рунге-Куты 4 порядкаДонецк 1998 год РЕФЕРАТ Дифференциальные Уравнения, Метод Рунге-Кутта, РК-4, Концентрация, Метод Эйлера, Задача Коши, Ряд Тейлора, Паскаль, Реакция, Интервал, Коэффициенты Дифференциального Уравнения. Листов : 28 Таблиц : 2 Графиков : 4 Решить систему дифференциальных уравнений методом Рунге-Кутты 4 порядка, расчитать записимость концентрации веществ в зависимости от времени, проанализировать полученную зависимость, удостовериться в действенности метода. Содержание: Введение 1. Постановка задачи…………………………………6 2. Суть метода…………………………………………8 3. Выбор метода реализации программы……………14 4. Блок – схема………………………………………...15 5. Программа…………………………………………..17 6. Идентификация переменных………………………19 7. Результаты…………………………………………..20 8. Обсуждение результатов…………………………...21 9. Инструкция к программе…………………………...23 10. Заключение………………………………………….27 Литература Введение Обыкновенные дифференциальные уравнения (ОДУ) широко используются для математического моделирования процессов и явлений в различных областях науки и техники. Переходные процессы в радиотехнике, кинетика химических реакций, динамика биологических популяций, движение космических объектов, модели экономического развития исследуются с помощью ОДУ. В дифференциальное уравнение n-го порядка в качестве неизвестных величин входят функция y(x) и ее первые n производных по аргументу x Единственные решения выделяют с помощью дополнительных условий, которым должны удовлетворять искомые решения. В зависимости от вида таких условий рассматривают три типа задач, для которых доказано существование и единственность решений. Первый тип – это задачи Коши, или задачи с начальными условиями. Для таких задач кроме исходного уравнения (1.1) в некоторой точке xo должны быть заданы начальные условия, т.е. значения функции y(x) и ее производных y(x 0 )=y 0 ’ , y ’ (x 0 )=y 10 , ... , y (n-1) (x 0 )=y n-1 , 0 . Для системы ОДУ типа (1.2) начальные условия задаются в виде Количество условий должно совпадать с порядком n уравнения или системы. Если решение задачи определяется в интервале x є [ x 0 ,x k ] , то такие условия могут быть заданы как на границах, так и внутри интервала. Минимальный порядок ОДУ, для которых может быть сформулирована граничная задача, равен двум. Третий тип задач для ОДУ – это задачи на собственные значения. Такие задачи отличаются тем, что кроме искомых функций y(x) и их производных в уравнения входят дополнительно m неизвестных параметров l 1 , l 2 , ¼ , хm , которые называются собственными значениями . Для единственности решения на интервале [x0 , xk] необходимо задать m+n граничных условий . В качестве примера можно назвать задачи определения собственных частот , коэффициентов диссипации , структуры электромагнитных полей и механических напряжений в колебательных системах , задачи нахождения фазовых коэффициентов , коэффициентов затухания , распределения напряженностей полей волновых процессов и т . д . К численному решению ОДУ приходится обращаться , когда не удается построить аналитическое решение задачи через известные функции . Хотя для некоторых задач численные методы оказываются более эффективными даже при наличии аналитических решений . Большинство методов решения ОДУ основано на задаче Коши , алгоритмы и программы для которой рассматриваются в дальнейшем . 1. Постановка задачи Многие процессы химической технологии описываются СДУ - начиная от кинетических исследований и заканчивая химическими технологическими процессами . В основу математических способов описания процессов положены СДУ и СЛАУ . Эти уравнения описывают материальные и тепловые балансы объектов химической технологии , а так же структуры потоков технических веществ в этих аппаратах . Для получения , распределения технологических параметров во времени и в пространстве (в пределах объекта) , необходимо произвести СДУ методом , которых дал бы высокую точность решения при минималььных затратах времени на решение , потому что ЭВМ должна работать в режиме реального времени и успевать за ходом технологического процесса . Если время на решение задачи большое , то управляющее воздействие , выработанное на ЭВМ может привести к отрицательным воздействиям . Методов решения существует очень много . В данной работе будет рассмотрен метод решения СДУ методом Рунге-Кутта 4 порядка. Для удобства работы на ЭВМ, необходимо данную кинетическую схему преобразовать в удобный для работы на компьютере вид. Для этого необходимо кинетическую схему процесса представить в виде уравнений. При рассмотрении кинетической схемы процесса необходимо учитывать коэффициенты скоростей реакций. Но, так как процесс протекает при изотермических условиях, коэффициенты скоростей реакций можно считать за константы скоростей химической реакции. Из приведенной ниже схемы мы можем составить ряд дифференциальных уравнений, учитывающих изотермичность процесса.
Исходя из полученных результатов, можно сделать вывод, что расчет произведен верно, так как, исходя из полученных значений скоростей реакций можно сделать вывод, что соблюдается баланс скоростей химической реакции. Вещество А на протяжении всего процесса расходуется на образование веществ В и С. Концентрации вещества А в начальный момент времени расходуется быстрее, чем концентрации его же в конце процесса. Это обусловлено тем, что скорость химической реакции зависит от концентрации реагирующего вещества. Производная имеет знак «минус». Это говорит о том, что вещество расходуется. Следовательно, чем выше концентрация вещества, вступающего в процесс, тем выше скорость его реагирования с другими веществами. Вещества В и С образуются пропорционально, так как, исходя из кинетической схемы процесса и значений констант скоростей химической реакции, видно, что образование этих веществ и расходование этих веществ, одинаково. Производная имеет знак «плюс». Это говорит о том, что вещество образуется. График . 4 Это видно также и по результатам расчета, на протяжении всего времени исследования процесса концентрации и скорости веществ В и С одинаковы. В этом можно убедиться по виду графической зависимости концентрации веществ В и С от времени. Можно сказать, что процесс протекает в сторону увеличения концентрации веществ В и С и уменьшения концентрации вещества А. Процесс будет протекать до момента установления равновесия, но в данном случае равновесие не установлено, так как вещества продолжают расходоваться и образовываться. На протяжении всего процесса ни одно из образующихся веществ не поменяло знак производной. Это говорит о том, что процесс протекает в одну сторону. 9. Инструкция к программм е Итак, программа состоит из 3 основных процедур: 1) Init - процедура инициалиации, включающую в себя ввод данных; 2) Run - процедура вычисления и обработки результатов, включает в себя вызов двух вспомогательных процедур Difur, RK-4, Stroka, первая из которых отвечает за вычисление, а последняя - за вывод результатов в файл в табличном виде; 3) Done - процедура подготовки к выходу из программы; и трех вспомогательных: a) Difur - процедура вычисления производных (изменение концентрации веществ за единикцу времени ) b) RK-4 - используя значения производных, вычисленных процедурой Difur, вычисляет последущие концентрации веществ методом Рунге-Кутта c) Stroka - процедура вывода результата в файл в табличном виде Рассмотрим все эти процедуры поподробнее: Процедура INIT:
Данная процедура выполняет функцию инициализации программных данных, считывание данных из файла in.dat, создание, открытие на запись файла out.rez и запись в него шапки таблицы результатов. Процедура RUN:
| ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
PROCEDURE DONE; BEGIN CLOSE(f1); WRITELN(f2,'___________________'); CLOSE(f2); WRITELN; END; |
Процедура DIFUR:
PROCEDURE Difur; BEGIN dC[1]:=C[3]*k2+C[2]*k4-C[1]*k1-C[1]*k3; dC[2]:=C[1]*k3-C[2]*k4; dC[3]:=C[1]*k1-C[3]*k2; END; |
Процедура STROKA:
PROCEDURE STROKA; BEGIN WRITE(f2,'¦',x:4:1,'¦',c[1]:7:3,'¦',c[2]:7:3,'¦',c[3]:7:3,'¦'); WRITE(f2,sum:3:0,'¦',dc[1]:7:3,'¦',dc[2]:7:3,'¦',dc[3]:7:3,'¦'); WRITELN(f2); END; |
PROCEDURE RK_4; BEGIN Difur; FOR i:=1 TO n DO BEGIN r1[i]:=dC[i]; C[i]:=cPR[i]+r1[i]*(dX/2); END; Difur; FOR i:=1 TO n DO BEGIN r2[i]:=dC[i]; C[i]:=cPr[i]+r2[i]*(dX/2); END; Difur; FOR i:=1 TO n DO BEGIN r3[i]:=dC[i]; C[i]:=cPR[I]+r3[i]*dX; END; Difur; FOR i:=1 TO n DO r4[i]:=dC[i]; FOR i:=1 TO n DO rSR[i]:=((r1[i]+r2[i])*(r2[i]+r3[i])*(r3[i]+r4[i]))/6; END; |
Программа представляет собой 2 файла – файл с исходным текстом на языке Паскаль smith.pas и исполняемый модуль smith.exe скомпилированный компилятором TNT Pascal 3.25 фирмы Layer`s Ins. Исполняемый модуль программы предназначен для запуска в операционных системах: MS Dos, Windows95, Windows NT, OS/2, а также в X-windows под Linux ( при наличии эмулятора ) Для нормальной работы программе необходимо 640 кb «нижней» памяти и 20 kb дискового пространства.
Согласитесь – требования минимальные, учитывая то, что сама программа абсолютно не требовательна к процессору. В процессе работы программа считывает данные из файла in.dat и записывает результаты работы в файл out.rez в табличном виде.
Исходный файл программма открывает стандартными средствами ОС, не проверяя его наличие перед работой, поэтому, если данный файл не будет доступен в каталоге, в котором расположена программа, компилятор выдаст сообщение об ошибке. Если Вы после запуска программы увидели что-то типа «Runtime error 202 at 0000:0A86» - это всего лишь значит, что программа не смогла найти файл с исходными данными в текущем каталоге. Если Вы забыли поместить его туда, скопируйте этот файл в каталог с программой и запустите исполняемый модуль еще раз. Если данный файл у Вас отсутствует, Вам прийдется сделать его самому. Для этого в любом текстовом редакторе наберите 3 выделенных строчки и сохраните созданный файл с именем in.dat 100 0 0 0.2 0.1 0.2 0.1 0 10 0.5 3 0.05 0 Создав файл и скопировав его к исполняемому модулю программы, запустите исполняемый модуль еще раз. В процессе работы программа будет выдавать сообщения об успешном окончании каждого блока. Если все прошло нормально, то на экране своего компьютера Вы увидите следуще сообщения: Step 1: Read data from file : in.dat - done. Step 2: Write header to file : out.rez - done. Step 3: Calculating data and writting results to file : out.rez - done. Step 4: Close all files and exiting... Первый шаг (step1) сообщает, что данные из файла in.dat были успешно прочитаны Второй – о том что программа успешно создала выходной файл out.rez и записала в него шапку таблицы с данными В третьем сообщении сказано, что данные успешно посчитаны и записаны в выходной файл out.rez Четвертое сообщение сообщает об окончании вычислений и завершении программы. После того, как программа отработает, Вы сможете познакомится с результатами, которые были вычислены и помещены в файл результатов out.rez. Просмотрев его любой программой просмотра текстовых файлов или вывев его на печать, вы получите таблицу c результатами . 10. З аключение. В результате выполнения расчета получена зависимость изменения концентрации вещества во времени. Из расчета следует, что на протяжении всего процесса вещество А расходовалось на образование В и С. Процесс не достиг конечного состояния (не достиг равновесия) Максимум концентрации вещества наблюдался при следующих значениях времени: при начальном значении времени max соответствовал веществу А; при значении времени, равном 10 часам, max соответствовал веществам B и С, однако, это не является максимумом концентрации веществ в процессе вообще, так как вещества B и С продолжают образовываться; В ходе выполнения работы был произведен расчет системы дифференциальных уравнений методом Рунге-Кутты четвертого порядка, произведен расчет кинетической схемы процесса при изотермических условиях при данных значениях концентраций и констант скоростей.