Математическое моделирование процессов гидродинамики

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

Рассмотрим состав этих условий.

• 1. Геометрические условия характеризуют форму и размеры тела, в котором протекает процесс гидродинамики. Например, если рассматривается круглая труба, то для определения распределения скорости и давления необходимо знать ее длину и диаметр. При других формах труб и каналов в геометрических условиях необходимо учитывать все их характерные особенности, связанные с геометрической конфигурацией объекта.
• 2. Физические условия характеризуют физические свойства среды (вязкость, плотность и др.).
• 3. Граничные условия характеризуют особенности динамического взаимодействия исследуемой среды на границах рассматриваемого объекта.
• 4. Временные, или начальные, условия характеризуют состояние среды в исходный (начальный) момент времени.

Перечисленные условия в совокупности определяют одно (конкретное) явление гидродинамики, поэтому они и названы условиями однозначности, или условиями единственности.

Для тела определенной геометрической формы с определенными (известными) физическими свойствами условия однозначности сводятся к заданию начального и граничного условий. Эти условия в совокупности называются краевыми условиями – начальное условие является временным краевым условием, а граничное условие – пространственным. Дифференциальные уравнения гидродинамики совместно с краевыми условиями составляют краевую задачу гидродинамики. Для установившегося (стационарного) процесса в задании начального условия нет необходимости, и в этом случае краевая задача будет состоять из уравнений гидродинамики и граничных условий.

Аналитические решения гиперболических уравнений в условиях гидравлического удара

Ниже будут рассмотрены конкретные примеры применения математического моделирования к решению сложных гидродинамических задач. И в частности, это задача о распределении скорости в динамическом пограничном слое, а также задача о распределении скорости при разгонном течении Куэтта с учетом релаксационных свойств жидкости (см. параграфы 5.7–5.9). В данном параграфе рассматриваются решения гиперболических уравнений в условиях гидравлического удара, который представляет собой скачок давления в какой-либо системе, заполненной жидкостью, вызванный крайне быстрым изменением скорости потока этой жидкости за очень малый промежуток времени.

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

(5.11)

(5.12)

где р – давление; х – продольная координата; t – время; – плотность; – скорость; – коэффициент гидравлического сопротивления (коэффициент трения), зависящий от скорости течения среды; – скорость звука в капельной упругой жидкости, текущей в трубе с упругими стенками, т.е. скорость распространения ударной волны – волны гидравлического возмущения (соотношение справедливо прии; k – модуль упругости жидкости;– толщина стенки трубы; d – диаметр трубы; Е – модуль упругости материала стенки трубы.

Решение нелинейных уравнений (5.11), (5.12) возможно лишь путем численного интегрирования. Для упрощения нелинейного уравнения (5.11) И. А. Парный [20] предложил способ линеаризации второго слагаемого в его правой части, приняв множительпостоянным и равным его среднему значению подлине трубы и времени:

Величина 2а для ламинарного режима, учитывая формулу Пуазейля

где– число Рейнольдса; v – кинематическая вязкость, приводится к виду

При турбулентном режиме течения второе слагаемое в правой части уравнения (5.11) осредняется в некотором интервале скоростис заменой кривой изменения функциинекоторым эквивалентным отрезком прямой. Тогда формула для величины 2а приводится к виду

,

где принимается .

С учетом рассмотренной линеаризации уравнение (5.11) принимает вид

(5.13)

Уравнения (5.12), (5.13) обычно решаются операционными методами либо методом контурного интегрирования, являющегося разновидностью операционных методов.

С целью упрощения получения аналитического решения уравнения (5.12), (5.13) в некоторых случаях сводятся к одному гиперболическому уравнению относительно давления или скорости. Дифференцируя уравнение (5.12) по переменной t, а уравнение (5.13) – по переменной х и сопоставляя полученные соотношения, находим

(5.14)

Уравнение (5.14) с учетом уравнения (5.12) принимает вид

(5.15)

Точно так же можно получить гиперболическое уравнение и для скорости

(5.16)

В качестве конкретного примера получения аналитического решения уравнений (5.15), (5.16) рассмотрим краевую задачу о распространении скачка давления в трубопроводе с неподвижной в исходном состоянии жидкостью. Допустим, что в сечении х = 0 произошло скачкообразное изменение давления, а сечение х= l перекрыто (скорость равна нулю). Требуется найти распространение давления по длине трубопровода во времени [20].

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

Математическая постановка задачи по определению давления в данном случае будет

(5.17)

где р0 – начальное давление в трубе; р1 – давление, приложенное в точке х = 0 и действующее в течение всего времени процесса вплоть до установления стационарного состояния (); l – длина трубопровода.

Из решения задачи (5.17) можно получить полную информацию о распределении скорости и давления по длине трубопровода во времени. Распределение скорости в сечении х= 0 при известном из решения задачи (5.17) давлении может быть найдено путем интегрирования уравнения (5.12).

Последовательность получения точного аналитического решения задачи (5.17) путем использования метода разделения переменных Бернулли – Фурье дана в [20]. Ниже будет рассмотрен метод получения точного аналитического решения путем совместного использования классического аналитического метода разделения переменных и ортогональных методов взвешенных невязок.

С целью упрощения математической постановки задачи и процесса получения аналитического решения, а также для нахождения решения наиболее общего вида введем следующие безразмерные переменные и параметры:

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

С учетом принятых обозначений задача (5.17) приводится к виду

(5.18)

Для удобства получения аналитического решения сделаем замену независимой переменной

Относительно переменной с, задача (5.18) примет вид

(5.19)

(5.20)

(5.21)

Решение задачи (5.19) – (5.21) принимается в виде

(5.22)

Подставляя соотношение (5.22) в уравнение (5.19), получаем

(5.23)

(5.24)

где – некоторая постоянная.

Подставляя соотношение (5.22) в граничные условия (5.21), получаем

(5.25)

Решение краевой задачи Штурма – Лиувилля (5.24), (5.25) принимается в виде

(5.26)

Собственные функцииввиду однородности уравнения (5.24) с точностью до постоянного множителя (который в данном случае можно принять равным единице) находятся из (5.26).

Очевидно, что соотношение (5.26) удовлетворяет граничным условиям (5.25). Подставляя соотношения (5.26) в уравнение (5.24), получаем формулу для определения собственных чисел:

Характеристическое уравнение для уравнения (5.23) будет

(5.27)

Уравнение (5.27) для каждого собственного числа имеет два корня

(5.28)

Если дискриминант , то из формулы (5.28) для каждого собственного числа будем иметь два действительных отрицательных корня и . С учетом найденных значенийирешение уравнения (5.23) для каждого собственного числа будет

, (5.29)

где , – неизвестные коэффициенты.

Подставляя выражения (5.26), (5.29) в формулу (5.22), находим

(5.30)

Каждое частное решение (5.30) точно удовлетворяет уравнению (5.19) и граничным условиям (5.21), но ни одно из них не удовлетворяет начальным условиям (5.20). Для их выполнения составим сумму частных решений

(5.31)

Неизвестные коэффициенты c1k и с2k определяются из начальных условий (5.20). Подставляя выражение (5.31) во второе равенство (5.20), получаем

(5.32)

Подставляя выражение (5.31) в первое равенство (5.20), с учетом равенства (5.32) получаем

(5.33)

Соотношение (5.33) представляет разложение единицы в ряд Фурье по собственным функциям краевой задачи Штурма – Лиувилля на отрезке |0; 1|. Умножим обе части соотношения (5.33) наи найдем интеграл в пределах отдо:

(5.34)

Соотношение (5.34) ввиду ортогональности косинусов принимает вид

(5.35)

Вычисляя интегралы в выражении (5.35), находим

После определения постоянных интегрирования точное аналитическое решение задачи (5.19) – (5.21) находится из равенства (5.31).

Если в соотношении (5.28) , то будем иметь следующие /1ва комплексных корня:, ,где;;

Частные решения уравнения (5.23) будут

Па основе частных решений запишем общий интеграл уравнения (5.23):

(5.36)

где , – неизвестные постоянные.

Соотношение (5.36) можно переписать следующим образом:

(5.37)

Используя формулы Эйлера, , соотношение (5.37) приведем к виду

(5.38)

Соотношение (5.38) с учетом обозначений; будет

(5.39)

Подставляя соотношения (5.26), (5.39) в решение (5.22) и составляя сумму частных решений, находим

(5.40)

Для определения постоянныхииспользуются начальные условия (5.20). Подставляя выражение (5.40) во второе равенство (5.20), получаем

Отсюда находим

Подставляя выражение (5.40) в первое равенство (5.20), находим

(5.41)

Соотношение (5.41) представляет разложение единицы в ряд Фурье по собственным функциям краевой задачи Штурма – Лиувилля на отрезке [0; 1]. Поступая так же, как при определении коэффициента(см. выше соотношение (5.33)), будем иметь

После определения постоянных интегрирования и точное аналитическое решение задачи (5.19) – (5.21) в замкнутом виде находится из выражения (5.40).

На рис. 5.1–5.7 даны результаты расчетов по формулам (5.31), (5.40) конкретной задачи о распределении давления в находящейся в стальном трубопроводе нефти при следующих условиях: v = 7,65•10-6 м2/с; d = 0,1м; λ = 0,01; ʋ0 = 0; ʋ1= 2 м/с; с = 998м/с; ρ = 840кг/м3; k = 1500МПа; E = 2-105 МПа; р0 = 10атм = 10•105 Па; р1 = 100 атм = 100 • 105 Па.

Расчеты выполнялись для следующих длин трубопровода: l = 500, 100 и 1 км.

Соответственно указанным длинам трубопровода число Ро, для ламинарного режима течения было 6,648•10-3; 0,1662; 1662,0277, для турбулентного режима – 0,8964 • 10_3; 0,02241; 224,1 соответственно.

Рис. 5.1. Перемещение фронта гидравлической волны по координате ξ, во времени Fo (For = 6,648 • 10-3)

Рис. 5.2. Распределение давления в трубопроводе при For = 6,648•10-3

Рис. 5.3. Распределение давления в трубопроводе при For = 10-7

Рис. 5.4. Распределение давления в трубопроводе при For = 0,1662

Рис. 5.5. Распределение давления в точке ξ = 0 во времени Fo при For= 0,1662 (n = 100, где п – число членов решений (5.31), (5.40))

Рис. 5.6. Распределение давления в трубопроводе при For= 0,1662 (п = 100)

Рис. 5.7. Распределение давления в трубопроводе при For = 0,1662 (п = 100)

Анализ результатов расчетов позволяет заключить, что при любых значениях числа For изменение давления характеризуется движением гидравлической волны, на фронте которой наблюдается скачок давления от его значения на фронте до величины давления невозмущенного потока. Область, находящаяся за пределами фронта гидравлической волны, оказывается невозмущенной, и давление здесь равно начальному давлению р0. Отмечается линейная закономерность движения фронта гидравлического возмущения по пространственной переменной во времени (см. рис. 5.1), что подтверждается исследованиями решения уравнения вида (5.17), выполненными другими авторами [1] применительно к определению температуры в теле с учетом конечной скорости распространения теплоты.

В зависимости от величины числа For наблюдается существенное отличие получаемых результатов. Так, при очень малых его значениях скачок давления имеет место лишь на некоторых начальных участках трубопровода. Например, при For =6,648-10-3 скачок давления наблюдается лишь в диапазоне 0 ≤ Fo ≤ 0,07 (см. рис. 5.2). Для всех Fo > 0,07 для этого значения For решение для безразмерного давления полностью совпадает с решением классического параболического уравнения теплопроводности, т.е. возрастание давления не сопровождается возникновением его скачков. Для очень малых значений For скачок давления наблюдается лишь на незначительном по длине участке трубопровода вблизи сечения х=l и при малых значениях времени. Например, при For = 10-7 скачки давления практически заканчиваются при Fo≈10-6. Фронт гидравлического возмущения для данного момента времени перемещается лишь на величину ξ = 0,003, что составляет 0,3% от всей длины трубопровода (см. рис. 5.3).

На всей остальной части трубопровода при Fo > 10-6 повышение давления происходит без скачков вплоть до достижения стационарного состояния, когда давление в сечении х=l становится равным давлению, заданному граничным условием первого рода в сечении х= 0. Решение задачи (5.19) – (5.21) на временном участке 10-6≤Fo≤∞ (при For = 10-7) полностью совпадает с классическим точным аналитическим решением параболического уравнения теплопроводности при симметричных граничных условиях первого рода [12] (см. рис. 5.3). Значения For, близкие к нулевым, могут быть в следующих случаях: большая длина трубопровода; высоковязкая жидкость; большое значение коэффициента гидравлического сопротивления.

При значениях числа For > 0,07 скачок давления в трубопроводе наблюдается вплоть до момента времени, когда фронт возмущения достигает координаты ξ = 0 (рис. 5.8– 5.10). Например, для For = 0,1662 фронт гидравлического возмущения достигает координаты ξ = 0 при Fo ≈ 0,5477 (см. рис. 5.4, 5.5). Давление на фронте в этот момент времени составляет Θ = 0,6. При дальнейшем увеличении времени наблюдается обратная волна гидравлического возмущения со скачком давления в сторону, противоположную скачку давления в прямой волне. По сути, это есть скачок давления, вызывающий гидравлический удар.

Рис. 5.8. Распределение давления в точке ξ = 0 во времени Fo при For= 0,1662 (n = 100)

Рис. 5.9. Распределение давления в трубопроводе при For= 224,1 (п = 100)

С увеличением времени наблюдается периодическое изменение давления в каждой точке трубопровода во времени вплоть до наступления стационарного состояния, при котором давление по всей длине трубопровода становится равным давлению p1=100 атм = 100 • 105 Пa, заданному на входе, т.е. в сечении х = 0 (см. рис. 5.8, 5.10).

Рис. 5.10. Распределение давления в точке ξ = 0 во времени Fo при For = 224,1 (п = 100)