CC BY
64
УДК 519.676
И.А. Соловьев, П.В. Зубков МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ СТОХАСТИЧЕСКИХ ЯВЛЕНИЙ ПЕРЕНОСА ТЕПЛОТЫ И МАССЫ
Стохастические дифференциальные уравнения предложены на основе концепции локальных средних значений. Даны постановки начально-краевых задач для плотности вероятности, которые описывают стохастические поля тепло-массопереноса.
Стохастические дифференциальные уравнения, дисперсия,
устойчивость решений.
I.A. Soloviev, P.V.Zubkov
MATHEMATICAL MODELLIG STOCHASTIC PHENOMENA OF HEAT AND MASS TRASFER
The stochastic difference equations by the conception of local mean values are suggested. Statements of initial-boundary problems for probability densities, which described stochastic heat and mass transfer, are given.
Stochastic differential equations, dispersion, stability of solutions.
Введение
Для тепловых и диффузионных полей фундаментальные детерминированные уравнения впервые получены в работах Ж.Б.Ж. Фурье, А. Фика, Ш. Соре, Л. Дюфура. Библиография по развитию детерминированной теории теплопроводности, диффузии и тепломассопе-реноса изложена в [1].
С точки зрения управления температурными полями в производственных условиях более перспективными представляются модели, которые основаны на решении задачи для функции плотности распределения вероятностей (в дальнейшем ФПРВ). Пионерские работы по стохастическому описанию температурных полей принадлежат Л. Больцману. Наиболее характерные трудности в традиционном способе описания случайных тепловых полей проявляются в модели, предложенной в работе Р.Ф. Куртайна [2] В ней дан аналог уравнения Эйнштейна-Фоккера-Планка-Колмогорова в функциональных производных. Предложенная Р.Ф. Куртайном модель широко не использовались из-за очевидных трудностей ее аналитической и численной реализации. Другой подход осуществлен в работах А.Г. Мадеры [3]. Отметим, что, несмотря на перспективность предложенной А.Г. Мадерой модели стохастической теплопроводности, эта модель имеет тот недостаток, что уравнение для подсчета корреляций имеет четвертый порядок. Кроме того, в работе [3] не установлена связь между предложенными уравнениями и соответствующим уравнение Колмогорова-Эйнштейна-Фоккера-Планка для случайных процессов. Библиография по развитию теории стохастических полей теплопроводности, диффузии и тепломассопереноса изложена в [4].
В настоящей работе предложен вариант изучения явлений переноса массы и теплоты на основе предложенных в работе [5] стохастических моделей, применение которых, прежде
всего, направлено на использование их в инженерных расчетах. В постановочной части хотелось бы обозначить актуальность проблемы устойчивости решений, решаемые при помощи исследования поведения дисперсии.
Метод локальных средних значений
В основе стохастического описания упомянутых явлений лежит феноменологический метод локальных средних. Поясним суть этого метода. Он исходит из математического описания физического эксперимента. Как правило, физики используют вначале немногочисленные данные, из которых затем нужно установить закон. В стохастике есть понятие среднего значения (математического ожидания), оно при ФПРВ П(г, х, С), где, например, г есть время, х - пространственная координата, а С - вещественная характеристика случайного явления, записывается так: | СП(г, х, С)йС . Локальное среднее значение определяется как
подынтегральное значение, то есть СП(г, х, С~)йС . На основе этого метода выведены стохастические уравнения теплопроводности, диффузии, тепло-массопереноса, а также волновые и поставлены начально-краевые задачи для этих уравнений [5]. Общая схема вывода одинакова для всех перечисленных задач, поэтому опишем эту схему на примере получения стохастического уравнения диффузии.
Схема вывода стохастических уравнений для случайных полей
Предположим, что между каждыми тремя соседними случайными диффузионными процессами, которые расположены вблизи друг от друга на одной прямой на расстоянии Н : х. - Н = х.-1, х., х. + Н = х.+1 за время т после г{ осуществляется обмен случайными значениями концентрациями С с известными вероятностями перехода р(гг + т; х. ± Н). Средние локальные средние значения в каждом узле этой сетки в момент времени + т соответственно равны: СП(£г +т, х3-1, С)ёС, СП(гг +т, х., С)ёС и СП(£г +т,х1+1,С)ёС . Поскольку локальное среднее значение в центральном узле х. в момент времени +т формируется в результате случайной связи только с ближайшими двумя соседними узлами х . ± Н и учитывает их состояние только в предыдущий момент времени , то локальное среднее значение во внутренних узлах сетки х. в момент времени + т формируется следующим образом. Во-первых, из двух соседних узлов х.±1 в центральный узел х. за время гг+1 - г{ =т может произойти перенос локального среднего значения. В результате такого переноса локальное среднее значение характеристики в х . может измениться на величину, равную
(Вт/ Н2)СП(г. +т,х ,С^С + (Вт/ Н2)СП(г. + т,х ,С^С, где В - коэффициент диффузии. А в узлы
I 1-1 I 1+1
х.±1 из центрального придут локальные средние значения, равные 2((Вг)/Н2)СП(г +т,х ,С)с1С. В
‘ 1
результате этого локальное среднее значение в центральном узле изменится на величину (1 -2• ((Вт)/к2))СП(г +т,х ,С)йС). В итоге получим разностное уравнение:
г 1
СП (г. + т,х ,С)йС « (Вт/ Н2)СП(т. +т,х ,С)йС +
1 1 1 3-1 (1)
+(1 - 2(( Вт) / Н2))СП(г +т, х., С )йС + (Вт / Н2)СП(и +т, х , С )йС + погрешность.
11 I з+1
Здесь не станем подробно описывать элементы процесса вывода разностного уравнения, связанного с матрицей вероятностей перехода из узла в узел и аппроксимации этих ве-
роятностей экспоненциальным законом. Подробно об этом написано в работе [5].
Возникающая при этом погрешность марковского поля компенсируется конечной (или бесконечной) суммой, которая определяется моментными функциями. После перехода к пределу при т ^ 0 и к ^ 0 получим следующее дифференциальное уравнение в частных производных для ФПРВ:
Э(П (г, х, С)) / Эг = пЭ2 (П(г, х, С)) / Эх2 + (погрешность) . (2)
С учетом того, что в каждой точке предполагалась непрерывность случайного поля, а также его марковость, естественно выбрать среди всего возможного множества записи погрешности ту, которая соответствует правой части уравнения Колмогорова-Эйнштейна-Фоккера-Планка [4] при пЭ2(П(г,х,С))/Эх2 = 0 и равенстве нулю члена, отвечающего за скачкообразность. В итоге получим следующее стохастическое уравнение диффузии:
ЭП(г,х,С) /Эг = ПЭ2П(г,х,С) /Эх2 -Э(/П) / ЭС + 0,5Вд2П / ЭС2, (3)
где /(г, х, С) - коэффициент сноса (функция источников (стоков)); В - коэффициент диффузии марковского поля.
Постановки начально-краевых задач для случайных полей
Полная постановка начально краевой задачи для описания диффузии (и переноса теплоты) такова:
ЭП(г,х,С) /Эг = ПЭ2П(г,х,С) /Эх2 -Э(/П)/ЭС+0,5ВЭ2П/ЭС2,г > 0,хе (0; I),Се (0;■+»); (4) П(0,х,С) = Пнач(х,С),хе (0; +~),Се (0; +~); П(г,0,С) = ПГ(г,С)/Эх2,г > 0,Се (0; +~); (5)
П(г,I, С) = Пг (г, С), г > 0, С е (0; +~); (6)
П(г, 0, -~) = П (г, 0,0) = П (0, х, 0) = П(0, х, +~) = П(г, 1,0) = П (г, I, +~) = 0. (7)
Здесь П - случайный коэффициент диффузии, для которого известна ФПРФ, в данной работе в отличие от работы [5], где П - известная константа, П имеет нормальное распределение: ф(%) = (1 / (ол/2п) ехр(-(% - ^0)2 / (2о2)) , где % - случайный фактор, не зависящий от С.
Стохастический аналог задачи классической Стефана [1] имеет вид:
С р (ЭП1 (г,х,Т) /Эг) = (Э2П1 (г,х,Т)/Э2х2) + 0,5В(Э2П1 (г,х,Т) /Э2Т2); (8)
г > 0, хе (0, s(г)),T е (0,+^);
С2Р2(ЭП 2(г, х,Т)/Эг) = Л2(Э2П2(г, х,Т)/Э2 х2)+0,5В(Э2П2(г, х,Т)/Э2Т2); (9)
г > 0, хе (5(г); X );Т е (0, +<^);
П1(г = 0,х,Т) = 01(х,Т), хе [0,я(г)], Те [0, +~); (10)
П2(г = 0, х,Т) = 02(х,Т), хе [5(0), X], Те [0, +~); (11)
П1(г,х = 0,Т) = Д,(г,Т), П2(г,х = X,Т) = Ах(г,Т), ^0, Те [0, +~); (12)
П1(г,х = 5(г),Т) = П2(г,х = 5(г),Т) = ф(г,Т), ^0, Те [0, +~) (13)
Д^(ЭП1(г, х=5 (г)+0,Т)/Эх)-Д^Э^^, х=5(г) - 0),Т)/Эх)=-Ьр(^у / 4г)Ч(г;Т); (14)
ТО
г>0,Те (0,+<х>); | Т^(г;Т)сГ=1.
0
+х> +то +то
IТП2(г, х = 0, Т)^Т = | ТАх (г,Т)^Т = Т2х (г), г > 0; | Тф(г,Т)ЛТ = Т^, г > 0. (15)
0 0 0
Функцию ^(г, Т) удобно представлять в виде ^(г, Т) -—Т (ду / дТ), где IУ(г;Т )ЯТ-1.
0
Здесь также все коэффициенты, перечисленные раннее, имеют стохастический характер, описываемый нормальным законом. Отметим, что все характеристики для температурного поля в твердой фазе обозначены нижним индексом 1, в жидкой - нижним индексом 2; Х1, I -1,2 , р, I -1,2 , Сг, I -1,2 - соответственно коэффициент теплопроводности, плотность, коэффициент теплоемкости.
Функция-преобразователь решения задачи о среднем значении в решение задачи о нахождении дисперсии
Все доказательства о нормировке ФПРВ на единицу, о среднем значении и дисперсии приведены в работе [5]. Здесь мы опишем, как можно рассчитать дисперсию случайного поля на основе решения соответствующей детерминированной задачи с использованием введенной нами функции-преобразователя.
Рассмотрим случай, когда влиянием коэффициента сноса / (г, х, С) и коэффициента диффузии В марковского поля можно пренебречь. Тогда для среднего значения задачи (4)-
(7) М(1) (г, х) - | СП (г,х, С)йС получается следующая постановка:
дМ(1)(г,х)/дг -ад2М(1)(г,х)/дх2, г > 0, хе (0,1); (16)
М(1)(0,х) -г(х), хе (0,1); М(1)(г,0) -в(г), г > 0;М(1)(г,I)-ц(г), г > 0. (17)
+ж>
Для моментной функции второго порядка М(2)(г,х) - | С2П(г,х,С)ёС постановка задачи имеет вид:
дМ(2)(г,х)/дг - ад2М(2)(г,х)/дх2, г > 0, хе (0,1); (18)
М(2)(0, х) - б(х) -^ + гЧх), х е (0,1); (19)
М(2) (г, 0) - Р(г) - о-2 + в2 (г), г > 0; (20)
М(2)(г,I) -Щ(г) -о2 +п2(г), г > 0. (21)
Здесь а - коэффициент температуропроводности.
Введем две функции р1 (С) и р2 (С), которые назовем функциями-преобразователями решения задачи для средних значений в решения задач для ФПРВ и моментной функции второго порядка. Эти функции должны удовлетворять следующим условиям:
I р (С)ёС -1,1 -1,2, | Ср1(С)ёС - 0, | Ср2(С)ёС -1, | С2р(С)ёС -а1,
—X —X —X —X
+х
| С2р2(С)^С -а2, при этом получается решение задачи для ФПРВ следующего вида:
П(г, х, С) = р (С) - (р (С) - р2 (С))М(1) (х, г).
Решение соответствующей задачи для моментной функции второго порядка:
М(2) (г, х) = а - (а1 - а2 )М(1) (х, г), соответственно дисперсия а2 (г, х) случайного температурного поля имеет вид:
а2(г, х) = а - (а\ -а2)М (1)(х, г) - (М (1)(г, х))2 при г > 0, х є (0; I); (22)
а2(0, х) = аі - (аі -а2Жх) -^2(х), х є (0,1); (23)
a2(t,0) = а1 — (а1 —a2)P(t) — в2(t),t > 0; (24)
a2 (t, l) = а — а — а2 )n(t) — r/2 (t), t > 0. (25)
Введенная функция преобразователь позволяет в указанных ограничениях использовать все найденные к настоящему моменту численные и аналитические решения задач о распространении не только концентрации, но и теплоты и волн.
Исследование устойчивости формы фронта с помощью анализа поведения дисперсии
Введем понятие устойчивости решения задачи для средних значений.
Определение. Решение задачи для описания средних значений случайного поля назовем асимптотически устойчивым в стохастическом смысле, если выполняется равенство:
lima2(t, x, у, z) = L = const <ж. (26)
При начальном значении дисперсии и ее граничных значениях равных нулю и L = О решение будем называть абсолютно устойчивым.
В последнем случае получается аналог асимптотической устойчивости по Ляпунову.
С помощью функции-преобразователя получены некоторые оценки дисперсии случайного теплового поля и дисперсии формы фронта фазового перехода на основе точного решения классической задачи Стефана [1] и стохастической ее постановки (8)-(15). Дисперсия температурного поля в твердой a12( x, t) (нижний индекс — 1) и жидкой а^( x, t) (нижний индекс — 2) фазах имеет вид:
of ‘ x, ґ ) = M10 + (m Іавл— M10) • Ф( xI (2^°і )I Ф(а I (2^ )) — (27)
—{T10 + (Tnmm — M10) • Ф( x I (2^) I Ф(а I (2^)) I2
*2 (x, t) = M 2нач — (M 2на ч — M Щ) • (1 — Ф( x I (2^07) I Ф(а I (2^/07)) —
—{^ч — (T2нaч — )I(1 — Ф(а I (^л/О^Ґ)) • (1 — Ф(x I (^7)) |2 .
Зона неустойчивости формы фронта фазового перехода, который в описании среднего ее значения представляется изотермической поверхностью, находится с помощью кинетического условия из уравнений:
^плавл — M10> Ф(Х/<2Ja^))/Ф‘a/(^^a17)) = Тплавл^(x,,), (29)
T2нач — (T2нaч — )(1—Ф( x I(^1/°27)) / ‘ 1-ф(а/(^%/027))) = (30)
- Т +о (х, г).
плавл 2
При нулевом коэффициенте диффузии марковского поля и нулевом коэффициенте сноса наблюдается постоянство зоны устойчивости формы фронта фазового перехода.
Репрезентативным представляется результат решения задачи абляции [1] шарового тела с учетом стохастичности поля температур (при ненулевом коэффициенте диффузии марковского поля и нулевом коэффициенте сноса), где закон движения границы раздела фаз имеет вид ^ - а4~г, а форма поверхности есть изотерма. Дисперсия подсчитывается с помощью следующего выражения:
<72(х,г) - Вг+М 2+(м|2) — М2)Ф(х/(2у[аг))/Ф(а/(24аг))—
—{м2+(М1 — м2)Ф(х /(24аг))/ Ф(а/ (^л/Ог)) }2.
Развитие формы поверхности фазового перехода — изотермической сферической поверхности представлено на рисунке.
Поскольку на фронте фазового перехода среднее значение температуры равно
+ж>
M(r>(t,x = ayft) = j ГП(г,x = a4t,T)dT, то зону устойчи-
0
вости формы фронта, который определяется как изотермическая поверхность, можно оценить, решая следующие два уравнения относительно s при данном t:
Ml + (M(1> - M®> • (Ф(* / (2yfat > / фа/ (2yfat >> = Mf> ± oi(t>.
При этом получим два значения s: S( (t> и
s2 (t>, которые и будут определять границы для фронта
фазового перехода в данный момент времени: s(t>е [S((t>;^2(t>] • Из приведенных формул видно, что
зона неустойчивости формы фронта растет со скоростью, пропорциональной Vt, причем коэффициент пропорциональности определяется величиной коэффициента диффузии марковского поля B . Если скорость перемещения фронта меньше или равна скорости роста дисперсии, то следует говорить о неустойчивости формы границы раздела фаз. Это означает, что фронта как такового не существует. Подобная картина определяется статистическим распределением и развитием во времени центров образования новой фазы, а также внешних случайных воздействий, которые и определяют величину B • Однако, если в течение времени t е [0, tmax], когда
проходит технологическая работа агрегата, из-за малости коэффициента B дисперсия медленно возрастает во времени, то решение можно считать пригодным для расчетов в указанном диапазоне времени. В противном случае система даст сбой, обусловленный неучетом стохастичности процесса.
ЛИТЕРАТУРА
1. Карташов Э.М. Аналитические методы в теплопроводности твердых тел. М.: Высшая Школа, 2001. С. 540.
2. Curtain R.F. Stochastical Partial Differential Equations / R.F. Curtain // Stochastical Nonlinear Systems. Berlin: Springer, New York: Heidelberg, (981. 20 p.
3. Мадера А.Г. Стохастическое моделирование процессов теплопередачи в твердых телах // Тепломассообмен-ММФ-92 : труды 2-го Минского международного форума по тепломассообмену. Т.9. Вычислительный эксперимент в задачах теплообмена. Ч.2. Минск: АНК “ИТМО им. А.В. Лыкова” АНБ, 1992. С. ПЫШ.
4. Гардинер К.В. Стохастические методы в естественных науках / К.В. Гардинер. М.: Мир, 1986. С. 528.
5. Соловьев И. А. Уравнения для стохастических полей / И. А. Соловьев. М.: ГУЗ, 2006. 143 с. Соловьев Игорь Алексеевич -
профессор кафедры высшей математики и физики Г осударственного университета по землеустройству, г. Москва
Зубков Павел Валерьевич -
кандидат физико-математических наук доцент кафедры «Высшая математика и физики» Государственного университета по землеустройству, г. Москва
Статья поступила в редакцию 2.08.11, принята к опубликованию 14.11.11
Развитие формы поверхности фазового перехода — изотермической сферической поверхности
