Соболь И.М. Численные методы Монте-Карло (1973) (1186217), страница 36
Текст из файла (страница 36)
Вероятность того, что частица при столкновении в Р' рассеется, равна 2 (г')/2 (г'). Обозначим 1=!г — г'1, ы=(г — г')/1, н в качестве координат точки г выберем сферические координаты с центром в (рис. 59), так что нг=(зо(оы — элемент «физического» объема. Для того чтобы рассеянная частица могла попасть в эле.
мент пг около точки г, направление ее движения после рассеяния должно оказаться в конусе еы около направления ю. Вероятность такого события выражается через иидикатрису рассеяния Р (Р, И',(2! и равна р ( г', (2'! в)ею. Есин это условие выполнеао, н части- "! Одиогрупповая (иля односкоростная) теория часто используется для расчета диффузии быстрых нейтронов, особенно в среде, состоящей из тяжелых атомов (25, 511.
$2! МОДЕЛИРОВАНИЕ СВОБОДНОГО ПРОБЕГА 223 па рассеялась по направлению ю, то вероятность столкновения в интервале (1, 1+Я), согласво (6), равна ! *и-г( — Хэ( + )»)». а Вакопсц, направление скорости й (как независимая величина в рпс. 59. фазовом пространстве) обязана совпадать с е и вероятность этого равна б(й — м)т(й. Следовательно, ! — ) елт К (Р', Р)г(Р= т, р (г', й'!М)пмк(г)» О о(б(й — ю)т(й. Введя сюцэ ЛР=т(гт(й=(тт(1 пас(й и сократив т(Р, нолучим форлтулу ) И(г'+из! Лэ К (Р', Р) .= — ', р (г', й',м)б(й — ю) Х (г') 2 (г) » ст Р (у) 2.2. Моделирование свободного пробега нейтрона. 2.2.1. В о д н о р о д н о й с р е д е. Так как во время пробега энергия нейтрона не меняется, то полное сечение постоянно: Х(Е, х) =Е(Е)ммЕ.
Из уравнения (6) следует Е(х) =1 — е-и* н метод обратных функпий позволяет записать явную формулу для расчета 2 $= — (1/Е) 1п Т. (б), 224 модГлиговмчиа естественных пРоцессов [гл б Легко вычислить, что средняя длина свободного пробега М5=!/Е. Если в качестве единицы длины выбрать 1/Х, то $= — 1п т средних длин свободного пробега.
2.2.2. В к у с о ч н о о д н о р о д н о й с р е д е (м е т о д о б р а т н ы х ф у н к ц и й). Мы рассмотрим случай, когда интересующая нас область пространства б состоиг из конечного числа однородных областей. Именно такова геометрия большинства задач, встречающихся на практике. В то >ке время очевидно, что всякую неоднородную среду можно, с любой точностью аппроксимировать такой кусочно однородной средой. 3 Поверхности, являющиеся границами областей, как правило, не сложнее, чем поверхности второго порядка: плоские, сферические, цилиндрические, конические.
Поэтому не представляет труда найти аналитически пересечение любого луча с любой из этих поверхностей. Рассмотрим простейший (по своей логике) алгоритм для расчета свободного пробега в такой среде, разработанный автором в 1955 г. В качестве иллюстрации выберем оба 4 ~К з ласть 4т, представляющую собой р с ао трехслойный конечный цилиндр с «крышкой» в форме полусферы. Сечение 6 плоскостью у=0 изображено на рис. 60, где цифры означают номера различных областей. В этом примере шесть граничных поверхностей: 1) х'+ у' = й1ц 4) г = 0; 2) х + уу = Рзы 5) г = Н; 3) х + у' = /тз; 6) х'+ у'+ (3 — Н)' = д,.
Пусть иейтрон вылетает из точки г, по направлению единичного вектора в. Тогда уравнение г'(5)=1 — т запишется в виде $ / ($): — г! ~ (гз + ыз) дз = — ! и у. (9) 'о Ф 21 МОДЕЛИРОБАННЕ СБОБОПНОГО ПРОБЕГА 225 Чтобы решить это уравнение, находим пересечение луча к=го+<аз со всеми граничными поверхностями, и все положительные значения з, соответствующие точкам пересечения, располагаем в порядке возрастания: О=а<о>(л<»(з<з>~...(з<„,>.
(10) Вычисляем длину отрезков луча 1«> =з«> — з<, >Р каждый из которых принадлежит одной области, и находим соответствующие этим областям значения Х, которые обозначим через а«>, <=1, 2,..., и> е). Обозначим через <» интегралы т»=1'(з<»>). Они легко вычисляются по формуле (» мы <<«>1«> '%' Нетрудно доказать, что если 7,( — 1и у 1,»>, то Е< а<»>+(1/с<<») ( — 1п 1' — <»); (11) если же — 1п у~(, то й;Р:з<„> н нейтрон вылетает нз обл ас < и <г. <т!»,и оп! стгл! йп Ь йп с %л Рис. 61. В самом деле, легко видеть, что з<»>($(з<„+» тогда и только тогда, когда <»(! (е) (l,+>.
Уравнение (9) ') Программа расчета должна содержать блок, позполиююнй определить, какой области приналлежпт любая заданнан точка. Выбран середину отрезка<«>, можно найти номер области, которой принадлежит (<О. 1б и. м. ссп ь 226 моделРгпонлнлге естсствеггг1ых процессов (гл а в этом случае можно переписать в виде (рпс. 61) ~л+сз<л+гг (е з(лг) 1п откуда сразу вытекает (11).
На рис. 62 изображен пример такого луча с указанием всех хи>. Легко заметать, однако, что пересечения згзг, згаг и.згаг здесь «лиш. ниеь. Чтобы исключить «лишние» пересечения, можно, например, изменить определение граничных поверхностей, добавив к иич неравенства (одно или несколько), выделяющие действительный участок границы. В примере, изображенном на рис. 62, граничные поверхности запишутся так: 1) х' -1- и' = ~з 0 ( л ( Н; 2) х +р )гзх, О( (Н' 3) х'+„,' И~~, 0(з(Н; 6) х' 4) =О, + р ( глгтз; 6) з= Н, ха+уз( )ф +Р'+(з — Н)'= Нзз, а> Н.
Рис. 63. Рпс. 62. каждый раз пришлось бы вычислять очень много пересечений, хотя в действительности кагкдый нейтрон пересекает лишь олпу-две поверхности. Для построения более экономного алгоритма в 1461 предлагается разбить граничные поверхности на элементарные граничные поверхности, каждая из которых разделяет две области. В нашем прилгере (рцс. 60) окагкется 1О таких поверхностей (в скоб.
Определив точку пересечения луча г=гл+ыз с какой-нибудь из граничных поверхностей и получив полохсительное значение з, следует подставить координаты точки пересечения в соответствуюгцие этой поверхности неравенства. Если координаты точки пересечения хотя бы одному из этих неравенств не удовлетворяют, то данная точка пересечения исключается (рис.
66). Если количество граничных поверхностей в задаче велиио, а фактические пробеги нейтронов малы, то такой алгоритм невыгоден: 4 г) ыодалнпонлннп спокон»ого пповсгл 227 ках указаны индексы повсрхпостн, т. е, номера областей, которые зта поверхност~ разделяет): 1) х +уз=)тр 0« и, (1; 2); 2) к'+ На = йтт, 0 < з( Н, (2; 3); 3) ха+У'=и, О« и, (О; 3); (О; Ил (О; 2); (О; 3)) (1; 4) (2: 4); 6) з=О 4) г=0, кз-чна(>>т и 5) х=-О, >4з) ( хз -';- уз ( Нт, 2' от (ха ( „з < >ст 7) а = Н, кз + рз ( >тт), 3) г Н >ст хт + Нг . >>т 2' 9) г = Н, )(з < ха+ и' < >(т, (3; 4); 1О) хз -1- рз -1 (г — Н)' = Нз, г > Н, (О; 4). ' Если, например, исходнан точка га принадлежит области 8, то надо найти пересечения луча т=гз+ыз только с теми злементарными граничными поверхностями, в индексах которых фигурирует 3.
Если область 3 выпуклая, то такая точка пересечения будет единственной; в общем случае выбираем точку пересечения с наименьшим г, превосходящим г)а), которое и назовем зп>. Второй индекс пересекаемой поверхности показывает, в какую область переходит нейтрон. И т. д. В примере, изображенном на рис. 63, г(1) соответствУет точке пересечения луча с поверхностью 2, а нейтрон переходит в область 2. Затем находим пересечения луча со всеми элементарными поверхностил)и, в индексе которых имеется 2, и отбираем пересечение, которому отвечает наименьшее г, превосходящее з)». Таким образом, получим значение з)з), соответствующее первому пересечению луча с поверхностью Д и номер следующей области — >. Конечно, условия У ~ †>п у< Уа)), входзщие в (11), следует проверять постепенно, по мере нахождении ( » ) » ) » т 12)' )2)' 12) т 2.2.3.
В произвольной среде ( и е т о д посто- я нного сечения). В последние годы широкое рас- пространение получил совсем другой метод моделирова- ния пробегов в сложной среде, предложенный, по-види- мому, Е. Р. Вудкоком. Кусочпая однородность среды при этом ие предполагается. Выберем произвольную постоянную )х- эпр Е, и обоз- начим через Хв разность Хо=а — Х)0. 'о'словиыся счи- 15' 22В МОДЕЛИРОВАНИЕ ЕСТЕСТВЕННЫХ ПРОНЕССОВ !ГЛ Е тать, что при столкновении нейтрона с ядрами, кроме реакций, входящих В 2, возможна еше одна — фиктивное столкновение, при котором ни энергия, ни направление движения нейтрона не меняются.
Сечение фиктивного столкновения будем считать равным г.м. Если (ср. пример гл 2, п. !.2.2) Х=Е,+л,+гн, то вероятности соответствуюших реакций в условной задаче равны у./а, Х,/а, Х,/а, а вероятность фиктивного столкновения равна г.о/а. Пробеги в этой задаче легко вычисляются по формуле (8) (! 2)' $= — (1/сс) !и Т, а тип столкновения разыгрывается с учетом всех четырех вазможностей.
Ниже доказано, что сумма таких пробегов до первого нефиктивного столкновения подчиняется тому же закону распределения (6), что истинный случайный пробег. Те о р е м з 1 ([112)). Рассмотрим нейтрон, вылетающий из точки»=0 по направлению оси Ох и подчинлющийся законам условной задачи Обозначим через $ координату первого нефиктивного столкновенич.
Тогда функция распределения Е выразсаггсл формулой (б) Д о к з з з т е л ь с т з о. Обозначим через т случайное количество фиктивных столкновений з интервале (О, $) и рзссмотрим случай т=( (1=0, 1, 2, ...). Выберем произвольные числа хь хг, ... ..., хг+г, удовлетворяющие неравенствам 0 < хт < х, < ... < хг < хе+1. Вероятность того, что фиктивные столкновения окзжутся з окрестностях точек хь ..., хь з первое иефиктизиое столкиозеиие — в окрестности точки хг+ ! равна (ср. (15) гл.
5) р(хд,..., х;+1) йх, ... ух!+, = = П аг ( ) 1 ') йх.~ф( ');< -а(х, — х ) [ ~„е (х;ц.~)) г+!— ~,р (»,) ... ~~~~ф (х,.) [а —,,Р (»! 1 з)1й»т ... йгг ь !. Проинтегрировав зту вероятность по всем возможным хь..., х,+1 э г) МОДЕЛИРОВАНИЕ СВОБОДНОГО ПРОБЕГА 229 таким, что хг+г < х, получим вероятность х х. г+! Х Р (я С х, т — 1) — ) ахг+л ) л(х .. ° ) и (хл, ." . хг~л1) ахл— о о о — ах .
г+! ) е г+' [а — ~ч~Ф (хг+1)1 4хг+г [ ~~Ф (х;) л)х, о й ... ~ ~Ф(х,) г(х,. 'о Все внутренние интегралы легко вычисляются Если ввести обозна-, к чение )ю ~х) ) ~ч~~Ф(з) г(з, то (ззллешлв хг+л иа У) можно бУдет й записать результат в виде х ))ф (У))' Р (в < х, т = л) = )г е аэ [а — УФ (у)], бу. а Отсюда следует, что — аг( l Ш1 Р (ь <к) =~~ Р ($ < х, э =1) = [е [и — ~ч„,ь (у)~ ЙУ !=о 'о — ах+1 (х) л 1 1 — е Ф = 1 — ехр [ — ) ~яр~ (з) г(з, о что и требовалось доказать. Так вак при увеличении и количество фиктивных столкновений возрастает, то обычно стараются выбрать минимально возможное значение, т.