Метод Монте-Карло для расчёта УКИ

Расчет теплообмена в высокотемпературных тепловых агрегатах связан, как правило, с определением угловых коэффициентов облучения. Это требует проведения четырех-, пяти- и шестикратного интегрирования для каждой пары объемных и поверхностных зон излучающей системы, что в большинстве случаев математического моделирования теплообмена в реальных агрегатах весьма затруднительно. Поэтому при решении задач радиационного теплообмена все чаще применяется метод Монте-Карло [5].

В настоящее время разработаны и успешно применяются алгоритмы и программное обеспечение для расчета методом Монте-Карло угловых коэффициентов поглощенного излучения в сложных по конфигурации заполненных оптически неоднородной средой излучающих системах, состоящих, в частности, из объемных зон, имеющих форму прямоугольного параллелепипеда или прямоугольной треугольной призмы, и поверхностных зон, представленных прямоугольниками или прямоугольными треугольниками. Ниже для простоты изложения приведено описание алгоритма расчета угловых коэффициентов поглощенного излучения в системах, состоящих из объемных зон, имеющих форму прямоугольного параллелепипеда, и прямоугольных поверхностных зон. Последние в зависимости от расположения относительно объемной зоны, гранями которой они являются, подразделяются на шесть типов. Алгоритм основан на имитации распространения излучения в поглощающей среде и предусматривает возможность укрупнения ряда объемных или поверхностных зон в зоны более сложной конфигурации с неоднородными радиационными характеристиками. Для определения угловых коэффициентов проводится серия статистических испытаний. Реализация каждого испытания сводится к следующим этапам.

1. Назначается излучающая зона.

1.1. Если зона объемная, то задается ее номер, по которому определяются координаты начала и конца объемной зоны: Х_н, Y_н, Z_н; Х_к, Y_к, Z_к (размеры зоны: ΔХ = Х_К— Х_Н, ΔY = Y_К— Y_Н; ΔZ = Z_K— Z_H), а также соответствующий этой объемной зоне коэффициент поглощения.

1.2. Если зона поверхностная, то определяется ее тип, координаты и размеры, а также номер той объемной зоны, с которой граничит рассматриваемая поверхностная.

2. Внутри излучающей зоны с помощью случайных чисел γ₁, γ₂, γ₃  (0; 1) выбирается произвольная точка, представляющая собой начало траектории луча.

2.1. Если зона объемная, то координаты точки излучения находятся по формулам

2.2. Если излучающая зона поверхностная, расчёт точек излучения производится по формулам

Значения коэффициентов для этих формул принимаются в зависимости от типа излучающей поверхностной зоны по таблице:

Таблица 3.1 Коэффициенты для метода Монте-Карло

3. С помощью случайных чисел γ₄   (— 1; 1), γ₅ (0; 1), γ₆ (0; 1) выбирается направление луча.

3.1. В случае излучения объемной зоны направляющие косинусы траектории луча при изотропном излучении находятся по следующим формулам:

3.2. Если точка излучения принадлежит поверхностной зоне, а характер излучения соответствует закону Ламберта, направляющие косинусы траектории луча для соответствующих типов поверхностных зон вычисляются в соответствии с таблицей 3.2, в которой W₁, W₂, W₃ соответственно равны:

Таблица 3.2 W₁, W₂, W₃ в зависимости от типа зоны.

4. Рассматривается распространение пучка фотонов в объемной зоне, в которой находится точка излучения, а также в тех объемных зонах, которые пересекаются лучом при его распространении к поглощающей поверхности.

4.1. Вычисляется длина пути луча в объемной зоне. Для этого с учетом знаков направляющих косинусов определяются типы тех трех граней, на одну из которых попадает луч. Затем вычисляется длина пути луча до точек пересечения с плоскостями, в которых лежат указанные грани:

где A, B, C, D — коэффициенты, определяемые в зависимости от типа грани по таблице

Таблица 3.3 Коэффициенты для длины пути луча

Действительная длина траектории порции энергии излучения R_j в рассматриваемой объемной зоне j соответствует минимальному из трех найденных значений R.

4.2. В соответствии с законом Бугера вычисляется оставшаяся и поглощенная части энергии пучка фотонов после прохождения им объемной зоны:

где Е₀ — энергия пучка фотонов в начале отрезка траектории в рассматриваемой объемной зоне; а_j — коэффициент поглощения в объемной зоне j, м^-1.

Одновременно (по мере проведения серии испытаний) производится суммирование поглощенной энергии в каждой из объемных зон.

4.3. В случае достижения лучом грани, являющейся поверхностной зоной, производится суммирование оставшейся энергии пучков фотонов, попавших на данную грань (также по мере проведения всей серии испытаний) и осуществляется переход к следующему испытанию.

4.4. Если поверхностная зона не достигнута, то определяются координаты конца траектории луча в данной объемной зоне j, являющиеся уже координатами точки излучения в той объемной зоне, в которую распространился луч:

4.5. По номеру типа грани, которую пересекла оставшаяся часть порции энергии излучения, определяется порядковый номер следующей по траектории движения объемной зоны. Для новой зоны при известных значениях координат начала распространения порции энергии и направляющих косинусов повторяются вышеизложенные вычисления, начиная с п. 4.1. После попадания остатка порции энергии на поверхностную зону, где энергия полностью поглощается, в зоне излучения определяются координаты новой случайной точки и направляющие косинусы для следующего статистического испытания.

5. После завершения для излучающей зоны i серии статистических испытаний в памяти ЭВМ накапливается информация о распределении поглощенной энергии между всеми зонами системы. Учитывая, что суммарное число порций излучаемой энергии равно числу статистических испытаний N, угловые коэффициенты определяются из соотношения

где ΔЕ^k_ij — уменьшение энергии пучка фотонов, выпущенного из зоны i в испытании k, в результате прохождения объемной зоны j. Если зона j является поверхностной, то ΔЕ^k_ij есть та энергия, которой обладал пучок фотонов при достижении зоны j; — первоначальная энергия пучка фотонов.

6. Переход от угловых коэффициентов между мелкими зонами к угловым коэффициентам между укрупненными зонами производится в два этапа. После определения угловых коэффициентов для первой мелкой излучающей зоны входящей в укрупненную зону i, вычисляются угловые коэффициенты между ней и укрупненными зонами модели, т. е. производится укрупнение по зонам поглощения. Допустим, укрупненная зона j включает в себя мелкие зоны j₁^’, j₂^’,.., j_n^’. Угловой коэффициент излучения из зоны i’ в зону j определится суммой

Искомые коэффициенты определяются в результате последующего укрупнения по зонам излучения. Если излучающая зона объемная, то

где m — число мелких зон i’, входящих в укрупненную зону i; V_i, — объем, м3, и коэффициент поглощения, м^-1, мелкой объемной зоны i’.

Если излучающая зона поверхностная, то

где F_i’ — площадь мелкой поверхностной зоны i’, м2.