метода Волощенко (551678), страница 6
Текст из файла (страница 6)
ЭФфективность этого алгоритма можно существенно повысить, если предусмотреть возможность численного решения с переменным шагом интегрирования. Методы численного интегрирования. Методы численного интегрирования основаны на полиномнальной аппроксимации. Если при определении значений переменных в узле номер и + 1 используются ранее вычисленные значения переменных в нескольких предыдущих узлах, то метод является многошаговым.
Порядок метода зависит от количества учитываемых предыдущих узлов. Расчетные формулы для основных методов численного интегрирования приведены в табл, 3.1. Заметим, что в случае системы дифференциальных уравнений расчетные формулы, приведенные в таблице, необходимо применять к каждому уравнению системы. Первые пять методов относятся к классу неявных методов, так как правая часть расчетных формул зависит от искомых переменных У„+ 1. Поэтому при использовании этих методов необходимо каким-либо итерационным методом решать систему разностных уравнений. Методы Адамса — Башфорта являются явными„значения переменных У„+1 могут быть вычислены по формулам (п.б и п.7 табл.
3.1)непосредственно. Метод прогноза и коррекции строится путем комбинирования явных и неявных методов (п.8 табл. 3.1.). При этом на каждом шаге по времени расчет проводится в два атапа: 1-й этап (прогноз): по явной Формуле вычисляются прогнозируемые значения переменных УУ в узле и + 1; и+1 2-й этап (коррекция): используя найденные значения Ус' и+1 вычисляются значения функций ~ (Ус„1, си г 1) и определяются новые скорректированные значения переменных У и+ 1 Достоинством методов прогноза и коррекции является то. что они позволяют отказаться от итерационного процесса решения неявных конечно-разностных уравнений на каждом шаге по времени.
32 Таблисса 3. с Поряд- ковый Расчетная формула Метод номе Трапеций (порядок 2) Гира (порядок 2) 18 9 2 б — и 2 И у -.— у + — у + — нЛу Гира (порядок 3) 5 8 НУ(у,С )+=НЛУи,уи) и 12 ссз-1' пз-1 У2 1 — — Н у(у,с ) и-!' и-1 23 1б !2 5 —,-НЛУи 2 си 2) В табл. 3.1: л — номер шага по времени; Н вЂ” величина шага по времени; ~(у, 1) — правая часть дифференциального уравнения, записанного в нормальной форме. 33 Адамса- Маултоиа (порядок 3) Адамса— Маултона (порядок 4) Адамса " Башфорта по ядок2) Адамса— Башфорта (порядок 3) Прогноза и коррекции (порядок 4) ! у =у + — Н1Яу С )+С(уп Си)1 у — — у + — Ну(у,с ) 4 ! 2 и+! 3'сс 3'и-! 3 "+!' с!+1 у = у +- Н (9/'(У,С ] у! ЧЛуи.Си) 1 и+1 сс 24 и+1' и+1 -5ЛУи с ) Луи 2'и 2)) 3 ! „„.. нЛу,с ) -нЛу,с ) пч-1 п 2 "' п 2 и-1' и — ! Прогноз У =Уп" (- ЛУп сп) У(уп !' и-1) и -37Лу,с 1 9ЛУ с И Н у „из- 24 и+ и+ и — !' и-1) и-2' и — 2 схемы.
Рвс. 3.4 Рис. 3.3 З.э.1. Анализ задания ния (3.10) в виде с(эг — = а1 Уг, (3.11) ~~32 д = аз Е(') аз у2 а2 у1 оиС 1 1 д1 С (3 10) 1. В. 1 =' =- — ' Е(1) — — ( — — и ,р г, гс г с Для применения многошаговых методов надо знать значения переменных в узлах номер и, а + 1, а + 2 и тэк далее, а из начальных условий известны значения переменных только в одном узле. Зту проблему обычно решают путем применения на первых шагах того или иного одношагового метода, например одного из методов Рунге — Кутта. Заметим, что в этом случае порядок одношагового метода должен быть не ниже порядка основного многошагового метода, либо расчет начальных узлов надо проводить с более мелким шагом Н.
Алгоритмы многошаговых методов также целесообразно оформлять в виде процедур расчета на заданном интервале. времени ())1 ). Для хранения значений переменных в предыдущих узлах (не во всех, а только тех, которые используются в расчетных формулах) необходимо предусматривать соответствующие массивы переменных (и сдвиг значений на один шаг после расчета очередных значений У„г). Если метод решения неявный и используется итерационный метод решения системы алгебраических уравнений, то должны быть предусмотрены средства контроля количества итерационных циклов и формирования соответствующего признака окончания вычислительного процесса. Способ вычисления функций, соответствующих правым частям дифференциальных уравнений, зависит от используемого итерационного метода.
Обычно выражения для правых частей уравнений выносят в отдельную подпрограмму, чтобы можно было легко вносить изменения. 3.4. Пример разработки программы расчета переходных процессов Задание: "Разработать алгоритм и программу расчета переходных процессов в схеме, приведенной на рис.
3.3. Исследовать влияние шага интегрирования Н на погрешность расчета переходного процесса". Математическая модель устройства имеет внд где и — напряжение на конденсаторе; (ь — ток через индуктивс ность; ЕЯ вЂ” входной сигнал; В, Е, С значения компонентов Закон изменения входного сигнала задан в виде Функции, приведенной на рис. 3.4, где обозначено: ТХ вЂ” длительность задержки входного сигнала; ТГ1 — длительность переднего фронта импульса; Т1 — длительность импульса; ТУ2 — длительность заднего фронта импульса; А — амплитуда входного сигнала.
Для решения дифференциальных уравнений использовать метод Гира второго порядка. Преобразование математических соотношений. Для построения эффективного алгоритма необходимо преобразовать математические формулы таким образом, чтобы количество математических операций в цикле по времени было минимальным. Введем новые переменные: д = ас, уз — — (с и запишем уравне- где а1 —— 1/С; аз — — 1/Е; аз —— В/Л . Очевидно, что эти коэффициенты можно вычислить один раз до цикла по времени.
Обозначим правые части уравнений Р = а1У1 Рй = аз Е(1) аз Уз ай У1 тогда систему уравнений (3.1Ц можно записать в виде ИУ1 — =Р ~й с(уг — =Р с(1 3 Расчетная формула в векторной форме для метода Гира имеет вид 4 1 2 Уп+1= 3 Уп 3 Уп — 1+ 3 1(Уз+1 ~ и+1) где и — номер шага по времени, Н вЂ” расстояние между точками с номерами и и и + 1; ) (У„1, 1 ~ — значение правой части дифференциального уравнения в момент времени („т 1. В данном случае мы имеем систему из двух уравнений поэтому соотношение (3.13) надо записать для каждого из уравнений. Обозначим 3' 3 4 1 2 3' 3 3' 3 3 Тогда для нашей системы уравнений получим У, „+,-Ь,У,„-Ь,У, и,+Ь,Р, п+,, (3.13) У2, и+ 1 Ь1 У2,п 2 Уй, п — 1+ ЬЗ Рй, и+ 1 * где У1 и+ 1, Уй и+ 1 — искомые значения У1 и Уй в точке и + 1; Р „, Рй и+ — значения правых частей уравнений в точке и + 1.
при преобразованиях могут быть допущены ошибки, что увеличивает сроки отладки программы; при изменении исходной системы уравнений необходимо ааново выпол- нять все преобразования„ Замечание. Нетрудно наметить, что искомые значения у( и+ 1, уй и+ 1 входят в левые и в правые части соотношений (3.13), что не позволяет вычислить их значения непосредственно по Формулам (3.13).
Конечно можно раскрыть выражения для правых частей и записать соотношения (3.13) в явном виде относительно р1 и + 1 и рй, + 1. Однако обычно так не делают по следующим причинам: часто вообще невоаможно выразить искомые переменные в явном виде, на- пример, если правые части содержат нелинейные члены. Поэтому обычно в явном виде искомые переменные не выражают, а решают систему уравнений (3.13) тем или иным итерационным методом.
Для повышения степени общности и наглядности алгоритма будем хранить переменные в массивах: Х(Ц=у( „,,; ХЬ((Ц=у „; ХЬ(1(Цму1 „; Р(Ц=Р; Х(2) = Уй, и+1: Хг((2) Уз, и ' ХьчЦ2) = Уй и 1 ( Р(2) = Р3 . Тогда расчетные формулы примут вид Х(Ц = Ь . ХЬ)(Ц вЂ” Ь ЪЪ)1(Ц вЂ” Ь . Р(Ц; (3.14) Х(2) = Ь1 . ХМ(2) — Ьй . ХМ1(2) — Ьз Р(2), где для нашей системы уравнений выражения для Р(1) и Р(2) будут иметь вид Р(Ц = а( Х(Ц; Р(2) = ай Е(1) — а3 Х(2) — ай Х(Ц . (3.16) Здесь предполагается, что в дальнейшем вычисление входного сигнала Е(1) будет реализовано в виде подпрограммы функции. Первые два члена в правой части системы (3.14) не зависят от Х, поэтому их можно вычислять один раз до итерационного цикла. Окончательно запишем систему (3.14) в виде Х(К) = Е(К) — Ь3 Р(К); К = 1, 2, (3.16) где Е(К) = Ь1 ХЬ((К) — Ьй ХЬ(1(К) Замечание.
Важным свойством соотношений (3.15) является то, что такая запись расчетных бюрмул позволяет легко иамеиить алгоритм при модификации исходной системы уравнений. при иаменении правых частей уравнений надо внести соответствующие ивменения в выражения для злемеитов массива р (выражения (3.15)).
а при увеличении количества дифференциальных уравнений достаточно расширитЬ массивы и авписатв правые части дополиительных уравнений в соответствующие влементы массива у. Для расчета переменных в узле и + 1 по формулам (3.14) необходимо знать значения переменных в двух предыдущих узлах по Т. Так как из начальных условий известны значения переменных 3.4.2.
Разработка алгоритмов Рнс. 3.б только для момента времени Т = О, поэтому воспользуемся методом Рунге — Кутта четвертого порядка для определения значений переменных на первом шаге, т.е. в момент времени Т = О + Н. Метод Гира является неявным, поэтому на каждом шаге систему уравнений надо решать итерационным методом. Для решения уравнений (3.14) воспользуемся методом простой итерации. Входные и выходные данные. В данной задаче входными данными являются: значения компонентов схемы: С, 1,, В„.
параметры входного сигнала: ТЕ, ТР1„'П, 'ГР2, А; начальные условия 1)Ы, Тг(; ° параметры вычислительного процесса: НХ вЂ” начальный шаг интегрирования, ТК вЂ” конечное время; РЬ вЂ” шаг вывода на печать; ЕРБ — погрешность расчета. Чтобы можно было проанализировать характер переходного процесса и зависимость погрешности решения от времени, будем выводить в виде таблицы (с шагом РЬ) следующие данные: Т— текущее время", 1)С вЂ” значение напряжения на емкости; ТР— ток через индуктивность и текущую погрешность ЕТ. Для этого зарезервируем соответствующие массивы.
Схема иерархии для данной задачи очевидна. На первом уровне необходимо реализовать следующие функции. 1. Ввод исходных данных. 2. Вычисление постоянных коэффициентов. 3. Решение уравнений методом Рунге — Кутта для момента времени О + Н. 4. Решение уравнений на интервале О...Тк. 6. Вывод результатов. На втором уровне функцию № 4 можно разбить на три функции: 4.1) вычисление входного сигнала; 4.2) решение уравнений методом Гира на интервале 1И.; 4.3) формирование массивов результатов. Алгоритмы метода Рунге — Кутта и решение методом Гира на интервале РР оформим в виде процедур„а блок вычисления входного сигнала в виде функции.