Fedorenko-RP-Vvedenie-v-vychislitelnuyu-fiziku (810773), страница 65
Текст из файла (страница 65)
Итак, на данном этапе рассмотрения мы отказались от претензий правильно рассчитывать эволюцию 4/5 всех гармоник, т.е. в рассчитываемом явлении они не должны играть существенной роли. Интересующее иас решение не должно претерпеть существеннмх изменений, если мы ограничимся отрезком ряда Фурье й(1, х) = ',Г со к л!! з!и /!пх. в ! Если это не так, то ста точек для расчета недостаточно. Но и это еще не все: точность численного решения задачи зависит от шага по времени т .
Выше было введено условие, при котором разностный коэффициент Фурье с" воспроизводит правильное знзчение сок г!и' с точностью, зависящей как от /1, так и от т. Проверим условие т/сгпг~1, полагая к = ЛЧ5: кпг(д!/5)г И 0 5 т е т И 1//Уг /!г (Здесь мы считаем, что 0.5 м1, так как е оз 0.61 яи 1 — 0.5.) Таким образом, да!ке при не очень высоких требованиях к точности, мы пришли почти к тому же соотношению между т и /!, которое следует из требования устойчивости для явной схемы! Можно взять схему второго порядка точности по тз и -и 1и ! — 2и +и+! 1и,-2и +и и+! и и+! и+! и+! и и г 2 Зг 2 ° иг + Для с~в получаем выражение с"=св ~1 — 2 — зш — у <1+2 — зш и о / к .
ггиа)// ! . гвяа) ПРИБЛИЖЕННЫЕ МЕТОДЫ ВЫЧИСЛИТЕЛЬНОЙ ФИЗИКИ 1ч, и з!з При том же условии йяй«1 имеем Сравним функции (1 — 0.5х)/(1+0.5х) и е ", При х = 1 погрешность достигает 10%, при х=0.5 она порядка 1%. В этом случае мы можем претендовать на приличную точность при тхзпз < 0.5, т.е. прн т ц лз. (Это может породить неверное впечатление, будто бы схема второго порядка точности не имеет преимуществ. Конечно, имеет, но они сказываются в точности расчета более гладких компонент решения, соответствующих меньшим х. Не надо также упускать из вида, что как в точном решении, так и в разностном негладкая компонента решения быстро стремится к нулю с ростом и) В каких же задачах применение неявных схем дает существенный выигрыш по сравнению с явными? Это — задачи с нелинейной теплопроводносгъю не решением типа «тепловой фронт».
В решениях таких задач можно выделить три характерных области (см. рис. 33). 1. Зона перед фронтом тепловой волны (или «фон», по которому распространяется тепловая волна) характеризуется очень малыми температурами иф. Коэффициент теплопроводности и»ф в уравнении (2) так мал, что практически не происходит никаких перетоков теплоты. Это утверждение имеет смысл лишь относительно тех характерных времен Т, которые нас интересуют в данной задаче. Оно означает, что Ти~~/1.з«1, где Š— характерное для задачи расстояние по х.
2. За фронтом тепловой волны температура и+ очень велика и реализуется почти изотермический ре:ким: температура почти не меняется по х, но может меняться по 1. Опять-таки дело в том, что очень велик коэффициент теплопроводности и+, точнее (в безразмерных терминах) ти»Р/1Р~1. Из этого соотношения следует, что за малое с точки зрения характерных времен задачи время т«Т в изотермической зоне успевает выравняться температура. Пусть в начальных данных в изотермической зоне имеется какой-то профиль температуры. Разложим его в ряд Фурье, учитывая, что характерный масштаб по х есть Ты ! (0 х) ..1 ~ с е!тлл!ь (Для иллюстрации примем модель линейной задачи: и, = и»и „.) Тогда через время т решение будет »»» » и(т х) =с +"''с е!"' »!се о 319 нзлннвйнов гг«внвннв теплов»оводностн 9гц Даже для самой гладкой и наиболее медленной первой гармоники (т = ж1) временной множитель ехр ( — пви«+т//.з).з:1, т.е.
фактически и(т, х) = с, Возмущения и неровности, наложенные на изотермический профиль, мгновенно (с точки зрения времени Т) выровнялись. В этой зоне и~+ играет роль большого параметра, и решение определяется «квазистационарным» уравнением 1и"и,1„=0, или и"и„= с(1), т.е. тепловой поток понти постоянен по х, но, вообще говоря, может меняться по времени, 3. И наконец, есть еще переходная зона — зона тепловой волны, в которой и переходит от и«к из и профиль и(1, х) носит характер, близкий к автомоделъному. Передний фронт тепловой волны рассчитывается при условии Куранта ти«ьв 1Р.
Нельзя заранее разделить всю область переменных (г, х) на зти три зоны. Мы хотим их рассчитывать по единой схеме, не вводя разных формул в разных зонах. Здесь проявляется решающее преипущество неявной схемы, которая выдерживает сильное изменение критерия Куранта (безразмерной величины тив/й*) и не теряет устойчивости. Теперь уместно вспомнить, что мы только что скомпрометировалн расчет с очень большим «курантом», показав, что он не обеспечивает точности в передаче временной эволюции коэффициентов Фурье. Это так, но в изотермической зоне точный темп временной эволюции разных гармоник нам и не нужен. Важно только, чтобы разностная схема правнлъно передавала качественный характер— почти полное их исчезновение за время т, и это она обеспечивает. Вспомним еще раз жесткие системы уравнений: уравнение теплопроводности является жесткой системой в бесконечномерном пространстве.
Именно решение задач с описанным выше качественным поведением решения и было основным стимулом, приведшим к активному использованию неявных схем для уравнений теплопроводносги и к «изобретению» метода прогонки. Метод потоковой прогонки. При расчете решений, содержащих изотермический участок с очень большим коэффициентом теплопроводности, вычислители встретились с характерной трудностью, преодоление которой, в частности, привело к созданию специального варианта прогонки, названного потоковым.
В чем же было дело? На этом участке, как уже отмечалось, поток и«и„с(г) был почти постоянным по х, а величина с = О(1). Стало быть, и„= с/и"— величина очень малая н в некоторых ситуациях настолъко малая, что ее невозможно правнлъно вычислить по разностной формуле тиПа (и — ич 1)/й, пРЯБлиженные методы еьяислительной Физики 1ч. и ззо В самом деле, величины и в машинном представлении заменяются на й и„(1+ е ), где е — погрешность машинного представления чисел. Фактически ЭВМ вычисляет й — й и — и -- = - --+О/ — "~ Л Л ~Ч В некоторых случаях мы сталкиваемся с ситуацией, когда погрешность О(НБ/л) много больше основной величины (и — и~,)/й и ничего хорошего ожидать не приходится. Поясним сказанное несколько иначе. Конечное значение теплового потока с получается, так сказать, раскрытием неопределенности «с =* О.
м», причем И„О, иа м. Но на ЭВМ с конечным числом разрядов в представлении числа значение разностной производной (й — й,)/л не может стремиться непрерывно к нулю: оно либо нуль (при Й вЂ” й 1), либо не меньше ~ и~ е/л. Первый случай реализуется при совпадении и и и, со всеми машинными знаками, второй — при различии их хотя бы в последнем знаке мантиссы машинного числа. Если модуль ~с~ миеи„е//а, поток и" и либо нуль, либо много больше с. В этом источник трудностей. Его можно преодолеть переходом к расчету с двойной точностью, но можно поступить иначе, применяя метод потоковой прогонки. Его основу составляет представление уравнения теплопроводности (1) в виде си- стемы аи ап — = — +Ц, аР ах ди н — =П.
дх Запишем разностные уравнения на сетке (х ) и Б = О (х„= лзй): и"+ — и" П вЂ” П вЂ” '~+Я я=1 2 ...,М вЂ” 1 х 6 (11) и"+' — и" +' вР! в И„,РНЗ а — — П Рнп т =О, 1, ..., М вЂ” 1. Прогоночное соотношение имеет вид где Р, А — прогоночные коэффициенты, которые должны быть определены. первые прогоночные коэффициенты РБ, ЯБ определяются из левого краевого условия.
Пусть оно имеет вед и(г, О) = р(г). тогда РБ = О, /1Б = р(д„+,). предоставим читателю вывести формулы для Рр, Яп в случае, когда поставлено общее краевое условие и„+ пи= р. зш НЕЛИНЕЙНОЕ УРАВНЕНИЕ ТЕПЛОПРОВОДНОСТИ з Зи Рекуррентное соотношение получается посае исключения и"+' из второго уравнения (!1): и" +! — Р П +и — Я вЂ” (ЬIИ~,Р! )П,, =О, т.е.
и" ~~! = 1Р— ИИ +Пз)П .Р!!2+ Р!т. В этом выражении исключим П +,!2 через П +2!2 и и" ++',, используя первое уравнение (11)! ВР! !! Е Пт+Пз = Пт+ЗГ — —, Ит+1+ —, "тн! — от+ 1. Разрешая полученное соотношение относительно и,"„"+',, имеем и"+' =Р И А+! тт1 +3!2+ Р!' и" — йц +! +~!1 ) Это и есть формулы потоковой прогонки. Получив значения Рт „йм ! и разрешив правые краевые условия, т,е. определив им!!, обратную прогонку реализуем, вычисляя поочеРедно Пт цз, и„","„П„зж и т.д. СтандаРтный шаг имеет вид Пт-П2= Пт+иг, "т +, "т "Й * т-1 Рт — 1Пт-1!2+ ~~т' Механизм преодоления трудностей, связанных с конечной разрядностью машинных чисел, тот же, что был указан в э 5: переход от уравнения высокого порядка к системе уравнений первого порядка.
Возможен и другой способ расчета области с очень большим коэффициентом теплопроводности, позволяющий обойтись стандартной прогонкой. Нужно лишь скорректировать формулу для и: и(и) =1и при и( и', (и') при им и ). Значение и' зависит от разрядности машинных чисел. Вычислительный эксперимент показывает, что точное значение х1и) не сущест- 11 — 1ВЗЗ 322 пгивлкженные методы вычислительной физики [ч. и венно, важно только то, что это очень большая величина.