Представлено решение сопряженной задачи о деформации и ресурсе конструкции реакционной камеры вихревой горелки, использующей мелкодисперсную смесь вязкого жидкого топлива и газообразного окислителя. На основе полученных численных результатов проведена оценка гидродинамического воздействия потока жидкости на деформацию элементов конструкции. Показано, что учет этого воздействия позволяет определить напряженно-деформированное состояние в элементах реакционной камеры при функционировании устройства, а также провести оптимизацию конструкции для достижения рациональных значений коэффициента запаса прочности.
Studying coupled processes of hydroelasticity and durability of a centrifugal injector for preparing high viscosity fuel.pdf В настоящее время особое внимание при разработке объектов энергетики уделяется требованиям в области экологичности, энергоэффективности и энергосбережения [1]. Несмотря на развитие технологий альтернативной энергетики, основная часть произведенной тепловой и электрической энергии приходится на долю традиционных углеводородных видов топлива. Одним из распространенных топлив для котельных агрегатов электростанций, промышленных и бытовых печей, судовых котельных установок являются различные виды нефтяного мазута [2 - 4]. Однако методы сжигания с использованием форсунок традиционных конструкций обладают рядом существенных недостатков [4], например высокой засо-ряемостью топливных каналов форсунок вследствие низких эксплуатационных характеристик мазута, значительным количеством вредных выбросов в окружающую среду, связанных с неполным сгоранием сырья. Таким образом, создание и разработка новых устройств для сжигания высоковязких видов котельного топлива с учетом повышенных требований по экологичности и энергоэффективности является актуальной задачей. В данной работе рассмотрена задача о создании равномерного заполнения реакционной камеры вихревой горелки мелкодисперсной смесью топлива и окислителя. Предлагается использовать систему диспергирования, основанную на комбинации пульсационного и аэродинамического методов распыления, когда диспергируемая среда при прохождении через щелевой канал форсунки в режиме пульсации попадает в высокоскоростной завихренный газовый поток окислителя [5, 6]. При этом площадь истечения щелевого сопла предполагается переменной и зависит от воздействия пульсаций давления в смеси с упругим откликом элементов конструкции форсуночной установки. Для разработки и исследования устройств с деформируемыми рабочими элементами необходимо использовать подход, описывающий гидроупругое взаимодействие твердых тел и жидкости [7 - 9]. Таким образом, возникает потребность объединения в рамках единой вычислительной модели подходов механики жидкости и газа (МЖГ) и механики деформируемого твердого тела (МДТТ). Построение таких физико-математических моделей осложняется тем, что система уравнений, основанная на фундаментальных законах сохранения, должна формулироваться в нестационарной постановке для описания одновременно протекающих динамических процессов. Следует отметить, что вследствие специфики пульсационной истории нагружения, для определения ресурса деформируемых элементов конструкции форсуночного устройства необходимо проводить расчет ресурса усталостной долговечности путем введения в систему дополнительных уравнений, полученных из экспериментальных данных по установлению характеристик по сопротивлению усталости при периодическом нагружении. Рассмотрим сопряженную задачу об определении рабочих параметров диспергирования паромазутной смеси для последующего сжигания в реакционной камере вихревой горелки и надежности деформируемых элементов конструкции на примере форсуночной установки с щелевым соплом (рис. 1). Типовая схема устройства включает в себя камеру распределения топлива 1, которая состоит из упругих элементов, представляющих собой две вогнутые пластины 2, закрепленные Рис. 1. Схема конструкции форсуночной установки на полом распределяющем стержне 3 и вихревую камеру смешения 4 с двумя тангенциально расположенными штуцерами подачи окислителя 5. Устройство снабжается уплотняющими элементами в виде резиновых 6 и медных колец 7. Паро-мазутная смесь через штуцер 8 подается в распределяющий стержень 3. При возникновении волнового импульса, создаваемого генератором пульсации давления, вогнутые пластины форсунки 2, образующие камеру распределения 1, испытывают кратковременную упругую деформацию, соприкасающиеся стенки расходятся и формируют щелевое сопло, через которое осуществляется подача диспергируемой смеси в камеру смешения с окислителем 4. При смешении порция топлива приобретет тангенциальную составляющую скорости, после чего поступает в реакционную камеру вихревой горелки 9. Компонентами топливной смеси являются топочный мазут марки М-100 ГОСТ 10585-99 и водяной пар в массовом соотношении 1:1. Упругие элементы конструкции 2 камеры распределения выполнены методом штамповки из конструкционной стали марки 30 ХГС ГОСТ 4543-71. Схема наложения сварных швов по ГОСТ 16037-80. Толщина пластин 1 мм, диаметр вогнутости e = 300 мм, диаметр большего основания камеры распределения с = 100 мм, диаметр меньшего основания 20 мм. Высота камеры смешения a = 200 мм, внешний диаметр камеры смешения b = 150 мм, внутренний диаметр камеры смешения d = 90 мм. Математическая модель Сопряженная задача определения ресурса элементов конструкции форсуночного устройства включает в себя следующие одновременно решаемые задачи: 1. Определение динамических параметров движения жидкости и газа (распределение полей давлений и скоростей) и их влияния на механическое поведение элементов конструкции форсунки. 2. Определение структурного отклика конструкции на периодическую динамическую нагрузку (распределение напряжений и деформаций) и его влияния на рабочие параметры, соответствующие области МЖГ. 3. Определение усталостной долговечности конструкции устройства при циклическом нагружении. Задача решается в переменных Лагранжа. Вводится декартова система координат, в которой ось OZ совпадает с осью полого стержня 3 (рис. 2). Начало координат расположено в точке пересечения оси OZ с плоскостью симметрии камеры распределения топлива 1. Рис. 2. Схема расположения системы координат Обобщенная система разрешающих уравнений для решения нестационарной сопряженной задачи имеет следующий вид [10]: + V(pu) = 0, (1) dt где p - плотность вязкой паромазутной смеси, V - дифференциальный оператор Гамильтона, и - вектор скорости; д -(pv) + V(pvv) = -VP + V(t), (2) dt где P - давление, t - тензор напряжений; Уравнение состояния для жидкой среды описывается зависимостью вида Р = Р (P). (3) Для математического описания турбулентных течений используется модель переноса касательных напряжений (Shear Stress Transport - SST) с соответствующим разрешением пограничного слоя расчетной сетки. Система уравнений описания механического поведения конструкции имеет следующий вид [11]: (4) d 2ui Зг;, j = Р s ГГ" dt2 +^77 5-£kk I , (6) где иг - компоненты вектора перемещений, а-- - компоненты тензора напряжений, ps - плотность материала конструкции; Ч = 2 + "],г ), (5) где е- - компоненты тензора деформаций; Определяющее соотношение записывается в виде E_Г _v_ 2 (1 + v )е,] +1 - 2v где E, v - модуль Юнга и коэффициент Пуассона материала устройства, 5- - символ Кронекера. При формулировании физико-математической модели, описывающей динамические процессы, в основе которых лежит взаимодействие газообразных и жидких сред с податливыми элементами конструкции, необходимо дополнительно учитывать возможность потери несущей способности устройства вследствие многоцикловой усталости [12, 13]. На этапе схематизации истории нагружения происходит переход от реальной истории нагружения к её схематизированному аналогу. Для перехода от текущего многоосного напряженного состояния к экспериментальному одноосному используется гипотеза, основанная на удельной энергии формоизменения, которая описываются следующим соотношением [14]: aeq =72^(а! - а2 ) +(а2 - а3 ) +(а! - ^ ) , (7) где aeq - эквивалентное напряжение, ai, a2, а3 - первое, второе и третье главные напряжения данного напряженного состояния. Для описания перехода от текущего ассиметричного цикла нагружений к симметричному используют зависимость, связывающую амплитуду напряжения эквивалентного по повреждаемости симметричного цикла нагружения aa с амплитудным и средним напряжением ассиметричного цикла нагружения, соответствующего рассматриваемой сопряженной задаче, следующего вида [12]: % = О a + Va0m , (8) где aa = max _-mL - амплитуда напряжений ассиметричного цикла, am = CTmax + a™n- среднее напряжение ассиметричного цикла, ya - коэффициент учета асимметрии цикла, определяемый по формуле [12]: 2a_1 _ О0 Va =-1-Ч a0 где a-1 - предел усталости при симметричном цикле нагружения (индекс «-1» соответствует коэффициенту асимметрии для симметричного цикла R = min = -1), amax a0 - предел усталости при пульсирующем (отнулевом) цикле нагружения (R=0). Для учета масштабных, технологических и геометрических факторов вводится коэффициент пересчета предела усталости Kt и вместо предела усталости a-1 используется предел усталости элемента конструкции aT-1 определяемый по формуле [15] aT-j = KT a_j; (9) KT =-;-1--, (10) (KjKd +1/ Kf - 1)/Kv где Ka - коэффициент концентрации напряжений (отношение предела усталости гладкого образца к пределу усталости образца с концентратором напряжений), определяется по формуле Ka = 1 + ( aa _ 1) q , где q - коэффициент чувствительности материала к концентрации напряжений; из номограмм [12, 13] определяются: о^ - теоретический коэффициент, Kd - коэффициент влияния размеров (масштабный фактор), Kf - коэффициент влияния качества обработки поверхности, Kv - коэффициент влияния поверхностного упрочнения. С учетом (9) амплитуда напряжения aa принимает вид aT 0 a 0a = K-~ . KT Тогда амплитуду напряжения эквивалентного по повреждаемости симметричного цикла нагружения можно записать в виде % = ^ + Vaam . (11) По результатам проведения схематизации истории нагружения производится оценка поврежденности на основе модифицированной гипотезы суммирования повреждений Седлачека [16]: s n О ^ _ 1 _ amax 1 N О В ' где s - число уровней напряжения; щ - число циклов наработки за время эксплуатации при напряжении ста ; N, - число циклов до разрушения при постоянном нагружении; о^ - максимальное напряжение в блоке; oB - предел прочности. В этом случае суммарное число циклов Nc до разрушения выражается формулой: R "max Nc _ (12) с о 1 _где R _ nij ^ nt - относительная длительность действия напряжения oa в блоке. Граничные и начальные условия Обозначим область пространства, занимаемую потоками жидкости и газа Q, а область, соответствующую элементам конструкции устройства, Q2 (рис. 3). Тогда граничные и начальные условия, замыкающие систему дифференциальных уравнений, описывающих нестационарный процесс (1) - (12), записываются в следующем виде. Для задачи, относящейся к движению вязкой среды: В точке x _ {xl, x2, x3}, принадлежащей области подачи в распределительную камеру форсунки распыления Гь задается закон изменения давления: P(x)_A_cos(2пT-jj,x е^ сЦ, (13) где A - амплитуда, T - период; Подача окислителя осуществляется с постоянной скоростью: Ua (Х) _ Ua,0 , (Х) _ Ut,0, X е Г2 C ^ (14) где ua - осевая компонента вектора скорости, ut - тангенциальная компонента вектора скорости; На выходной границе Г3: P (x)_ PC, x еГ3 cQ1, (15) где Pc - давление на входе в камеру сгорания. На внешних стенках установки Г4 задается условие прилипания потока: и (x)_ 0, x еГ 4 сЦ. (16) Для задачи о напряженно-деформированном состоянии в элементах конструкций: Основание форсунки распыления жестко закрепляется на границе Г5: u (x)_ 0, x еГ5 cQ2. (17) Граница Г6 свободна от внешних воздействий: щOj (x) = 0, x еГб CQ2. (18) На граничной поверхности областей Г7: щOj (x) = P(x), x еГ7 CQ1 nQ2; (19) и(x) = 0, x еГ7 cQ1 nQ2. (20) При этом распределение давлений по координатам P (x) получено из решения системы уравнений (1) - (3) на границе Г7. Начальные условия записываются в виде P (x,0) = F1 (x), x eQ1, и (x,0) = F2 (x), x е Q1, T(x,0) = T0, x е Q1, U (x, 0) = 0, x е Q2 , где функции F (x) и F2 (x) определяются из решения стационарной задачи, когда P(x) = P0, x е Q1, T0 = 300 К = const. Для решения системы уравнений вместе с граничными и начальными условиями в работе применялся метод конечных элементов, реализованный в универсальной программной системе ANSYS v.15.0. Для задания выражений, связанных с сопряжением расчетных областей и перестройкой сеточного разбиения области Q1 в результате изменения положения стенок упругих элементов использовались пользовательские функции, написанные на языке APDL и CEL. Рис. 3. Схема расчетной области Г\/ х £2 Qi При выполнении конечно-элементного анализа по достижении невязок значений порядка 10-6 критерий сходимости на временном шаге интегрирования считается выполненным. Временной шаг интегрирования разрешающей системы уравнений определяется из условия Д t _ min(At1, A t2), где At1 и At2 - шаги интегрирования, удовлетворяющие критерию сходимости решения для областей Q1 и В результате исследования сеточной сходимости установлено, что в случае пространственной дискретизации расчетных областей задач МДТТ и МЖГ сеткой конечных элементов в размере 20 000 и 150 000 соответственно дальнейшее измельчение конечно-элементного разбиения не оказывало существенного влияния на результаты расчета. Для верификации программного комплекса использовались задачи, описанные в разделе «Verification Manual» [17]. Данный раздел содержит набор тестовых задач, в которых представлено сравнение результатов конечно-элементного анализа и аналитического расчета. Результаты и обсуждение Для оценки целесообразности решения сопряженной задачи и применения единой вычислительной модели (сопряженный подход - СП) численное решение системы уравнений (1) - (12) вместе с граничными (13) - (20) и начальными условиями сравнивалось с решением задачи теории упругости (ТУ) с упрощенными условиями нагружения. При решении задачи ТУ влияние течения потока жидкости на конструкцию упрощенно учитывается путем задания нагрузки в виде давления, изменяющегося по условию (13) с амплитудой 2105 Па и периодом 2 с, равномерно распределенного на внутренних поверхностях элементов Г7, т.е. решения системы уравнений (4) - (12) с граничными условиями (17) - (20). На рис. 4 показаны расчетные распределения динамического (а), статического (б) и полного (в) давления в момент времени 1,063 с в сечении камеры при заданном законе изменения давления (13) с амплитудой 2105 Па и периодом 2 с. Максимальное значение полного давления наблюдается в зоне распределяющего стержня и составляет 2,04 105 Па. Возникающие зоны разрежения в окрестности камеры распыления формируются в вихревой камере смешения за счет расходящегося потока. Указанные составляющие давления связаны соотношением Ptot _ Pstat + Pdyn _ Pstat + ^P (u ' U) , где Ptot - полное давление, Pstat - статическое давление, Pdyn - динамическое давление. Распределения осевой (а), радиальной (б) и тангенциальной (в) компонент вектора скорости в сечении камеры в момент времени 1,063 с показан на рис. 5. Максимальное значение осевой компоненты скорости составляет 170 м/с. Для нормативных условий функционирования устройства характерно возникновение геометрически обусловленной области застойной зоны, расположенной над камерой распределения. Рост осевой компоненты вектора скорости наблюдается в областях, прилегающих к входу в вихревую горелку, где тангенциальная и радиальная компоненты вектора скорости составляют 66,7 и 12,8 м/с соответственно. 0,87 0,77 0,67 0,58 0,48 0,39 0,29 0,19 0,10 0 Па 170,89 147,70 124,50 101,31 78,11 54,92 31,72 8,53 -14,67 -37,86 м/с а 1,98 1,74 1,51 1,27 1,03 0,79 0,55 0,32 0,08 0 Па 2,04 1,79 1,55 1,31 1,06 0,82 0,57 0,33 0,08 0 Па Рис. 4. Контуры распределения динамического (а), статического (б) и полного (в) давления б 76,71 63,93 51,15 38,37 25,59 12,81 0,03 -12,74 -25,52 -38,30 м/с 100 88,89 77,78 66,66 55,55 44,44 33,33 22,22 11,11 0 м/с Рис. 5. Контуры распределения осевой (а), радиальной (б) и тангенциальной (в) компонент вектора скорости На рис. 6 показаны зависимости эквивалентных напряжений (а) и перемещений (б) элементов конструкции от положения на координатной оси ОХ для моментов времени (0,25; 0,5; 1,063 с) рабочего процесса. Штрихпунктирными линиями показаны результаты решения упрощенной задачи в рамках ТУ. Максимальные значения эквивалентного напряжения на всем этапе нагруже-ния находятся в упругой области функционирования материала и расположены на расстоянии 13,5 мм по оси Х, что соответствует области крепления упругих элементов 2 к распределяющему стержню 3. Наибольшие значения перемещений наблюдаются в области щелевого канала форсунки. Разница между результатами, полученными в рамках СП и ТУ, достигает 8,5 % для значений эквивалентных напряжений (265 МПа - СП; 288 МПа - ТУ) и 16,5 % - для значений перемещений (0,132 мм - СП; 0,154 - ТУ). Расчет напряжений и фактических коэффициентов запаса производился по ГОСТ 52857.1-2007 «Нормы и методы расчета на прочность». Для учета влияния многоцикловой нагрузки в качестве предельного для данного материала напряжения используется предел длительной прочности. Коэффициент запаса по пределу длительной прочности в соответствии с нормативным документом принимался равным 1,5. Фактический коэффициент запаса составляет 1,7 (получен в рамках СП) и 1,6 (при решении задачи в рамках ТУ). я я * ts Л я й я ^ о я н я ^ ч й я я я и 300 п 3 200 100- 0 10 20 30 Ось Х, мм 40 50 m 0,16п 3 б И 0,08- 0,04 20 40 Рис. 6. Распределение эквивалентных напряжений (а) и перемещений (б) по оси Х для некоторых моментов времени (1 - 0,25; 2 - 0,5; 3 - 1,063 с) 0,12 0 10 30 Ось Х, мм 50 Зависимости максимальных эквивалентных напряжений (а) и перемещений (б) элементов конструкции от времени представлены на рис. 7. Функции распределения указанных параметров имеют вид волновой функции и достигают максимума в момент максимального раскрытия щелевого сопла (1 с - результаты, полученные при решении задачи в рамках ТУ; 1,063 с - СП). Смещение по времени на 0,63 с обусловлено влиянием динамики течения диспергируемой жидкости на элементы конструкции устройства. При использовании единой вычислительной модели в рамках СП рассчитанные значения напряжений и перемещений на 10 и 19 % меньше, чем при решении задачи только в рамках ТУ соответственно. 300 0,5 1 Время, с 1,5 2 я | 200 ts а я й я ^ § я 100 я я я и m 0 С 0,04 У ~ N / б / \ / / \\ ' / V / / / / / / ' / ' / \ / / \ / / \ / / Л д / / \\ i / \\ / / Л Л/ 1 1 ■ 1 1 1 1 -1-■-1-■-1- 0,5 1 Время, с 1,5 2 0,16 S 0,12 я | 0,08 0 Рис. 7. Зависимость максимальных эквивалентных напряжений (а) и перемещений (б) от времени Ввиду того, что пульсация давления диспергируемой жидкости происходит периодическим образом, за цикл работы устройства принимается промежуток, соответствующий открытию и закрытию щелевого сопла. Результаты исследования циклической долговечности устройства для подходов в рамках ТУ и СП представлены на рис. 8. 0 2,5107 т-•-1 5-107 7,5-107 Циклы Рис. 8. Зависимость циклической долговечности устройства от реализующейся максимальной амплитуды интенсивности напряжений в элементах конструкции Зона минимального запаса циклической прочности устройства совпадает с местом концентрации максимальных эквивалентных напряжений, которое находится в области крепления упругих элементов 2 к распределяющему стержню 3. Минимальное значение циклической долговечности, полученное в результате решения задачи в рамках ТУ, составляет 1,3 млн циклов, что на 35 % меньше по сравнению с результатами при СП (2 млн циклов). Такое снижение циклической долговечности обусловлено тем, что максимальные эквивалентные напряжения, полученные в рамках ТУ, превышают аналогичные параметры, рассчитанные при использовании СП. Заключение На примере задачи о деформации и ресурсе конструкции реакционной камеры вихревой горелки, использующей мелкодисперсную смесь вязкого жидкого топлива и газообразного окислителя, показано преимущество применения методики решения сопряженных задач МДТТ и МЖГ, связанных с взаимодействием жидкости, газа и элементов конструкции. Предложена физико-математическая модель, позволяющая определять реализующееся напряженно-деформированное состояние и течение жидких сред в рамках решения динамической задачи гидроупругости. Решение сопряженной задачи с использованием единой вычислительной модели отличается от решения аналогичной задачи в упрощенной упругой постановке при задании нагрузки в виде давления, равномерно распределенного на внутренних поверхностях элементов по эквивалентным напряжениям на 10 %, по перемещениям на 19 %, по циклической долговечности на 35 %. Установлено, что эквивалентные напряжения и перемещения элементов пластин сопла описываются волновыми функциями и достигают максимальных значений в момент максимального раскрытия щелевого сопла. Проведено исследование рабочего процесса, прочности и долговечности оригинальной конструкции форсуночного устройства, использующего принципы пульсационного и аэродинамического методов распылений котельного топлива. Использование сопряженного подхода к моделированию связанных процессов в устройствах реакционных вихревых камер дает возможность снизить затраты на проектирование за счет более точного прогнозирования напряженно-деформированного состояния и ресурса разрабатываемых устройств. Предложенная методика решения сопряженных задач может быть использована при разработке новых устройств, работающих в условиях вибраций и циклических нагрузок в энергетической, нефтеперерабатывающей, химической и других смежных отраслях машиностроения.
Парфирьева Е.Н., Пантелеева Ю.В. Перспективы развития мирового нефтегазохимического комплекса // Вестник Казанского технологического университета. 2012. Т. 15. № 12. С. 177-181.
Воликов А.И. Сжигание газового и жидкого топлива в котлах малой мощности. Л.: Ленхиммаш, 1989. 159 с.
Зверева Э.Р., Ганина Л.В. Повышение технико-экономических и экологических показателей мазутных хозяйств // Энергетика Татарстана. 2007. № 2. С. 62-66.
Забродин А.Г., Алибеков С.Я., Маряшев А.В., Сальманов Р.С., Филимонов С.С. Анализ физико-механических свойств мазута и устройство для его эффективной подготовки и сжигания // Вестник Казанского технологического университета. 2013. Т. 16. № 5. С. 226-230.
Архипов В.А., Бондарчук С.С., Евсевлеев М.Я. и др. Экспериментальное исследование диспергирования жидкости эжекционными форсунками // Инженерно-физический журнал. 2013. Т. 86. № 6. С. 1229-1236.
Коровина Н.В., Кудряшова О.Б., Антонникова А.А., Ворожцов Б.И. Распыление жидкости при импульсном воздействии // Известия вузов. Физика. 2013. Т. 56. № 9/3. С. 169-172.
Ковеня В.М. Некоторые тенденции развития математического моделирования // Вычислительные технологии. 2002. Т. 7. № 2. С. 59-73.
Данилов А. М., Гарькина И. А. Математическое моделирование сложных систем: состояние, перспективы, пример реализации // Вестник гражданских инженеров. 2012. № 2. С. 333-337.
Русаков С.В., Шуваев Н.В. Численное моделирование аэроупругого взаимодействия компрессорной лопатки с дозвуковым потоком воздуха в трёхмерной постановке // Вычислительная механика сплошных сред. 2013. Т. 6. № 3. С. 300-308.
Лойцянский Л.Г. Механика жидкости и газа: учеб. для вузов. 7-е изд., испр. М.: Дрофа, 2003. 840 с.
Мэйз Дж. Теория и задачи механики сплошных сред. М.: Мир, 1974. 310 с.
Когаев В.П. Расчеты на прочность при напряжениях, переменных во времени. М.: Машиностроение, 1977. 232 с.
Трощенко В.Т., Сосновский Л.А. Сопротивление усталости металлов и сплавов. Киев: Наукова думка, 1987. 1238 с.
Феодосьев В.И. Сопротивление материалов. М.: Наука, 1967. 552 с.
Когаев В.П., Махутов Н.А., Гусенков А.П. Расчеты деталей машин и конструкций на прочность и долговечность. М.: Машиностроение, 1985. 563 с.
Прошковец Й., Вайтишек Я. Расчет долговечности элементов машин, нагружаемых переменными колебательными силами // Проблемы прочности. 1980. № 8. С. 21-28.
ANSYS Workbench Verification Manual. SAS IP Inc., 2013. 236 p.