У даній роботі зроблена спроба побудови квазітрехмерной математичної моделі, яка описує режим вібраційного горіння в вертикальної камері горіння регенеративного повітронагрівача доменної печі, розглядаючи її як розподілену динамічну систему.

Анотація наукової статті з фізики, автор наукової роботи - Басок Борис Іванович, Гоцуленко Володимир Володимирович


Mathematical modeling of vibrating combustion mode in a distributed combustion chamber of a blast furnace air heater

This work attempts the construction of a quasi-three-dimensional mathematical model which describes the vibrating combustion mode in a vertical combustion chamber of a regenerative air heater in a blast furnace, considering it as a distributed dynamic system.


Область наук:

  • фізика

  • Рік видавництва: 2018


    Журнал: Математичне моделювання та чисельні методи


    Наукова стаття на тему 'МАТЕМАТИЧНЕ МОДЕЛЮВАННЯ РЕЖІМАВІБРАЦІОННОГО ГОРІННЯ В розподілі КАМЕРЕГОРЕНІЯ воздухонагревателі доменних печі'

    Текст наукової роботи на тему «МАТЕМАТИЧНЕ МОДЕЛЮВАННЯ РЕЖІМАВІБРАЦІОННОГО ГОРІННЯ В розподілі КАМЕРЕГОРЕНІЯ воздухонагревателі доменних печі»

    ?УДК 669.162.23

    DOI: 10.18698 / 2309-3684-2018-4-2540

    Математичне моделювання режиму вібраційного горіння в розподіленої камері горіння повітронагрівача доменної печі

    © Б.І. Басок, В.В. Гоцуленко

    Інститут технічної теплофізики НАН України, Київ, 03057, Україна

    У даній роботі зроблена спроба побудови квазітрехмерной математичної моделі, яка описує режим вібраційного горіння в вертикальної камері горіння регенеративного повітронагрівача доменної печі, розглядаючи її як розподілену динамічну систему.

    Ключові слова: автоколивання, тепловий опір, напірна характеристика теплопідводу, нестійкість

    Вступ. Багато елементів теплоенергетичного обладнання є потенційно нестабільними, і при зміні режимів їх роботи виникають автоколивання, що створюють різні аварійні ситуації. Зокрема, в металургійному виробництві для нагріву доменного дуття використовують регенеративні повітронагрівачі (Каупер). Відомо, що з підвищенням їх теплового навантаження створюється значний економічний ефект. Однак з підвищенням теплового навантаження зростають коливання тиску, через які виникають зусилля близько 20 т, що діють на конструкцію 5-7 разів на секунду [1].

    Введення феноменологічного запізнювання процесу горіння зіграло визначну роль в теорії вібраційного горіння і зумовило механізм нестійкості горіння [2]. Рішення проблеми порушення автоколебаний теплопідводу (конвективним або при згорянні палива) в роботі [2] звелося до обгрунтування іншого механізму, що породжує висхідну гілку на напірної характеристики теплопідводу F (Q). При теплопідводу до потоку в області переходу ламінарного режиму в турбулентний на графіку залежності гідравлічних втрат по довжині hf (Q) утворюється спадна гілка,

    де виконується нерівність dht (Q) / dQ < 0 [3]. В області теплопод-вода з-за зниження щільності потоку змінюється його швидкість. Виник в результаті цього місцеве тепловий опір h (Q), як обгрунтовано [3], також має спадну гілку. Таким чином, низхідні гілки опорів ht (Q) і hr (Q) породжують висхідну гілку на напірної характеристики теплопідводу F (Q) і сос-

    тавляют механізм збудження як термоакустичних автоколебаний, так і вібраційного горіння.

    Механізм, обумовлений тепловим опором, діє локально в області теплопідводу, і при його описі можна використовувати математичні моделі з зосередженими параметрами. Однак, так як Каупер має значну висоту (близько 30-50 м), то в даній роботі для вивчення гідродинаміки в камері горіння Каупера розглядається її дискретно-розподілена модель (рис. 1). У нижній частині камери горіння в зоні теплопідводу розташований коливальний контур з зосередженими параметрами, де діють механізми запізнювання згоряння Л. Крокко і негативного теплового опору, збуджуючи автоколебания. Дані коливання через кордон 1-1 передаються у верхню частину камери горіння - динамічну систему з розподіленими параметрами.

    Мал. 1. Схема умовного поділу камери горіння:

    0-1-1-0 - коливальний контур з зосередженими параметрами, всередині якого розташована зона теплопідводу; 1-2-2-1 - коливальний контур з розподіленими параметрами

    З появою ефективних алгоритмів чисельного інтегрування [4] і методів CFD (computational fluid dynamics) моделювання став можливим чисельний аналіз багатьох завдань гідрогазодинаміки. Однак розрахунок автоколебаний за прямими рівнянням гідрогазодинаміки навіть в цьому випадку пов'язане зі значними

    труднощами. При інтегруванні рівнянь гідродинаміки і теп-ломассопереноса стосовно розподіленої частини камери горіння Каупера ми використовуємо метод розщеплення [5] і його кінцево-різницеву реалізацію, грунтуючись на результатах робіт [6-7].

    Динаміка в зосередженої зоні теплопідводу камери горіння. В області теплопідводу (див. Рис. 1) виникає тепловий опір, яке для довільного політропної підведення теплоти з показником політропи п можна описати наступною залежністю [3]:

    к = п | ^ 1 р0с, (Т1 - Те) + Ро

    1 - п

    2

    ГТ Л Те

    до т1 У

    п-1

    З "V

    до 51 У

    Система рівнянь, що визначає рух в контурі з дискретними параметрами, має такий вигляд [1-3]:

    = Р (2Р ^ - й (, -т) -, (Р),

    (1)

    де Ьа, Са - акустичні параметри коливального контуру; Р (й) = Рпод (й) + А (й) Кт (й) - напірна характеристика теплопідводу; Рпод (й) - напірна характеристика нагнітача, що подає повітря в камеру горіння; й = - об'ємний витрата продуктів згоряння на виході з зосередженого коливального контуру; ф (Р) - функція, що визначає характеристику мережі кс (й) = кдрй2,

    підключеної до коливального контуру; т - запізнювання згоряння палива; I - мінлива часу.

    Наведемо систему рівнянь (1) до безрозмірного вигляду, вважаючи параметри стаціонарного режиму

    де до

    ін

    й = § ​​^ Р = кдр§,

    Рпод (5) + А (§) -до т (§)

    §2 .

    Вводимо безрозмірні змінні: г1 = й / ї * ° тй, z2 = Р / Р ° ° ТР,? ' = Т ^ і т '= т {т. У нових змінних система (1) матиме вигляд

    2

    0

    1

    і

    з

    Т1Т2 т

    тгт1

    т

    (1г

    = Т2 ^ під (т1) + т2 ^ (т) - т2 Ч (т1) - 22;

    =? 1 (ґ -X ') - Г22 (? Ін т2тГ2) "/ 2.

    (2)

    Під дією механізмів запізнювання згоряння палива Л. Крокко і негативного теплового опору [3] в контурі з зосередженими параметрами самозбуджується автоколебания. Позначимо через г * (г) і г * (г) - періодичне рішення системи (2), яке описує ці коливання. Дані функції далі використовують при формуванні граничних умов для коливального контуру з розподіленими параметрами (див. Рис. 1).

    Динаміка в розподіленої частини камери горіння. У безрозмірних змінних рівняння руху (Нав'є - Стокса), нерозривності і теплопереносу в наближенні Буссінеска стосовно вертикальної камері горіння мають вигляд [6-7]:

    - + (УУ) V = - (р) + У2 V + ВГ 0ег

    Ег

    Ш1у У (х, г) = 0,

    ее ^

    - + yo1У

    Ег

    ЄУ - Рг-1 gгad (е)

    : Від'їзд "

    (3)

    (4)

    (5)

    де ВГ =% у [Т1 - Т2] 13 / у - число Грасгофа; Рг = у / а - число Пранд

    тля.

    Розглянуту задачу будемо вирішувати чисельно за допомогою явної триетапної схеми розщеплення, діскретізіруя тимчасову вісь

    г

    п + 1

    | Гп + Аг (п + 1) Аг:

    У = уп-х (- (УПУ) У + У2Уп + вгее,);

    е = єп -X div {епУп - Рг-1gгad (єп)} + Ху пе ,;

    Арп = У УI х; У п + 1 = У-х gгad (рп); єп + 1 = єп -хШ1у {еУп - Рг-1gгad (е)} + хУпег

    (6)

    (7)

    (8)

    Дана схема являє собою комбінацію схеми О.М. Біло-Церковський розщеплення по фізичним факторам для рівнянь гідродинаміки [7] і перерахункових різницевої схеми Н.І. Микита-ко [8] для рівняння конвективного переносу концентрації газу і об'єднує їх переваги. Турбулентний характер руху далі буде враховуватися введенням ефективного коефіцієнта в'язкості [6].

    Запишемо векторне рівняння (3) в циліндричній системі координат в припущенні циліндричної симетрії (т. Е. Е / Еф = 0). Компоненти швидкості середовища в цій системі мають такий вигляд:

    V = іет + Мег.

    Отже,

    ЕІ / ЕГ + (VV) і = -ер / ЕГТ + Ді - і / г2;

    Ем / ЕГ + (VV) м = -ер / ЕГ + Дм + ВГ 0,

    де дії операторів VV і V2 визначаються формулами

    (VV) = ІЕ / Ет + МЕ / ЕГ;

    V2 = IЕ Г т 11 + .Е1

    т Ет до Ет) ЕГ

    2

    Рівняння (5) в циліндричній системі координат буде мати наступний вигляд:

    Е0 1 Е Г Г "-1 Е0Ц Е I." -1 Е01 - + -! Т \ 0і - Рг 1 + - Рг 1 - > = м.

    Ег т Ет {До Ет) \ Ег {Ег]

    Таким чином, в циліндричній системі координат схема розщеплення (6) - (8) має вигляд

    і = і + т

    1 Е (2 ^

    ---1 ти

    т Ет

    {Ті2) |

    Е (їм) Е Г1 Е (ти) ЕГ Ет до т Ет

    1 Е 2и + -

    Ег2

    мм = м + т

    1 Е Г Ем 1 Е2

    Ем 1 Е, ч 1 Е {Ем 1 Е "м ^

    ------ (тім) + - \ т- I + - + Се 0

    Ег т Ет 1 т Ет До Ет) Ег2 е

    (9)

    (10)

    0 - 0п 1 Е

    -+--!

    т т Ет

    0піп

    1 Е0п

    Рге / Ет

    + - \ 0пмп Ег

    1 Е0п рге ^ Ег

    = м

    (11)

    Г

    1 _е_

    г Ег

    г

    V

    п + 1

    Ер!

    Ег

    Е2 рп -1 11 е

    +

    Ег '

    = х

    г Ег

    V Ег У

    +

    Ей1

    ел;

    и- = й-хЕрп / ЕГ, 1 + 1 = 1-хЕрп / ЕГ;

    єп + 1 -еп 1 е

    -+--<!

    х г Ег

    е йп -

    еел

    РГГ Ег

    У

    + А] е 1 »- ± _ Е?

    Ег I Рг ^ ег

    п

    = й

    (12)

    (13)

    (14)

    Для отримання кінцево-різницевого алгоритму на основі розглянутого вище методу розщеплення скористаємося рознесеною (шахової) сіткою [6], наведеної на рис. 2.

    Мал. 2. Схема даної розрахункової області і шаблон для розрахунку швидкості й ;. ,, швидкості, тиску р; | І температури е; ,

    I, J;, J J;, J

    Отримуємо наступні кінцево-різницеві рівняння:

    Аг

    й; 1 = й; 1 + х

    (Йй) ,,, - (йй)

    Аг

    /

    1

    л

    /

    + |

    1 +

    V 2; у

    1

    й; - + 1,1-

    1

    1 -

    V 2; у

    й2-1,1

    (; + 3/2) Аг

    ((; + 3/2) й, - + 1,1 - (; +12) й, -.)-

    (15)

    (; +12) Аг '

    ((; +12) й, -. 1 - (; -12) щ, -1,1) +

    + 7 ^ 2 (й; \ 1 + 1 2й; 1 + й; \ 1 -1)

    Аг

    г

    2

    1

    мм. ., = М- |'|

    г,] г,]

    т (. ..2

    2Дг

    (/ + 1 -] '- 1) |

    + т !

    (Г +1/2) Дт 1

    (Г +12) (їм). ,, - - (г -1/2) (їм).

    г -1, /

    +

    (Г +12) ДТ2

    (Г + 3/2) (мг + 1,] '- м /) - (г + V2) (м, / - м-1, /)

    (16)

    +

    +

    ДГ

    1 г - | 0 "] +1 +0"]

    -I м- 1 - 2м- • '+ м- |> 1 I + ---- ВГ

    2 I г,] +1 уу1,] ^ уу1,] -1 ^ 2

    (Їм) л / = 4 (и- + і; '] + 1) (м- + м; +1,]); (17)

    0;, / = 0п - + м /

    т Е

    (; +12) Дт Ет

    0 "і"

    1 Е0 "

    К / Ет)

    ; > ]

    - Д- 1 Е0 "

    Ег

    Рг / ЕГ

    ;, /

    (Г - 1,5) ДТ2

    (Г-1) (Р "+1, / - Р" /) - (г-2) (р "- - р і, -)

    V (РП / + 1 - 2 РП / + РП / -1) =

    +

    Дг2

    1 {(Т-г ^ кг-1)] - (г-2) 1, / |

    +

    ДГ

    г, / г, / -1

    (18)

    (19)

    п + 1 ~ т

    ь] "г,] Дт

    і ^ = і '; - (РП + м / - РП /); <+ = - (РП / + 1 - РП /); (20)

    0 "+1 = 0" + тм "

    г,] г,] г,]

    т Е

    (Г +12) Дт Ет

    0 ип

    1 Е0

    Рге / Ет)

    г,]

    Е I м п 1 Е0

    - т-! 0мп --

    Ег

    Рге / ЕГ

    (21)

    г> /

    де г '= г +12; / =] +1/2.

    т

    т

    1

    1

    т

    Чисельне моделювання турбулентного характеру руху. При моделюванні обліку турбулентного режиму руху в рівняннях руху скористаємося результатами роботи [6]. Кінематична в'язкість газу і нагрітого повітря порівняно невелика. Рух в камері горіння має турбулентний характер, так як число Рейнольдса Яе ~ 1 / V досить значно, а саме воно в першу чергу визначає режими течії (ламінарний або турбулентний) та динаміку переходу між ними.

    Таким чином, адекватність математичної моделі в значній мірі залежить від правильного обліку в даних умовах ефектів турбулентності. Як зародження [6] турбулентності, так і нестійкість рішення різницевого рівняння мають одну і ту ж природу. Складові, що вводяться в різницеві рівняння для забезпечення їх стійкості або автоматично з'являються при деяких способах апроксимації, певним чином моделюють вплив усереднених в масштабах осередку розрахункової сітки флуктуації фізичних величин, т. Е. Турбулентність.

    З виникненням турбулентності завдяки зародженню мікровіхрей різко посилюється перемішування рідини, і ефективна в'язкість ті /, що розраховується як коефіцієнт пропорційності

    між тензором напружень і тензором швидкостей деформацій для усередненого макрорухи, зростає. В цьому випадку вона має дві складові (т ^ = т + тт), що володіють різною фізичною природою: перша забезпечує молекулярний перенос імпульсу, друга - конвективний перенос імпульсу мікровіхрямі. Ясно, що в разі ламінарного руху ефективний коефіцієнт динамічної в'язкості збігається з коефіцієнтом молекулярної в'язкості т. Широке коло моделей турбулентності ґрунтується на моделюванні коефіцієнта турбулентної в'язкості. Найпростішими моделями є алгебраїчні, коли коефіцієнт тт моделюється деяким

    алгебраїчним виразом, пов'язаних з різними (усередненими) характеристиками руху середовища.

    Для обліку турбулентності в різницеві рівняння необхідно вводити додаткові параметри. При реальному русі такі параметри «з'являються» автоматично з утворенням турбулентності. У розрахунок їх необхідно закладати заздалегідь, а відсутність цих параметрів як раз і призводить до розходження чисельного методу. «Фізична» причина такого роду розходження і навіть її необхідність стає очевидною, якщо врахувати, що дрібномасштабні турбулентні вихори відбирають енергію у великомасштабного (по відношенню до вихором) руху. Якщо ж цю спад енергії не враховувати, потоки розганяються до нескінченних швидкостей, що і стає

    проявом розбіжність розрахункової схеми. Наведені міркування вказують також на те, що турбулентність необхідно враховувати саме в дисипативний члені, що описує спад енергії макрорухи (наприклад, моделюючи коефіцієнт в'язкості).

    Кордон стійкості розрахункової схеми руху рідини на сітці з характерним розміром розрахункової осередку А [6] відповідає початку зародження турбулентних вихорів масштабом менше А (які неможливо описати на даній сітці). Починаючи з цього моменту мікровіхрі, відбираючи енергію у макрорухи, виражаються в збільшенні ефективної в'язкості рідини. Її величина відповідає кордоні стійкості розрахункової схеми.

    Таким чином, параметри ефективної в'язкості необхідно підбирати таким чином, щоб розрахункова схема з такою в'язкістю перебувала на межі стійкості. Це відповідає врахуванню при даних параметрах руху (швидкості і розмірах області, в яких воно відбувається) вихорів, що мають розміри менше розмірів осередку розрахункової сітки.

    Припущення про рівномірної турбулізації рідини по всьому об'єму, яку висувають в моделі з постійною ефективної в'язкістю, в багатьох випадках є занадто грубим. Розглянемо модель турбулентності, асоційовану з сітковим числом Рейнольдса

    ReA = VaAv- !. Тут VA - швидкість в межах даного осередку, V т = тт / Р - турбулентна складова кінематичної в'язкості. Оскільки саме число Рейнольдса визначає перехід в турбулентний режим руху рідини, можна припустити, що сіткове число Рейнольдса у всьому об'ємі рідини постійно. Це призводить до наступної моделі турбулентної в'язкості, лінійно залежить від швидкості [6]:

    V т = VaA ReA. (22)

    З урахуванням (22) ефективні кінематична в'язкість, числа Грасгофа і Прандтля визначаються формулами

    Vef = V + VaDReA;

    Grf = gy (T - r2) I3 (V + Va ^ a) "1; (23)

    Pref = 1 (v + VaD ReA). a

    У модель турбулентності (23) входять три параметра: ReA, V і А. При фіксованій сітці (А = const) в області, де VA ® 0, ефективна в'язкість збігається з молекулярної, т. Е. Vef ® V. В розрахункових

    осередках, в яких велике значення сіткової швидкості Уд, турбулентна в'язкість може у багато разів перевищувати молекулярну [7].

    При записи різницевих рівнянь гідродинаміки на рознесеною сітці необхідно на гранях розрахункової осередку задавати значення в'язкості УР, - і у /, -, числа Грасгофа О /., = 2gу (Г1 - Т2) х

    хI3 (УР, / + У /, /) і Прандтля Рг /.-- = УР, - / а. Значення коефіцієнтів в'язкості в кутах розрахункової сітки УР - 'і у / - > визначають як середнє арифметичне з сусідніх коефіцієнтів УР, - і у /, -: У - = (у - + Ур, - + 1) Д Уг, - = (у / - + Уг + 1, - у 2. Сіткове число Рей-

    нольдса вибирають за допомогою численних експериментів. Число Яед (при фіксованій розрахункової осередку Д) підбирають таким чином, щоб забезпечити розрахунковий процес на кордоні стійкості. Велика кількість чисельних експериментів, проведених для різних течій і сіток (з різними значеннями Д), дозволяє зробити висновок, що при Яед »2 забезпечується стійкість розрахунку гідродинамічних характеристик і це значення знаходиться в області зриву стійкості [6].

    Розглянута вище модель відповідає перенесенню вихору масштабу Д в напрямку швидкості V. Однак ця модель не передбачає умов зародження вихору. Вираз для турбулентної в'язкості, яке враховує умови зародження турбулентних вихорів, було запропоновано Л. Прандтлем:

    у т = I2 Еу |, (24)

    де I - довжина перемішування; у - змінна, параметрізірую-щая вісь, перпендикулярну напрямку потоку.

    Запис на сітці вираження для ефективної в'язкості повинна враховувати її тензорний характер, т. Е. Той факт, що в горизонтальному і вертикальному напрямках в'язкості задаються в різних точках:

    Ур, - = У +12

    у г, - = у +

    м • '- м-' 1

    ;, - -, - -1

    'Дт;

    (25)

    ДГ,

    і; - - і; '- 1, -

    де = (м -, - + у2; Щ = (ur;] +) / 2.

    Довжини перемішування вважають пропорційними сторонам розрахункової осередки: I г = ЬРДт і Іг = ЬгДг, де, Ьг - параметри моделі.

    Отже, остаточно отримуємо:

    ,/ = П + Р2АгЬ + 1, / + / -1 - Ч / - Ч / (26)

    Ч / = п + Р2Аг | мГ; / + 1 + їм>/ + 1-й (// - "М>// 2. (27)

    Дана модель турбулентної в'язкості - двопараметричного. Параметри і Ь підбирають з умови наближення до порога

    нестійкості розрахункової схеми.

    Об'єднуючи розглянуті вище моделі турбулентності, отримаємо трипараметричної модель, яка враховує як перенесення вихору уздовж потоку, так і умови його зародження. В цьому випадку вираз для ефективної в'язкості має вигляд

    V = п + УА А Кеа +12 | Еу / Ц • (28)

    Розстановка початкових і крайових умов (див. Рис. 2). На кордоні АВ маємо наступні граничні умови (г = 0; иг):

    "Ло = 0, Чо '=, Рг, 0 =? 2 * (Ь), 0г, про = 0 *;

    на кордоні СБ (/ = 0; N):

    , / = 0, ^, / '= 0 (умова прилипання); на кордоні АС (/ = 0;):

    "0", / = 0 (умова непротеканія), / = Ч / (умова вільного ковзання); на кордоні СБ (г = 0; иг):

    "Ґ, = 0, ^ - 1 = 0.

    Граничні умови для тиску в вузлах, розташованих на кордоні Г = [ВБ] і [СБ] і [АС], визначають за допомогою рівнянь руху в комбінації з граничними умовами для швидкості і температури. Вони являють собою умови Неймана Ер / Еп | Г = у, де п - нормаль до кордону Г, а у - відома величина, що отримується проектуванням рівнянь руху на кордон Г. Для температу-

    ри 0 на кордоні Г '= [ББ] і [АС] застосовуємо умова повної теплоізоляції Е0 / Ея | г = 0.

    Початкові умови припускаємо однорідними і обумовленими стаціонарним режимом в першому коливальному контурі з зосередженими параметрами (r = 0; иг, 1 = 0; N г):

    при п = 0 <, = 0, = Х / ^ 0П1 = 0, р ^ = (X).

    Рішення рівняння Пуассона для тиску. У розглянутій явною схемою розщеплення (9) - (14) основні труднощі укладена в рішенні рівняння Пуассона (12) на кожному часовому шарі. При записи схеми розщеплення на рознесеною сітці це призводить до вирішення системи кінцево-різницевих рівнянь алгебри (19). Існує досить велика кількість методів, як точних, так і наближених, вирішення цієї системи [9, 10]. Розглянемо метод послідовної верхньої релаксації (ПВР), попередньо записавши рівняння (19) у вигляді

    Рп = АПРП + 1, J + БПРП-1, J + СПРП, ^ + Щ.м + ^;

    ап = - -

    11

    г х

    2г - 3

    (Г - 1,5) Аг2 Аг2

    -1

    х

    х

    {(^ Дт ^ -1) '1 - (г - 2) -1,1] + ^ [1 - 1 -1]

    А

    г -1

    (Г - 1,5) Аг2

    Бп = г-2 1 (г - 1,5) Аг2

    сп _ бп _ 1

    4 = Б = а?

    +-

    (Г -1,5) Аг2 Аг2

    2г - 3 2

    -1

    -1

    -1

    (Г -1,5) Аг2 Аг2

    Метод ПВР [10] визначається наступною итерационной схемою

    (5 = 17):

    "П, 1 (л" \ _п, 1-1. "\ ЛП п, 1-1. Т>п п ,, 1. / ~ Т _п, 1-1. пп "п, 1. / -лп1

    рп = (1 - Цр) рп + Цр {Ач рг + 1+ щр ^ і + + Бпр «,) - 1+},

    2

    де Др - параметр релаксації, що обумовлює швидкість сходи-

    мости ітераційного процесу; 5 - номер зупинки обчислень. Параметр 5 * знаходять з умови

    тах

    I,}

    р г - р

    де е - задана точність.

    Вважаючи = 1, метод ПВР перетворюється в метод Зейделя [10].

    При реалізації обчислень на ПЕОМ, враховуючи їх сучасні обчислювальні можливості, для чисельного рішення рівняння Пуассона (12) більш простим у застосуванні є метод простої ітерації, що дає задовільні рішення при істотно

    більшій кількості ітерацій (5 = 1; 5

    "П, 5" І, 5 1. Г ,<І_І, 5-1. Г) П "П, 5 - 1. / -1П "П5-1. Г »П" Я, 5-1 .

    Р {/ = РП + ®р {АцР + і + ВцР1-й + сцРй + 1 + вцРй-1},

    де ю Р - параметр збіжності методу простої ітерації.

    Результати чисельного моделювання. У розглянутій дискретно-розподіленої моделі камери горіння (див. Рис. 1) є два пов'язаних між собою коливальних контура. У зоні теплопідводу через негативного опору і під дією

    0,2 {

    Мал. 3. Характер деформації граничних циклів і відповідних їм автоколебаний з ростом х запізнювання згоряння палива: а - X = 0; б - X = 0, 25 - релаксаційні коливання; в - X = 0; г - X = 0,002 - коливання, близькі до гармонійних

    *

    інших механізмів нестійкості, наприклад запізнювання згоряння палива, самозбуджується автоколебания, що моделюються зосередженої динамічною системою. На рис. 3 наведено форми автоколивань вібраційного горіння і відповідні їм граничні цикли, самозбудні в зоні теплопідводу (в коливальному контурі з зосередженими параметрами).

    Через відповідне граничне умова Діріхле нестаціонарності передаються в коливальний контур з розподіленими параметрами (рис. 4), де динаміка описується рівняннями Нав'є - Стокса (3-5).

    Мал. 4. Поле швидкостей руху продуктів згоряння в розподіленої частини камери горіння на етапі її заповнення

    Висновок. Розроблено квазітрехмерная дискретно-розподілена математична модель, яка дозволяє виконати чисельні розрахунки гідродинамічних процесів, що протікають в камері горіння регенеративного повітронагрівача доменних печей. Чисельні експерименти, проведені на основі отриманої математичної моделі, дозволяють моделювати режим вібраційного горіння в камері горіння повітронагрівача як динамічній системі з розподіленими параметрами.

    література

    [1] Басок Б.І., Гоцуленко В.В. Розрахунок параметрів автоколивань у вертикальній камері горіння повітронагрівача доменної печі при нестійкому горінні. Теплоенергетика, 2015 року, № 1, с. 59-64.

    [2] Гоцуленко В.В. Математичне моделювання особливостей феномена Рійк. Математичне моделювання, 2004, т. 16, № 9, с. 23-28.

    [3] Басок Б.І., Гоцуленко В.В. Термогидродинамические нестійкість потоку теплоносія. Київ, ТОВ ВД «КАЛИТА», 2015 року, 412 з.

    [4] Дімітріенко Ю.І., Коряков М.Н., Захаров А.А. Застосування методу Яків для чисельного рішення рівнянь газової динаміки на неструктурованих сітках. Математичне моделювання та чисельні методи, 2015 року, № 4, с. 75-91.

    [5] Марчук Г.І. Методи розщеплення. Москва, Наука, 1988, 264 с.

    [6] Самохвалов С. Є. Теплофізичні процеси в багатофазних середовищах: теоретичні основи комп'ютерного моделювання. Київ, Вид-во Інституту системних досліджень. Мін. обр. Україна, 1994, 172 с.

    [7] Білоцерківський О.М. Чисельне моделювання в механіці суцільних середовищ. Москва, Наука, 1984, 520 с.

    [8] Нікітенко М.І. Парні і зворотні завдання тепломассопереноса. Київ, Наукова думка, 1988, 240 с.

    [9] Роуч П. Обчислювальна гідродинаміка. Москва, Мир, 1980, 616 с.

    [10] Поттер Д. Обчислювальні методи в фізиці. Москва, Мир, 1975, 392 с.

    Стаття надійшла до редакції 10.10.2018

    Посилання на цю статтю просимо оформляти наступним чином: Басок Б.І., Гоцуленко В.В. Математичне моделювання режиму вібраційного горіння в розподіленої камері горіння повітронагрівача в доменній печі. Математичне моделювання та чисельні методи, 2018, № 4, с. 25-40.

    Басок Борис Іванович - д-р техн. наук, проф., чл.-кор. НАН України, завідувач відділом теплофізичних основ енергоощадних технологій інституту технічної теплофізики НАН України. e-mail: Ця електронна адреса захищена від спам-ботів. Вам потрібно увімкнути JavaScript, щоб побачити її.

    Гоцуленко Володимир Володимирович - д-р техн. наук, провідний науковий співробітник відділу теплофізичних основ енергоощадних технологій Інституту технічної теплофізики НАН України. e-mail: Ця електронна адреса захищена від спам-ботів. Вам потрібно увімкнути JavaScript, щоб побачити її.

    Mathematical modeling of vibrating combustion mode in a distributed combustion chamber of a blast furnace air heater

    © B.I. Basok, V.V. Gotsulenko

    Institute of Engineering Thermophysics of National Academy of Sciences of Ukraine,

    Kyiv, 03057, Ukraine, 03057

    This work attempts the construction of a quasi-three-dimensional mathematical model which describes the vibrating combustion mode in a vertical combustion chamber of a regenerative air heater in a blast furnace, considering it as a distributed dynamic system.

    Keywords: self-oscillations, thermal resistance, pressure characteristic of heat supply, instability

    REFERENCES

    [1] Basok B.I., Gotsulenko V.V. Teploenergetika - Thermal Engineering, 2015-го, vol. 62, то. 1, pp. 58-63.

    [2] Gotsulenko V.V. Matematicheskoe modelirovanie - Mathematical Models and Computer Simulations, vol. 16, no. 9, 2004, pp. 23-28.

    [3] Basok B.I., Gotsulenko V.V. Termogidrodinamicheskaya neustoychivost potoka teplonositelya [Thermohydrodynamic instability of coolant flow]. Kiev, TOV VD KALITA Publ., 2015-го, 412 p.

    E.H. EacoK, ​​B.B. ro ^ nerno

    [4] Dimitrienko Yu.I., Koryakov M.N., Zakharov A.A. Matematicheskoe modelirovanie i chislennye metody - Mathematical Modeling and Computational Methods, no. 4 (8), 2015-го, pp. 75-91.

    [5] Marchuk G. I. Metody rasshchepleniya [Methods of splitting]. Moscow, Nauka Publ., 1988, 264 p.

    [6] Samokhvalov S.E. Teplofizicheskiye protsessy v mnogofaznykh sredakh: teore-ticheskiye osnovy kompyuternogo modelirovaniya [Thermophysical processes in multiphase media: theoretical foundations of computer simulation]. Kiyev, Institut sistemnykh issledovaniy Publ., Min. obr. Ukraine, 1994, 172 p.

    [7] Belotserkovskiy O.M. Chislennoye modelirovaniye v mekhanike sploshnykh sred [Numerical modeling in continuum mechanics]. Moscow, Nauka Publ., 1984, 520 p.

    [8] Nikitenko N.I. Sopryazhennyye i obratnyye zadachi teplomassoperenosa [Conjugate and inverse problems of heat and mass transfer]. Kiev, Naukova dumka Publ., 1988, 240 p.

    [9] Rouche P.J. Computational Fluid Dynamics. Hermosa, 1972. [In Russ .: Rouche P.J. Vychislitelnaya gidrodinamik. Moscow, Nauka Publ., 1980, 616 p.].

    [10] Potter D. Computational Physics. John Wiley, 1973. [In Russ .: Potter D. Vychis-litelnyye metody v fizike. Moscow, Nauka, 1975, 392 p.].

    Basok B.I., Dr. Sc. (Eng.), Professor, corresponding member of National Academy of Sciences of Ukraine, Head of Department of Thermophysical Bases of Energy-Saving Technologies, Institute of Engineering Thermal Physics of National Academy of Sciences of Ukraine. e-mail: Ця електронна адреса захищена від спам-ботів. Вам потрібно увімкнути JavaScript, щоб побачити її.

    Gotsulenko V.V., Dr. Sc. (Eng.), Senior Research Fellow, Department of Thermophysi-cal Bases of Energy-Saving Technologies, Institute of Engineering Thermal Physics of National Academy of Sciences of Ukraine. e-mail: Ця електронна адреса захищена від спам-ботів. Вам потрібно увімкнути JavaScript, щоб побачити її.


    Ключові слова: автоколиваннями /Теплове ОПІР /Напірної характеристики теплопідводу /нестійкість /SELF-OSCILLATIONS /THERMAL RESISTANCE /PRESSURE CHARACTERISTIC OF HEAT SUPPLY /INSTABILITY

    Завантажити оригінал статті:

    Завантажити