метода Волощенко, страница 5

Погрешность численного решения дифференциальных уравнений (ДУ) сильно зависит от шага интегрирования Н и скорости изменения переменных. Для каждого численного метода решения ДУ н конкретной системы ДУ существует некоторое оптимальное значение шага Ншс при Н > Нсп затраты машинного времени на решение уменьшаются, но возрастает погрешность решения; при Н < Нш погрешность практически не изменяется, но возрастают затраты времени, поэтому обычно выбирают Н < Нш. Значение Нш можно определить экспериментально путем решения задачи для разных значений Н. Обычно начинают с достаточно малого значения Н и увеличивают его в 1,б — б раз для каждого решения до тех пор, пока не начнет возрастать погрешность.

Погрешность решения можно определить как модуль разности значений переменных для двух последовательных решений и одинаковых моментов времени. Различают три типа погрешностей: текущая относительная погрешность (ЕМТ) вычисления перемен- ной номер т в момент времени С, определяемая из соотношенияс ,я+1 ут„с ут, с (3.1) ут, с где р — номер решения; ус,, у с — значения переменвои р+ 1 номер т в момент времени с для решений номер р и р + 1 соответственно с ° погрешность решения (ЕТ) в момент времени с, определяемая как максимальное значение из погрешностей для всех переменных, т.е.

(3. 2) максимальная погрешность решения (ЕЕ), определяемая как максимальное значение из погрешностей решения ЕТ для всех шагов по времени, то есть (З.З) ЕЛ = снах ЕТ л где и — номер шага по времени. 3.2. Решение дифференциальных уравнений на ЭВМ При анализе переходных процессов в электрических цепях обычно решается задача Коши для системы дифференциальных уравнений: ус=(С(С УС у21 "~ут)' У2 с2 ~С* У1 * У2 ~ '" Ут)" (3. 4) Ут =)~(С, У,, У2, ...,Ут); с начальными условиями У1( О) = У10: У2(го) = У20 ' -. Уш(го) = У,но где у1, у2 ., у — токи или напряжения в компонентах электрической цепи, à — время. Систему (3.4) можно записать в векторной форме: У=У(у, г) (3. 5) где У вЂ” вектор производных, à — вектор правых частей системы (3.4).

Существуют различные методы Еи 7 численного решения задачи Коши. В Н и — 1 основе их лежит представление производной в виде алгебраического выра- Г жения и переход от системы диффе- 2 Иг Ил+1 ренциальыых уравнений вида (3.4) к системе конечно-разыостных алгебраических уравнений для дискретных моментов времени, отстоящих друг от друга на величину Н, называемую шагом интегрирования (рис. 3.1).

Систему конечно-разностных уравнений можно записать в векторыой форме: 0 1 Рис. 3.1 (3. 6) У =а +Ч' +1= оул Е~У Ул-1 -' у —.) ° либо в форме +1= по у + Ч2 у +1 ° ун ° Ул — 1 -* Ул-1) (3.7) где ӄŠ— искомые значения переменных в узле и + 1; У„, У„Е, ..., У„г — известные значения перемеыных в узлах и, и — 1, ..., и — Е соответственно. Формула (3.6) называется явной, так как искомые переменные (у„1) входят только в левую часть уравнений, поэтому у„+ 1 моясно непосредственно определить из (3.6).

Формула (3.7) называется неявной, так как у„+ 1 входит как в левую, так и в правую части формулы, поэтому выражение (3.7) представляет собой систему уравнений (в общем случае нелинейную), которую необходимо решать тем или иным итерационным методом относительно у„ Так как значения переменных для момента времени г = г известны из начальных условий то, используя соотношение (3.6) или (3.7), можно последовательно определить значения переменных для моментов времени 11, 12,, гл . Задача расчега переходных процессов обычно формулируется следующим образом: построить зависимости переменных у, У2, ..., у от времени на заданном интервале. Обычно шаг Н выбирают так, чтобы получить минимальную погрешность. Поэтому количество узлов на интервале [О, гл ) может быть очень большим.

В то же время для анализа получаемых зависимостей обычно достаточно 30 †1 узлов, поэтому в массивы выходных данных целесообразно заносить значения переменных не для каждого узла по времени ц а только для узлов, отстоящих друг от друга на величину РЬ. Величину ВЬ обычно задаег пользователь. Отсюда следует, что алгоритм должен содержать два циклических процесса: внутренний — для расчета на интервале П1 с шагом Н и внешыий — - для расчета ыа интервале [О, Тк) с шагом П1. Обычно внутренний цикл оформляют в виде процедуры. При численном решении дифференциальных уравнений количество шагов по времени может быть очень большим, поэтому на этапе анализа задания все расчетные формулы необходимо преобразовать к виду, содержащему миыимальное количество арифметических операций в циклах.

3 лечение. Длл повышения вФфективности алгоритма рекомендуется расчет на интервале ПЬ проводить с переменным шагом Н. Для автоматического изменения шага Н можно испольэовать сгмдуюший алгоритм. На каждом шаге проводим два решения: с шагом Е1 и шагом Н/2. Затем проверяем, если погрешность меньше ег, то увеличиваем Н в два раза. если же погрешность больше е2, то уменьшаем Е1 ые числа. причем е1 меньше е2).

в лва раза (адесь прелполагается что е1 и е2 — маа 3.3. Численные методы решения дифференциальных ураинений Болыпинство численных методов решения дифференциальных уравнений базируется ыа двух основных подходах: разложении в ряд Тейлора (их называют методами Рунге — Кутта) и полиыомиальыой аппроксимации (методы численного интегрирования). Методы Рунге — Кутта. Существует несколько методов Рунге— Кутта разного порядка точности, при этом все они являются одно- шаговыми и явными Пусть имеем систему из т ДУ первого порядка, записанных в нормальной Форме (3.5).

с начальными условиями у(О) = у В методе Рунге — Кутта четвертого порядка значения переменных в узле определяются из соотношения У) и+ 1 = Уй и + ~К) о + 2К) 1 + 2К) 2 + К) 3 ) б-, (3.8 ,Н где к. о, к 1, ..., к. 3 — вектоРы коэффициентов, Рассчитывае мые по формулам; К.

=РН у„с ( п,У1,п'Уз,л~ -*Ут,п)' у,1 ("л 2'У1п ~2~10'Узп+~ЗК20' — * Утл 2 тс) * (3.9) й 2 ~ ~1п ~2 ' У1,п + ~2 К1,1 г -"™ Ут,п ~ 32 Кт,1) ) у,з (1л У1 л+ 2 1,2~ '"' 'Ут,л+ 2 т,2)» где й2 — — Н/2. Из соотношений (3.8) и (3.9) видно, что можно последовательно вычислять коэффициенты К- о, ..., К. в порядке их запис 10'"' )3 иси в ( .9) и одновременно накапливать значения переменных на оче- 13 редном шаге Н в соответствии с соотношением (3.8). Все зти вычисления целесообразно оформить в виде процедуры решения уравнений на интервале от Т до Т + ПЬ . Т ак как козффиЦиенты К. С, ..., К 3 вычислЯютсЯ по оДним и тем же формулам для всех уравнений, целесообразно расчет этих коэффициентов и накопление их суммы на каждом шаге по времени организовать в виде циклического процесса.

Для упрощения процедуры и повышения степени общности алгоритма целесообразно вычисление функций (правых частей уравнений) вынести в отдельную процедуру. Один из возможных вариантов схемы алгоритма процедуры решения ДУ методом Рунге — Кутта четвертого порядка (назовем ее МВК) приведен на рис. 3.2. Обращение к процедуре имеет зи ВК(М,Н,Т,ВЬ,Е, г' ), где М вЂ” количество дифференциальных д уравнений, Н вЂ” шаг интегрирования, Т вЂ” начальный момент вре- 30 мени, ПЬ вЂ” интервал интегрирования, 2 — массив переменных для момента времени Т; У вЂ” массив значений переменных для момент времени Т+ ПЬ.

Рис. 3.2 В алгоритме используется вспомогательный массив Я для на. копления суммы коэффициентов ке, ...„кз, а также переменные Х и Х1 для определения текущего времени. Вычисление правых частей уравнений осуществляется с помощью процедуры ВРР(Х1,У,Р), 31 где Х1 — текущее время; У вЂ” массив переменных; Р— массив значений правых частей уравнений. В каждом цикле по времени (рис. 3.2) сначала накапливается сумма коэффициентов сгс, ..., Нз в массиве Б, а затем вычисляются аначения переменных (массив У).