Пасконов В.М., Полежаев В.И., Чудов Л.А. Численное моделирование процессов тепло- и массообмена (1185910), страница 30
Текст из файла (страница 30)
ниже и. 6.8.3). Иэ сказанного следует, что требования к вычислительным методам для решения уравнений Навье — Стокса должны различаться в аависимости от рассматриваемого диапазона чисел Рейнольдса и тех целей, которые ставятся при численном моделировании. Общие требования к вычислительным методам можно сформулировать следующим образом: 1) Вычислительная устойчивость. 2) Точность расчета основных характеристик, приемлемая для соответствующих приложений. -3) Экономичность; минимальный объем оперативной памяти; простота реализации. Первое требование заключается в том, чтобы весь вычислительный процесс в целом был устойчив. Оно относится как к самой разностной схеме, так и к методу решения соответствующей системы алгебраических уравнений, Основные определения были даны выше в гл.
2 — 4. Для разностных схем, аппроксимирующих уравнения Навье — Стокса, причин неустойчивости, однако, больше, чем для простых модельных уравнений, рассмотренных в упомянутых главах, причем в ряде случаев явления вычислительной неустойчивости трудно отличить от возможного сложного поведения решений. Второе требование означает необходимость высокой пространственно-временной разрешимости, которой можно 173 в принципе достигнуть, либо применяя схемы пе слишком высокого порядка точности, реализуемые на подробных пространственно-временных сетках, либо существенно повышая порядок точности схем.
Для уравнений Навье— Стокса особенно важным является построение разностных схем, аппроксимирующих общие нестационарные уравнения (к позволяющих в частном случае'определять стационарные решения, если таковые существуют). При этом практика показывает, что для расчета весьма широкого класса течений достаточно использования схем первого порядка точности по времени. В отличие,от течений не- вязкой жидкости, при этом характерны более высокие требования к пространственной аппроксимации решения (пограничные слои, основные и вторичные течения и т. д.).
Наиболее удобными являются раэностные схемы второго порядка точности по пространственной координате на неравномерной сетке, сгущающейся в зоне больших градиентов. 'Третье требование на самом деле может состоять иэ двух (или даже трех) требований: минимального числа операций на временном слое, миимального объема оперативной памяти ЭВМ и минимальных затрат труда программиста на реалиэацию программы.
Перечисленные требования в известной мере условны, так как значение каждого иэ них зависит от ряда дополнительных факторов, таких, например, как режим течения по числу Рейнольдса, тип ЭВМ, квалификация исполнителя, ограничения на время для получения результата, серийность расчетов и т. д. Эти требования, кроме того, противоречивы, так как одновременное и полное их выполнение практически невозможно, что требует компромиссных решений..
'Ф э 6.2. Раэностные схемы для уравнений Навье — Стокса. Предварительное рассмотрение Ф 6.2Л. Простейшая разностная схема для двумерных уравнений. Для того чтобы скорее подвести читателя к вопросу о конструировании конечно-разностных схем для уравнений Навье — Стокса, рассмотрим сначала одну . иэ простейших схем численного интегрирования. Будем в качестве исходной испольэовать систему двумерных уравнений Навье — Стокса для однородной иэотермической жидкости в переменных вихрь, функция тока, со- йтй стоящую из иестационарного уравнения вихря и стацио. нарного уравнения Пуассона для функции тока: дв дф де дф да 1 ( д в д и ') у дд др дз дх ду Ке.~ дзз дрз / д'1р д'р — + — = ю~ (6.2.2) дзз ду Рассматривается течение в замкнутой квадратной области при граничных условиях на твердой границе ф=О, —,.„=О.
дф (6.2.3) Начальные условия заданы в виде ~р(х, у, 0) = ф'(х, у), ю(х, у, 0) = в'(х, у). Для аппроксимации дифференциальных уравнений разпостными введем пространственно-временную сетку с координатами х,=й, у,=)1,. т„=ит, где Ь, 1 — шаги сетки по координатам х, у соответствен- ной т — шаг по времени; 1=0, 1, ..., М вЂ” 1; 1=0, 1,,... ...,М вЂ” 1;л=0,1,...,К. Введем следующее обозначение: — <р (й, у1, пт) = рп;.
Производные по пространственным переменным будем аппроксимировать центральными разностями, например, дч т1+1д %~-1д д Е Ч~+1д ть)+ Чр-к) 3 — 2 дх . 2й ' дзз й~ де Фар+1 'Рйд-г д <р %;,)+~ — 7~я+%;,) — 2 ду 2) ' дуз ра Производную по времени заменим разностным отношением «вперед» в виде де %(,д %ад дю Запишем, используя указанные аппроксимации, следующую явную схему для уравнения вихря (6.2 1): <у+1 и щ~ ~~а Я ей Е'д З1ь) р Е1 +ь) Е(-ь) а Ейд+г Ей)-1 2й М 22 й (6.2.4) 175 Здесь .»го п "».»- 2) .
оМ=— д»г »со= о ° да~ Запишем вторую производную от функции тока граничном узле следующим образом: 1»,о — 2$»»+ Ф»,о 2Ф»»д»,о+ Фьо ьо ьо ьо * = — + Рассмотрим раз постную запись граничного дфдп О в виде 3»о» о — 4»р» +»гьо = О. в при- условия Подставляя выражение»васо из последней формулы в предыдущую, получим, используя второе граничное условие »р»,, = О, выражение для вихря на границе в виде о»,",+' = 2»Р",/Ьо. (6.2.5) Значения поля вихря во всей области й в соответствии со схемой (6.2.4) и граничным условием (6.2.5) могут определяться различными способами: вдоль линий о 176 По этой схеме по известным в момент времени г, значениям полей функции тока (скорости) и вихря $»,„ »о»,, внутри расчетной области П, включая ее границу, можно определить значения вихря в области Й, исключая ее границу, в следующий момент времени г +, — — ~„ + т.
Связи, определяемые схемой, имеют локальный характер, так как для определения величины о»»+ требуется знать значения вихря на слое п в пяти точках: «»»лч»о»+,,;, »о» и;» »о»,»+г» о»»,»-» При определении вихря с помощью уравнения (6.2.4) требуется использовать те или иные условии для вихря на границе. Заметим, что условиями задачи вихрь на границе не задан, а заданы граничные условия для функции тока (которые, вообще говоря, относятся ко всей системе (6.2Л), (6.2.2)). Граничные условия для вихря можно получить, например, из уравнения для функции тока, считая его справедливым вплоть до границы; тогда получим, например, для границы р =сопзФ: ° совзФ, вдоль линий 7' сопела или последовательно по отдельным участкам, начиная от той или иной грашщы области, что представляет определенные преимущества при реализации алгоритма в виде программы для ЭВМ.
Перейдем теперь к решению уравнения Пуассона для функции тока (6.2.2). В отличие от уравнения для вихря, это уравнение стационарно. Это значит, что для получения решения системы (6.2.1), (6.2.2) на одном временном слое нужно решить стациопарное уравнение (6.2.2), где Я+1 правая часть — вихрь эв1Л вЂ” определена ранее. Для этого мы применим простейший явный итерационный метод (см.
$2.5). Его можно сформулировать по аналогии с решением нестационарного уравнения, если ввести фиктив. пое время а следующим образом: — = Ав(в — 1э. двг до Обозначая через г индекс внутреннего итерационного цикла, запишем схему для решения этого уравнения на временном слое п + 4 в виде 7и+1 в+1 ри в 1,'в 0 йа+1,в и Ра+1в+ в(1и+1,в+1 „Ри+1,в 2Ри+1в+ в(а+1в+1 1+1,1 вя 1-1,1 + 1Л 1-1 1 В 1Л-1 ь~ ь~ — е11~,+'. (6,2,6) Аппроксимация оператора Лапласа здесь строится таким образом, чтобы искомое значениевР1,1 ' можно было опи+1,в+1 ределить, не прибегая к решению системы алгебраических уравнений.
После преобразований (6.2.6) можно представить в виде „ри+1,в+1 и+1,в „~ 1 ( и+1,в,и+1,в+1 + и+1,в + и+1,в+1 + ив и+1,в) „и+1,в~ где ав — итерационный параметр, определяемый через сеточные параметры о, Ь (вхв= 4о/й*). Граничным условием при расчете по формуле (6.2.7) является условие вр- О, задаваемое на границе области (другое иэ граничных условий дв(/дп= О уже было использовано при получении граничного условии для вихря (6.2.5)). 1о в.
м, паокоиов и ар. 177 Расчет поля функции тока по формуле (6.2.7) проводится до получения стационарного решения. Зто значит, что внутренний итерационный цикл с параметром г должен заканчиваться при определенном условии, которое характеризует достижение стационарного режима. При плавном изменении ф,з в процессе итераций можно использовать самое простое из этих условий, состоящее в том, что разность значений функции тока в двух соседних итерациях г и з + 1 не превосходит некоторой заданной величины шах ~ ф+) ' — ф ) ~ ( е. (6.2.8) При выполнении условия (6.2.8). расчет уравнения Пуассона по формуле (6.2.7) прекращается, и мы имеем поле вихря и поле функции тока, удовлетворяющие разностным аналогам уравнений (6.2.1), (6.2.2) на временном слое и+1.
Для получения решенйя в следующий момент рассмотренная выше процедура повторяется, с той лишь разницей, что в качестве начальных значений теперь используются найденные величины полей ф,;, а;,; . . Повторим в краткой форме еще раз всю последовательность расчета полей вихря и функции тока при переходе от слоя к слою. 1) Значения а~;, ~р,".,) предполагаются известными внутри расчетной области 1). 2) Значении вихря на границе расчетной области О определяются по формуле (6.2.5).
3) Определяется поле вихри внутри области П на слое и+1 по формуле (6.2.4) при граничных условиях (6.2.5). 4) Определяется поле фу)акции тока при граничном условии ф = 0 и при использовании найденных значений вихря ю,"т' путем итераций по формуле (6.2.7) до тех пор, пока не будет выполнено условие (6.2.8). Злементарный анализ устойчивости (гл.