Диссертация (1143492), страница 13
Текст из файла (страница 13)
[261];DCT квадратуры, Koch, et al. [160]−1.20.213.17−0.320.570.52−0.160−0.91−1.36−2.723.360.97−1.071.160.550.740.05−2.70−2.563.054.31−0.57−0.142.93−6.233.063.41−1.92−0.43−4.640.560.910.480.44−2.160.61−1.443.300.330.35−1.12−0.85−0.342.37−0.57−0.28−0.42−0.71−3.592.43−0.52−2.85−0.3764гралы вычисляются аналогично тому, как это выполнялось для триангуляции 1-го типапри помощи формул (3.1)–(3.7), но теперь — элементарные кусочно-квазилинейныефункции второго типа (конечные кусочно-квазилинейные элементы второго типа насфере).
Интегрирование по сфере также сводится к интегрированию по сферическимтреугольникам, аналогично тому как это описано в параграфе 3.A.2.Приведем в качестве примера квадратуры, полученные при помощи кусочно-квази(2)линейной интерполяции 2-го типа. Эти квадратуры обозначаются PQLA , = 2, 3(Piecewise QuasiLinear Angular, триангуляция 2-го типа), где значения = 1, 2, 3 соответствуют различному расположению узловых точек.(2)1Квадратура PQLA2 соответствует оригинальной триангуляции 2-го типа из работы [38]:(1) ≈ 0.00620996; (2) ≈ 0.02649607; (3) ≈ 0.02468077;(5)(4) = 0; 1 ≈ 0.2357023, 1 ≈ 0.01894174; (6) = 0.(2)2В квадратуре PQLA2 узловые точки совпадают с узловыми точками в квадратуре Лебедева L11:(1) ≈ 0.01053022; (2) ≈ 0.02212533; (3) ≈ 0.02249134;(5)(4) = 0; 1 ≈ 0.3015113, 1 ≈ 0.02041330; (6) = 0,(2)3Квадратура PQLA2 соответствует триангуляции, в которой относительная разница между наибольшей и наименьшей площадями сферических треугольников минимальна:(1) ≈ 0.01314671; (2) ≈ 0.01995346; (3) ≈ 0.02113951;(5)(4) = 0; 1 ≈ 0.3333216, 1 ≈ 0.02115905; (6) = 0.(2)Таблица 3.4 сравнивает относительную точность квадратур PQLA при вычис∫︀лении моментов >0 d по полусфере > 0.
Таблица демонстрирует, что расположение квадратурных узлов существенно влияет на точность квадратуры.Таблица 3.4. Относительные ошибки при вычислении моментовсфере > 0.Квадратура50505050110Точные значения >0 dпо полу-Относительные ошибки (в %) для моментов (, , )Числоузлов(0, 0, 1)(2)1PQLA2(2)2PQLA2(2)2PQLA2(2)3PQLA2(2)2PQLA3∫︀−2.27−0.66−0.690−0.42(0, 0, 3)(0, 2, 1)1.62−6.16−0.58−0.76−0.56−0.80−0.2−0.64/2/4−1.821.88(0, 0, 4)(0, 2, 2)3.05−4.60−1.12−1.081.701.63−3.555.30−0.42/5(0, 0, 5)4.49−1.66−1.57(0, 4, 1)(0, 2, 3)(2, 2, 1)−4.32−4.13−15.78−1.121.47−4.53−1.041.55−4.48−0.58−1.044.780.58−5.11−1.20.57−1.362/15/3/8/12/244.81С увеличением расстояние между гранями многогранника и сферой уменьшается, и разница между плоскими и сферическими треугольниками становится прене-65брежимо мала.
В этом случае интегрирование по сфере в интеграле (3.3) может бытьаппроксимировано интегрированием по многограннику, и формула (3.3) принимает вид ≈112 (3.9)где — сумма площадей сферических треугольников, имеющих узловую точку вершиной. Чтобы понять, какую ошибку вносит аппроксимация (3.9), мы строим квад(2)2(2)2ратуры псевдо-PQLA2 и псевдо-PQLA3 , в которых веса вычисляются по форму(2)2(2)2ле (3.9), эти квадратуры обозначаются PQLA2and PQLA3 , соответственно:(2)2PQLA2 :(5)(1) ≈ 0.010628; (2) ≈ 0.022103; (3) ≈ 0.022464; 1 ≈ 0.020410;(2)2PQLA3 :(1) ≈ 0.00376575; (2) ≈ 0.00996536; (3) = 0; (4) ≈ 0.00962567;(5)(5)(5)1 ≈ 0.00802249; 2 ≈ 0.00973055; 3 ≈ 0.01002474; (6) = 0.(2)2(2)2Сравнение квадратур PQLA2 и PQLA2и таблица 3.4 показывают, что приближенная формула (3.9) аппроксимирует формулу (3.3) с высокой точностью уже при(2)2 = 2.
Поэтому можно ожидать, что разница между квадратурой PQLA3 , осно(2)2ванной на приближенных весах (3.9), и квадратурой PQLA3 , основанной на точныхвесах (3.3), будет весьма мала. Кроме того, стоит отметить, что разница между весами(2)2квадратуры PQLA3и квадратуры Лебедева L17 не превышает 1-2%.(2)2(2)2Таблица 3.3 показывает, что квадратуры PQLA2 и PQLA3весьма точны исравнимы с лучшими квадратурами, имеющими близкое число узловых точек.(2)1(1)В заключение заметим, что плохая точность квадратур PQLA и PQLA2 указывает на весьма неоднородное расположение узловых точек в оригинальных квадратурахпервого и второго типа в статье Соболева [38].3.1.5Дискретизация граничного условияИзлучение, падающее внутрь области с поверхности в точке задается граничнымусловием^ + (, ) = ext (, ) + (, )∫︁′· >0| · ′ |(, ′ ) d ′ , · < 0,(3.10)где — внешняя нормаль к поверхности в точке , ext — интенсивность внешнего излучения (то есть, либо собственное излучение непрозрачной границы, либо доля излучения, прошедшего сквозь (полу-)прозрачную границу извне), , — коэффициентызеркального и диффузного отражения, соответственно,^ = − 2( · )66— направление луча, зеркально отраженного в направлении , см.
рис. 3.4.nSDΩΩ̂Рисунок 3.4. Зеркальное отражение.После угловой дискретизации граничное условие (3.10) должно быть записано вузловых точках , · < 0 (узловые точки нумеруются последовательно так, что|| = 1, . . . , и − = − , то есть число узловых точек равно 2 ). С использованием представления интенсивности в виде (3.1) эта процедура осуществляется элементарно. Подставляя представление (3.1) в граничное условие (3.10), полагая = , иучитывая, что (, ) = (), мы получаем граничное условие () = ext (, ) + ∑︁^ ) + () (||=1 () (),||=1где1 () = −и∑︁∫︁′· >0 · < 0,(3.11) · ′ ( ′ ) d ′^ = − 2( · ).Коэффициенты всегда удовлетворяют тождеству∑︁ = 1,(3.12)||=1которое выражает закон сохранения энергии.
Заметим, что коэффициенты в гранич∑︀ном условии (3.11) играют ту же роль, что и сумма 3=1 при вычислении потокаизлучения на границе по квадратурной формуле в оригинальном методе дискретныхординат ( и , = 1, 2, 3, — координаты нормали и направляющие косинусы направления , соответственно).Граничное условие (3.11) это дискретный аналог граничного условия (3.10). Кажется, что первое содержит формальное противоречие. Действительно, суммированиев (3.11) выполняется по всем направлениям , в то время как в условии (3.10) учитываются только направления, удовлетворяющие условию · > 0 (падающее награницу излучение). В действительности (бо́льшая) часть направлений, для которых · < 0, не дает никакого вклада в правую часть условия (3.11) поскольку соот-67^ ) и () равны нулю.
Тем не менее некоторые изветствующие коэффициенты (этих коэффициентов, соответствующих таким направлениям, отличны от нуля: это направления, имеющие по крайней мере одно смежное направление ′ , удовлетворяющееусловию · ′ > 0), то есть близкие к касательным направлениям. В результате праваячасть условия (3.11) может зависеть от интенсивностей (), для которых · < 0,что с физической точки зрения бессмысленно. Чтобы преодолеть эту бессмысленность,можно исключить соответствующие узловые точки из суммирования, однако такой прием может привести к значительным ошибкам.
Лучший способ состоит в экстраполяцииинтенсивностей, в направлениях, близких к касательным. Простейшая экстраполяцияможет быть реализована следующим образом. Пусть — направление, для которого^ ) или () отличны от нуля. Тогда интенсивность · < 0 и коэффициенты ( () в правой части условия (3.11) может быть экстраполирована усреднением интенсивностей ′ () по всем направлениям ′ , смежным с и удовлетворяющимусловию · ′ > 0. Более общая линейная экстраполяция также может быть легкореализована.В заключение заметим, что для вычисления интеграла в условии (3.10) можно использовать вместо представления (3.1) приближенное представление подынтегральнойфункции∑︁′′ · (, ) ≈ · () ( ′ ),(3.13)||=1тогда вместо последней суммы в условии (3.11) будет− ∑︁(3.14)∫︁(3.15)||=1где1 () =3.2 · () (),· ′ <0 ( ′ ) d ′ .Уравнение переноса излучения и граничные условия в случае осевой симметрииЭтот раздел основан на результатах статей [228, 230].3.2.1Уравнение переносаРассмотрим перенос теплового излучения в трехмерной осесимметричной области, содержащей серую (то есть без спектральной зависимости) поглощающую, излучающую и рассеивающую среду.
В декартовых координатах перенос излучения описывается68стационарным уравнением переноса излучения [29, 32, 143, 196] · ∇ + ≡ + + + = s + b ,(3.16)где ≡ (, ) (или в другой записи (, , )) — интенсивность излучения, = ( cos , sin , )— координата точки в пространстве (здесь сразу удобно перейти к цилиндрическимкоординатам), = (sin cos , sin sin , cos )— направление излучения, и s — коэффициенты поглощения и рассеяния, соответственно, = + s— коэффициент ослабления,1 ≡ () =4∫︁(, ) d,S2— интеграл рассеяния, d = sin d d (для простоты рассеяние предполагается изотропным),2 4b ≡ b () =— интенсивность излучения черного тела, ≡ () — температура, — показательпреломления среды, — постоянная Стефана-Больцмана.Если задача для уравнения переноса осесимметрична (осью симметрии являетсяось ), интенсивность излучения зависит только от переменных , , и | − |.Поэтому достаточно искать значения интенсивности только в точках(, , = 0) ≡ ( cos , sin , , , = 0), ∈ [0, ](углы ∈ (, 2) исключены из соображений симметрии), то есть искать решение уравнения переноса (3.16) в половине { > 0} трехмерной цилиндрической области и толькодля направлений = (sin , 0, cos ),(3.17)см.