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

Работа выполнена при финансовой поддержке гранта РНФ №23-21-00269,

Как известно из теории переходных процессов в сетях энергосистемы, токи короткого замыкания достигают сверхкритических значений, а с учетом закона Джоуля — Ленца не сложно дать оценку теплового воздействия указанных токов, которые могут привести к катастрофическим последствиям в энергосистеме. Для защиты от данных явлений в энергосистеме предусматриваются силовые коммутационные выключатели, задача которых максимально в короткое время разорвать поврежденную цепь.

Однако, даже если контакты выключателя разомкнуты, ток короткого замыкания может восстановить электрическую дугу за счет явления плазмодугового разряда между контактами. На этот случай в выключателе, в зависимости от типа и конструктивных особенностей, имеется механизм и средства гашения упомянутой дуги, позволяющие резко сбросить температуру на оси дуги с 7000 K до примерно 2000 K и обеспечивающие надежное отключение.

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

Для устранения этого парадоксального явления Дж. К. Максвелл в рамках теории газодинамики, А. В. Лыков в процессах тепломассообмена, К. Каттанео и П. Вернотте при исследовании теплопроводности из молекулярно-кинетических представлений, используя гипотезу о конечности времени соударения молекул и длины свободного пробега, получили новый обобщенный закон тепломассопереноса, в котором фигурирует дополнительный элемент, так называемое время релаксации, т. е. время установления термодинамического равновесия между потоком и градиентом потенциала. При этом выводится уравнение теплопроводности гиперболического типа [1, 2].

• 1    Постановка начально-краевой задачи

Исследуя данное явление, сформулируем математическую модель упомянутого переходного процесса в виде гиперболо-параболического уравнения в частных производных второго порядка с нелинейным источником тепла и переменным коэффициентом теплопроводности:

∂ ² u        ∂u

k ^{( x,t )} dtp + a ^{( x,t )} dt

∂ ∂x

^{X ( x,t )} dx

+ c(x, t)u + q(u),

в области G = [0, L] x [ — T i ,T 2 ] с положительными T i и T 2 - Здесь при t <  0 k(x, t) = 0, а при t >  0 k(x, t) > 0, k(x, t) — коэффициент тепловой релаксации; X(x, t) — коэффициент теплопроводности; c(x, t) — коэффициент теплоотдачи на боковой поверхности; q(u) — нелинейный внутренний источник тепла.

Начально-краевая задача заключается в поиске решения уравнения (1) при выполнении следующих граничных условий третьего рода и заданного начального профиля температуры u o (x):

'TX(x,t) dx

+ a o,L (x,t) u(x,t)

= Г qo ^{( t )} ;

:=o,l    (qL ^{( t}V

u(x, -T i ) = u o (x),                           (3)

где α 0_,L и q 0_,L — коэффициенты теплоотдачи и плотности теплового потока поверхностных источников на левом и правом концах стержня [3].

Для численного расчета поставленной задачи зададим конкретные значения коэффициентов уравнения (1):

k(x, t) = k при t > 0, a(x, t) = a > 0, X(x, t) = x + 2, c(x, t) = c < 0, q(u) = u(x, t)2, тогда уравнение (1) принимает следующий вид:

k^Ugr¹ + a^udx^ = dx [ (x + 2)5^ 1 + cu(x,t) + u(x,t) ² . (4)

• 2    Численное решение поставленной задачи

Проведем дискретизацию по пространственной переменной { x i J Nli , X i = ih, h = L/N и по временной: {T j } M= ^{M 1 + M} 2 , при t <  0 T i = T 1 /M 1 , t > 0 T 2 = T 2 /M 2. ³

Составим конечно-разностное уравнение для произвольного внутреннего элементарного объема ячейки [x i - o.5 , x i+0.5 ] на промежутке времени [t j — i ,t j ] с использованием интегроинтерполяционного метода

• [3] . В таком случае неявная разностная схема будет выглядеть следующим образом:

, ( U j⁺¹ - 2U j + U j ^-1 ) (U^ - U j )                      _Т +^ ,_{т 2} .2

k--------Т2+a-----Т----- = Wi+O^-Wi-O^+cUi +(Ui ) , где xi±o.5 = (xi ± h), Wi±o.5 — потоки тепла, где основанием для расчета является обычный закон Фурье:

Uj+1 - Uj+1 i+1

Wi+0.5 ≈ λi+0.5

Uj+1 - Uj+1 ii

Wi-0.5 ≈ λi-0.5

Эффективные теплопроводности λ i ± 0 _. 5 определим следующим образом:

λ i ± 0 . 5

2 X(x i )X(x i ± i ) X(x i ) + X(x i ± i ).

Используя тот же метод конечных объемов, составляем разностные уравнения для краевых элементарных объемов [3], учитывая при этом тепловые потоки на границе области в окружающую среду, выражения для которых вытекают из краевых условий (2):

j +1      j +1                                               j +1 j

- X 0.5 -1-- h- ⁰- + a U = q o + h ^   + (U 0 ) ² - aU°-^ ---⁰

j +1 j +1

j+1     j cUjj+1 + (UN)2 - aUN Т UN

U N - U N - 1      j +1       h

^X N - 0.5----- h+^a L U N = ^q L ⁺ 2

Полученная система является частным случаем системы канонического вида [3] для метода прогонки:

A i U j + + B i Ui^¹ + C i U j+11 = F i , (i = 2,N - 1).

Соберем в соответствии с каноническими уравнениями:

Г A o . 5 I _п_ ch I ah 1 T rj+¹ Г A o . 5 1 T rj + l _ I h(TT^j \2 i q^hTT^j

["h" ^{+ a} 0 - T ⁺ 2tJ U 0    - [ T" J U 1 = q 0 ⁺ 2 (U 0 ) ⁺ 2t U 0

λ i- 0 . 5      j +1       k a    λ i +0 . 5     λ i- 0 . 5           j +1       λ i +0 . 5      j +1

• - [“h^ J U i-1 ^{+    +} T ⁺ h ^{2 +} h ² c\U i    - [-"h 2-J U i+l =

= (U j ) ² + Tk U j - T 2 U j ^-1 + a U j

• I" A n -0.5"! ^j++l 1  Г A n -0.5  .           ch . ah J j++l         , hfTTj \2 ।  ahjTj

• _^- -- h --JUN - 1 ⁺ [-- h ^{+ a} L ^- "2" ⁺ 2 T_| U N = q L ⁺ 2 ^{( U} N ) ⁺ 2T U N .

• 3 Формирование начальных условий для гиперболического этапа

Первоначально, уравнение (4) соответствует параболическому типу, т. е. при t <  0 k(x, t) = 0, где начальное условие для уравнения (4) имеет вид (3). Однако при 0 < t < T 2 k(x, t) > 0 и уравнение (4) преобразуется в гиперболический вид, для которого требуется найти первое начальное условие u(x, 0) из предыдущего расчета для параболического этапа при t = 0:

u(x i , 0) = u M ¹ , 1 <  i <  N.

Далее требуется определить второе начальное условие для гиперболического уравнения. Для этого указанное условие зададим с повышением порядка аппроксимации по времени до второго, в отличие от первого порядка для формулы:

du(xj, 0) = U i^{M 1} - U M ^{1 - 1} ∂t               τ 2         ^.

Для этого запишем ряд Тейлора на точном решении u(x, t) по времени в окрестности t = 0:

u(x„T2) = u(x„ 0) + Т2 ди(Х12> + Т4 д^10 + О(Т2), dt 2 dt где вторую производную выведем из уравнения (4):

d ² u(x, t) ∂t ²

1 d k ∂x

(x + 2)^du^)] - aSu^ + C (x_it) + 1 u ² (x,t) = ∂x k ∂x k        k

1 du(x,t)

k

∂x

+ (x + 2) ^d 2 ^u 1

∂x ²

- a^dUr-t ⁾ +C"(x,t)+1 uhxA k ∂x k        k

Найдем отсюда неизвестное дискретное значение:

M 1 +1 - M^M1 _T ^{du ( x} i , ⁰⁾ T 2 i         i ^{+ 2} dt ⁺ 2

1 du(x i , 0)    (x i + 2) d ² u(x i , 0)

k ∂x        k      ∂x ²

adu(x i , 0)     c ,         1

-T----- "⁺ ^^{u ( x} i ^{, 0) +} Tu ²⁽ x i , ^{0) + O ( t 3 )}

k ∂t      k          k

т. е.:

M 1     M 1 - 1

U M 1 +1 = U M 1 + T 2 U i - U i---

τ ²

+ —

^U i^M+ ¹ 1 ^{- U} i^M - ¹ 1

2hk   ⁺

(X i +2) (U Ml - 2U M 1 + uMp

⁺ k           h ²

M 1     M 1 - 1

- a ^{( U} i   _TUi   ⁾ + cu M ¹ + W¹ ) ² + o- .

k       τ 1           k        k

В результате получаем:

du(x i , 0) = U i^{M 1 +1} - U MM ¹

∂t             τ 2

M1M

Ui  - Ui+ T ti2

^U i^M+ ¹ 1 ^{- U} i^M - ¹ 1

2hk   ⁺

(X i +2) (U iM1 - 2U M ¹ + U MI) _

⁺ k            h ²

M 1    M 1 - 1

- a ⁽ i ^{- 1} —)+cu ^M 1 ⁺ 1 ⁽ u M ¹ ) 2 + ° ^{< T} 2² ).

k         τ 1              k ⁱ k ⁱ

Отсюда получаем второе начальное условие для расчета гиперболического этапа со вторым порядком аппроксимации:

du(x i , 0) = U M 1 - U M 1 ^{- 1} + T 2 dt               T 1           2

(x i + 2) (U M1 - 2U M 1 + U M1 )

⁺ k           h ²

-

^U i^M+ ¹ 1 ^{- U} i^M - ¹ 1

2hk   ⁺

M 1    M 1 - 1

a ⁽ U i ----Ui----) + cU M i + 1(U M i ) ²

k         τ 1             k ⁱ k ⁱ

Заключение

После определения всех требуемых компонентов поиск решения уравнения (4) с начально-краевыми условиями (2-3) производится методом прогонки по неявной консервативной разностной схеме с учетом переменного коэффициента теплопроводности, нелинейного источника тепла и бокового теплоотвода. Полученные температурные поля хорошо согласуются с имеющимися экспериментальными данными. Приведенный выше алгоритм может быть легко реализован на различных языках программирования, в частности, с использованием инженернопрограммных средств пакетов Mathcad и Comsol Multyphysics. В дальнейшем планируется применить метод итераций для завершения данной нелинейной схемы и распространить ее для двухмерного пространственного случая.