Диссертация (1145283), страница 38
Текст из файла (страница 38)
Старкова, например [134].4.5.2. Численное решение тепловой задачи в эйлеровыхкоординатах в подвижной изменяющейся областиАналитический подход к решению нелинейной начально-краевой задачи (4.77)–(4.80) в движущейся и изменяющейся по сложному закону областипредставляется мало перспективным.
Как известно [135], [136], решение даже251линейной краевой задачи конвективной теплопроводности в области с движущимися границами аналитически в общем случае затруднительно. Болееэффективным при решениии подобных задач оказывается численный подход.Для приложений обычно представляют интерес оценка общих теплопотерь впроцессе расширения слоя и рекомендации по их возможному уменьшению.Анализ гидродинамической части задачи позволяет утверждать, что для всехдопустимых режимов характерно наличие двух этапов процесса: сравнительно короткого первого этапа [0, t∗ ], на котором слой нельзя считать тонким,и второго этапа (t∗ , tk ], в несколько раз дольше первого, на котором слой сдостаточной точностью можно считать тонким. Здесь t∗ — время окончанияпервого этапа, найденное из решения гидродинамической части задачи (соотношение (4.75)).
В связи с этим целесообразно разделить и решение тепловойзадачи на две части. На первом этапе на интервале времени [0, t∗ ] строитсячисленное решение, на втором этапе на интервале (t∗ , tk ] находится прибли-женное аналитическое решение задачи о поведении во времени средней послою (который можно считать тонким) температуры Tc (t). Приведем алгоритмы расчета обоих этапов.Решение задачи на интервале [0, t∗ ]Достаточно эффективным для численного решения задачи (4.77)–(4.80)оказалось использование двуслойной явной схемы на подвижной равномерной сетке. Расчет перехода от n-го слоя по времени к (n+1)-у осуществляетсяв два этапа.
На первом этапе рассчитываются новая конфигурация оболочки и изменение в ней температуры, соответствующее только конвективномупереносу тепла с известной скоростью u(r, t). На втором этапе найденныйпрофиль температуры интерполируется на новую равномерную сетку в изменившейся области и решается уравнение теплопроводности без конвективногослагаемого с граничными условиями (4.79), ( 4.80).Первый этап. Здесь ставится задача Коши для уравнения∂ T̃∂ T̃+ u(r, t̃)= 0,∂r∂ t̃T̃ (r, 0) = T (r, tn ) = T n (r),(4.81)(8.82)252t̃ = t − tn ,t̃ ∈ [0, τ ],τ — шаг по времени, T n (r) — распределение температуры в слое в моментвремени tn . Для малых τ допустимо считатьu(r, t̃) ≃ u(r, tn ) =A,r2A = R2 (tn )Ṙ(tn ).ОбозначимA.(4.83)r2Аналитическое решение задачи Коши (4.81)–(4.82) известно [136], оно приũ(r) =u(r) = ũ(r) (4.83) имеет видp3T̃ (r, t̃) = T ( r3 − 3At̃).nУпростим выражение аргумента в правой части, воспользовавшись малостью величины At̃/r3 на интервале τ .
Для исследуемых режимов это приводит к приближенному равенству:T̃ (r, t̃) ≃ T n (r − t̃ũ(r)),(4.84)в котором скорость ũ(r) определена формулой (4.83). Формула (4.84) даетаналитическое выражение для непрерывной функции T n (r).В дискретном варианте предлагается поступить следующим образом.Сначала находим новую равномерную сетку в изменившейся области. Пустьri — равномерная сетка в момент времени tn :ri = R(tn ) + (i − 1)∆,i = 1, . . . , N + 1,где ∆ — шаг по пространственной координате. Жидкость, находящаяся вмомент времени tn , в точке ri , в момент времени tn + τ окажется в точке r̃i∼r̃i = ri + ui τ,(4.85)а величина скорости ui , в соответствии с формулой (4.83), будет равна:ũi =A.r̃i2(4.86)Из формулы (4.86) видно, что скорость ũi с ростом номера узла i падает,поэтому равномерная в момент времени tn сетка переходит в неравномерную253сетку r̃i , сгущающуюся к внешней границе оболочки.
Кроме того, изменяетсяи общая длина расчетной области, равная толщине H(t) жидкой оболочки.Новая толщина слоя H(t) определяется в момент времени tn+1 = tn + τ поформуле:H(tn+1 ) = R̂(tn+1 ) − R(tn+1 ).(4.87)Введем в момент времени tn+1 , наряду с рассчитанной неравномернойсеткой r̃i (4.85), новую равномерную сетку r1i , определив общее количествоновых узлов Nn равенством:Nn = ближайшее целое от (H(tn+1 )/△),и положивr1i = R(tn+1 ) + (i − 1)△,i ∈ [1, Nn+1 ].(4.88)В момент времени tn значения температуры T n (ri ) в узлах старой равномерной сетки считаются известными. Температура в узлах новой равномерной сетки r1i уменьшившейся области находитс следующим образом:• припишем неравномерным узлам r̃i (i = 1, . .
. , Nn+1 ) (4.85) сгустившейсясетки прежние значения температуры Tin равномерной старой сетки ri ,а именно положимT̂i r̃i = Tin ri ,i = 1, . . . , N ;• для каждого узла j новой равномерной сетки r1j (4.88) в уменьшившейсяобласти найдем узел r̃i∗ (ближайший к узлу j) из сетки (4.85) и положимT̃˜jn = T̂i r̃∗ , j = 1, . . . , Nn+1 ;i• найденный массив температуры T̃˜jn (j = 1, 2, . . . , Nn+1 ) примем «началь-ным» на новой сетке (4.88), и по нему на втором этапе рассчитаем температуру Tjn+1 в момент времени tn+1 .Второй этап. На втором этапе ищем численное решение теплового урав-нения (4.77) без конвективного слагаемого в левой части в изменившейся области толщиной H(tn+1 ) (4.87)ν ∂∂T= 2∂tr ∂r∂Tr2.∂r(4.89)254Рисунок 4.14 – Новая конфигурация оболочки и изменение в нейтемпературы на первом этапе.На равномерной сетке (4.88) уравнение теплопроводности (4.89) аппроксимируется с помощью явной разностной схемы. В качестве начального массива T n используются значения температуры T̃˜n , рассчитанные в узлах r1jji(4.88).
Граничные условия (4.79), (4.80) аппроксимируются с помощью введения дополнительных фиктивных узлов и последующего их исключения всоответствии с разностным аналогом уравнения теплопроводности (4.89). Алгоритм расчета нового массива температуры Tin+1 на втором этапе по явнойсхеме записываются следующим образом:Tin+1 = Tin +τν nn(Ti−1 − 2Tin + Ti+1),2△i = 1, ..., M − 1,M = Nn+1 ,2τ ν n2aτ νn(T−T)+(T∗ − T0n ),102△△2τ ν n2τ νb n2τ νγn 4nn+1n(TM) − (T ∗ )4 .= TM+ 2 (TM(TM − T ∗ ) −TM−1 − TM ) −△△△В каждый момент времени значение средней температуры Tc определяT0n+1 = T0n +ется равенствамиZR̂(t)3Tc (t) =T (r, t)r2 dr,K(4.90)R(t)K = R̂(t)3 − R(t)3 ,и рассчитывается по найденному массиву Tin+1 численно, например по формуле Симпсона.
Как показали расчеты необходимо пользоваться максимальноточным алгоритмом численного интегрирования.255Таким образом, переход от n-го слоя по времени к (n + 1)-у осуществляется в два этапа. На рисунке 4.15 представлены результаты расчета попредложенному алгоритму профиля температуры в момент времени t = 5.Уже через 5 единиц безразмерного времени толщина слоя становится малойвеличиной, более чем на порядок меньшей исходной толщины и возникаетвозможность перехода ко второму этапу решения тепловой части задачи.Рисунок 4.15 – Профиль температуры в жидком слое в моментбезразмерного времени t = 5.4.5.3.
Модель поведения средней по слою температуры жидкостиКак отмечалось ранее, значение времени перехода t∗ ко второй части за-дачи находится независимо от тепловых процессов из решения гидродинамической части (соотношение (4.75)). На протяжении всего интервала времени[t∗ , tk ] слой можно считать тонким, поэтому достаточно следить только заизменением средней по слою температуры Tc .Утверждение 6. ◮ Динамика средней по слою температуры Tc (t) наинтервале (t∗ , tk ] в первом приближении по δ(t) = H(t)/R(t) описываетсяобыкновенным дифференциальным уравнениемm8 (t) m4 (t) + Tc (t)m1 (t) + Tc (t)3 m2 (t) + Tc (t)4 m3 (t)+Ṫc (t) = −m5 (t) + Tc (t)3 m6 (t) + Tc (t)2 m7 (t)+Tc (t)m9 (t) + m10 (t),в котором функции mi (t) (i = 1, ..., 10) представляют собой аналитическиезависимости от R(t), R̂(t) и параметров задачи K, a, b, ν, γ, T∗ , T ∗ .
◭256Доказательство. ⊲ Найдем приближенное уравнение для динамики средней по слою температуры, пользуясь малостью функционального параметраδ(t) на данном интервале времениδ(t) = H(t)/R(t),(4.91)H(t) — толщина слоя жидкости. При t ≥ t∗ , как отмечалось, величина δ(t) ≪1.Выведем уравнение динамики средней температуры Tc (t) (4.90). Для этого запишем тепловое уравнение (4.77) с учетом выражения для поля скорости(4.83) и проинтегрируем его по r в пределах от R(t) до R̂(t).ZR̂(t)R(t)∂T (r, t) 2r dr + A∂tZZR̂(t)R(t)∂T (r, t)dr = ν∂rZR̂(t)R(t)∂∂r∂T(r,t)r2dr.∂rR̂(t)∂T (r, t) 2r dr = −R(t)2 Ṙ(t) TR̂ (t) − TR (t) +∂tR(t)∂T(r,t)∂T(r,t)2 .−νR(t)+ν R̂(t)2∂r R̂(t)∂r R(t)(4.92)Воспользуемся этим равенством для вычисления скорости изменения Ṫc (t)средней температуры Tc (t)3 ddṪc (t) = Tc (t) =dtK dtZR̂(t)T (r, t)r2 dr.R(t)В соответствии с формулой взятия производной от интеграла, зависящего от параметра, имеемddtZR̂(t)T (r, t)r2 dr =R(t)ZR̂(t)R(t)∂T (r, t) 21r dr +∂t3!dR̂(t)3dR(t)3TR̂ (t) −TR (t) .dtdt(4.93)Первое слагаемое в правой части задается равенством (4.92).
Кроме того,учтем равенства:R̂(t)3 = R(t)3 + K,K = const =⇒dR̂(t)3dR(t)3== 3R(t)2 Ṙ(t).dtdt257В результате приходим к следующему уравнению, моделирующему поведение средней температуры Ṫc (t) расширяющегося слоя жидкости:3ν bR̂(t)2 TR̂ (t) − T ∗ +Ṫc (t) = −K24∗ 42+γ R̂(t) TR̂ (t) − (T ) + aR(t) (TR (t) − T∗ ) .(4.94)Найдем выражения для температуры на внутренней и внешней поверх-ностях слоя TR (t), TR̂ (t) в первом приближении по малой величине δ(t) =H(t)/R(t).
Аппроксимируем реальный профиль T (r, t) температуры параболической зависимостью следующего вида:r2aT (r, t) = ϕ(t) + ψ(t)(r − R(t)) +(TR − T∗ ).2R(t)2(4.95)Вид зависимости выбирался, исходя из численного решения тепловойзадачи на начальном этапе. Также нетрудно убедиться в том, что при этомграничное условие на внутренней поверхности (4.79) удовлетворяется автоматически, а именно из (4.95) следует:2ra∂T (r, t)= 2ψ(t)(r − R(t)) +(TR (t) − T∗ ).∂r2R(t)∂T (r, t) 2R(t)a(TR (t) − T∗ ) = a(TR (t) − T∗ ).=2ψ(t)(R(t)−R(t))+∂r r=R(t)2R(t)Выразим ϕ(t) через Tc (t).
Для этого подставим зависимость T (r, t) (4.95) ввыражение для Tc (t) (4.90), проинтегрируем и воспользуемся малостью H(t),представив R̂(t) в виде:H(t)= R(t)(1 + δ(t))R̂(t) = R(t) + H(t) = R(t) 1 +R(t))(ZZ R̂(t)R̂(t)3ψ(t)(r − R(t))2 r2 dr +ϕ(t)r2 dr +Tc (t) =KR(t)R(t)Z3 R̂(t) r4 a+(TR (t) − T∗ )dr.(4.96)K R(t) 2R(t)R R̂(t)Покажем, что интеграл R(t) (r − R(t))2 r2 dr в первом приближении по δ(t)равен нулю.ZR̂(t)R(t)R̂(t)5 − R(t)5(r − 2rR(t) + R(t) )r dr =−522225833R̂(t)4 − R(t)42 R̂(t) − R(t)+ R(t)=−2R(t)43111((1 + δ(t))5 − 1) − ((1 + δ(t))4 − 1) + ((1 + δ(t))3 − 1) .= R(t)5523(Здесь использовано равенствоR̂(t)R(t)= 1 + δ(t).) Учтем, что в первом прибли-жении по δ(t) справедливо приближенное равенство:(1 + δ(t))k ≃ 1 + kδ(t),позволяющее записатьZ R̂(t)11122 2(r − 2rR(t) + R(t) )r dr ∼(5δ(t)) − (4δ(t)) + (3δ(t)) = 0.=523R(t)Таким образом, средняя температура Tc (t) (4.96) в первом приближениипо δ(t) имеет вид:Tc (t) ≃ ϕ(t) +3aR(t)4 δ(t)(TR (t) − T∗ ) =⇒2K3aR(t)4 δ(t)(TR (t) − T∗ ).2KНайдем выражение для TR (t) из зависимости для T (r, t) (4.95):ϕ(t) ≃ Tc (t) −TR (t) =R(t)aT∗2.R(t)a1−2ϕ(t) −(4.97)(4.98)Подставим это соотношение в равенство (4.97) и найдем вид функции ϕ(t):3aR(t)a+δ(t)R(t)4 T∗Tc (t) 1 −2K2ϕ(t) =.R(t)a3a41−+δ(t)R(t)2K2Можно упростить функцию ϕ(t), используя приближенное представлениеδ(t)R(t)4 ∼ H(t),=K3δ(t)и записать окончательное приближенное выражение для ϕ(t) в следующемвиде:ϕ(t) ≃ Tc (t)(1 − aR(t)/2) + aR(t)T∗ /2.(4.99)259Найдем функцию ψ(t), входящую в аппроксимацию T (r, t) (4.95).
Дляэтого предварительно докажем, что из приближенного равенства (4.99) следует в первом приближении по δ(t), приближенное равенство:TR(t) ∼= Tc (t).(4.100)Из (4.95) при r = R имеем (4.98).Подставим сюда соотношение для ϕ(t) (4.99), в результате получим:TR(t) =Tc (t)(1 − aR(t)/2) + aR(t)T∗ /2 − aR(t)T∗ /2 ∼= Tc (t).1 − aR(t)/2(4.101)Выражение для ψ(t) найдем из граничного условия на внешней поверхности слоя R̂(t) (4.80). Для этого запишем выражение для TR̂ (t) из аппроксимации (4.95) и подставим найденное приближенное соотношение для ϕ(t)(4.99). После ряда преобразований получим следующее выражение для TR̂ (t)в первом приближении по δ(t):TR̂(t) ∼= Tc + aR(t)δ(Tc (t) − T∗ ) + ψ(t)H(t)2Запишем это соотношение, учитывая равенство H(t) = δ(t)R(t), в виде:TR̂(t) ≃ Tc (t)(1 + ε(t)) + ψ(t)H(t)2 ,ε(t) = aH(t) (1 − T∗ /Tc (t)) .(4.102)В граничные условия на внешней поверхности слоя R̂(t) (4.80) входитслагаемое вида TR̂4 (t).
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
