3.1. Общий вид системы (Математическое моделирование задач газодинамики и тепло-массообмена. Молчанов А.М. 2013)

Глава 3. Численные методы решения уравненийдвижения вязкой жидкости.Системадифференциальныхописывающаядвижениевязкойуравненийвньютоновскойчастныхпроизводных,жидкостиназываетсяуравнениями Навье-Стокса.Существует большое количество численных методов решения уравненийНавье-Стокса. Автор не ставил задачей дать обзор и сравнительное описаниеэтих методов.В этой главе представлены:- некоторые особенности, присущие уравнениям Навье-Стокса с точкизрения их численного решения;- основные требования к численному решению уравнений Навье-Стокса.Построен численный метод, который был практически реализован в видепрограмм на языках FORTRAN-90 и C++ (программа UNIVERSE-CFD).1.

Общий вид системы1.1. Особенности численного решенияВ главе 1 были представлена система основных дифференциальныеуравнений динамики вязкой жидкости в различных системах координат. Онавключает уравнение неразрывности, уравнения количества движения,уравнение энергии.По своей форме уравнения количества движения и уравнение энергииочень похожи на дифференциальное уравнение теплопроводности систочником тепла.Это сходство наводит на мысль, что для решения этих уравнений можноиспользовать такой же метод, как для решения дифференциальногоуравнения теплопроводности. Причем, каждое из этих уравнений можнорешать отдельно, последовательно друг за другом.Действительно, если задано поле давления, такой подход дает хорошиерезультаты.Вэтомслучаечлен,содержащийградиентдавления,рассматривается как заданный источник.А что же делать, если поле давления не задано и его требуется определитьна основе решения системы уравнений?Первое, что приходит в голову, это использование итеративного процесса.Например:1) задать в первом приближении поле давления;2) решить уравнения количества движения и энергии, используя это поледавления как заданный источник, и получить распределение скорости;3) определить поле плотности, решая уравнение неразрывности;4) по уравнению состояния определить поле давления и вернуться ко 2-мупункту данного алгоритмаК сожалению, этот алгоритм редко сходится.

Приходится прибегать кразличным ухищрениям, например, использование шахматной сетки,дополнительных уравнений для поправок скоростей и давления и т.д.Подробно эти подходы разобраны в книге Патанкара [1].Существует подход, основанный на исключении давления из системыуравнений и переходе к системе переменных функция тока – вихрь. Но этотподход имеет ограниченную область применения (см. [1]).Практика последних лет показывает, что наиболее эффективным являетсяподход, основанный не на последовательном решении каждого уравнения,входящего в систему уравнений Навье-Стокса, а на решении одногоуравнения с одной неизвестной величиной.

Только эта неизвестная величинаявляетсянескалярнойфункцией,авекторной.Соответственно,коэффициенты, входящие в полученное уравнение тоже не скаляры, аматрицы.Вэтомслучаеотпадаетнеобходимостьсначаланаходитьоднунеизвестную величину (например, скорость ), потом другую (например,давление), а потом искать способы их как-то связать.Все неизвестные рассматриваются как единое целое и определяются врезультате решения уравнения одновременно.Такой подход имеет глубокий физический смысл: в реальном газовомпотоке все газодинамические и термодинамические параметры тесно связанымежду собой; любое изменение одного параметра вызывает неизбежноеизменение других.

Таким образом, эти параметры является единым целым, ирасщепление при численном решении на отдельные части противоречитфизическому смыслу процесса.Этот метод требует существенное усложнение математического аппарата,но достигаемая эффективность с лихвой оправдывает это усложнение.1.2. Векторная форма системы уравнений Навье-СтоксаДля декартовой системы координат ( x, y, z ) все уравнения, входящие всистему уравнений Навье-Стокса, можно записать в следующей общейформе∂U ∂F ∂G ∂H+++=0∂t∂x ∂y∂z(1.1)1) В уравнении неразрывности через U , F , G, H обозначены:U1 = ρ , F1 = ρ u , G1 = ρ v, H1 = ρ w(1.2)2) В проекциях уравнения количества движения на оси x,y,z черезU , F , G, H обозначены соответственно:U 2 = ρ u , F2 = ρ u 2 + p − τ xx , G2 = ρ uv − τ xy , H 2 = ρ uw − τ xz ,(1.3)U 3 = ρ v, F3 = ρ vu − τ yx , G3 = ρ v 2 + p − τ yy , H 3 = ρ vw − τ yz ,(1.4)U 4 = ρ w, F4 = ρ wu − τ zx , G4 = ρ wv − τ zy , H 4 = ρ w2 + p − τ zz(1.5)3) В уравнении энергии через U , F , G, H обозначены:pU 5 = ρ E , F5 = ρ u  E +  − uτ xx − vτ yx − wτ zx + qx ,ρpG5 = ρ v  E +  − uτ xy − vτ yy − wτ zy + q y ,ρ(1.6)pH 5 = ρ w  E +  − uτ xz − vτ yz − wτ zz + qzρЗдесь использованы следующие обозначения:ρ - плотность; u,v,w - компоненты вектора скорости в проекциях на осиx,y,zсоответственно; p - давление; E - полная внутренняя энергия;компоненты тензора вязких напряжений и вектора плотности тепловогопотока выражаются по формулам: ∂u 2 ∂um  4 ∂u 2 ∂v 2 ∂w −−− = µ, 3 ∂x 3 ∂y 3 ∂z  ∂x 3 ∂xm τ xx = µ  2 ∂v 2 ∂um  4 ∂v 2 ∂u 2 ∂w −−− = µ, 3 ∂y 3 ∂x 3 ∂z  ∂y 3 ∂xm τ yy = µ  2 ∂w 2 ∂um  4 ∂w 2 ∂u 2 ∂v −−− = µ, 3 ∂z 3 ∂x 3 ∂y  ∂z 3 ∂xm (1.7)τ zz = µ  2 ∂u ∂v  ∂u ∂w +  , τ xz = τ zx = µ  +, ∂z ∂x  ∂y ∂x τ xy = τ yx = µ  ∂v ∂w + ∂z ∂y τ yz = τ zy = µ qx = −λгде∂um ∂u ∂v ∂w=++∂xm ∂x ∂y ∂z∂T∂T∂T, q y = −λ, q z = −λ∂x∂y∂z- дивергенция скорости,(1.8)- коэффициентµдинамической вязкости, λ - коэффициент теплопроводности.Удельная внутренняя энергия e связана с полной внутренней энергией Eсоотношением:e=E−ВыраженияU1 , U 2 , U 3 , U 4 , U 5некоторого 5-мерного вектора U1 2u + v 2 + w2 )(2можнорассматривать(1.9)каккомпонентыρ  ρu  U = ρv  , ρ w ρ E (1.10)Аналогично с использованием формул (1.2)-(1.6) строятся векторы F , G, H ρu 2 ρ u + p − τ xx ,F = ρ vu − τ yx ρ wu − τ zx ρ uH − uτ xx − vτ yx − wτ zx + qx (1.11) ρv ρ uv − τ xy,2G =  ρ v + p − τ yy ρ wv − τ zy ρ vH − uτ xy − vτ yy − wτ zy + q y (1.12)ρw ρ uw − τ xzH = ρ vw − τ yz ρ w2 + p − τzz ρ wH − uτ xz − vτ yz − wτ zz + qz (1.13)Ввод этих векторов позволяет рассматривать уравнение (1.1) каквекторную форму основных уравнений динамики вязкой жидкости.

Такаяформа позволяет существенно сократить объем записей при построениичисленного метода.Для идеального газа:e = CV T , γ =CP, R = CP − CV = (γ − 1) CV ,CV(1.14)p = ρ RT = (γ − 1) CV T ρ = (γ − 1) eρH =E+pρ,(1.15)где CP , CV - удельные теплоемкости при постоянном давлении и объемесоответственно, γ - показатель адиабаты, R - газовая постоянная..