Марчук Г.И., Лебедев В.И. Численные методы в теории переноса нейтронов

Подождите немного. Документ загружается.

Qu = c

+ с&в + c

CtQ

(9.4.11)

где 2c

= 1; для схемы (9.4.9) считаем, что c

= 1/4.

При аппроксимации (9.4.8) возможны три трехточечные разност-

ные схемы: ^

a) RT:m = U

+ Й

, т

= 0; тогда D

= 0 и

Q^—m

/i /р

—^з

+ Q

(S—

—m

)

6) TD : m = li

= 0; тогда Rij = 0 и

^i+^2+

(2—

—m

—Qi—tn

ij = —

*^"

Qi—m

+ (2 —m

—m

) h Q

+(2 —т

—m

)

в) RD :m

=Лй

1э

= 0; тогда 7^ = 0 и

Qi — Q

—m

h . n _ ^2—

*fi = —

I>i# =

Qi+(2—mi—m

)/i

-\-(Jb—m

—m

(9.4.12)

(9.4.13)

(9.4.14)

Каждая из схем включает в себя две произвольные постоянные.

Их можно выбирать, руководствуясь дополнительными соображения-

ми.

Например, если в ЯГ-схеме взять т

= т

= 0, то она будет с по-

ложительными коэффициентами, а если в RT-cxeue взять т

= т

= 2/2, в TD-схеме — т

= т

= 2/2, в #£>-схеме — т

= т

= 2/2, то в формуле для ошибки локальной аппроксимации исчезнут

члены вида h -^-, ft -^- .

Аналогично коэффициенты в схеме RTD второго порядка аппрокси-

мации будут иметь вид:

и==

_Й

—(2—2m

)h)(Q

+ 2m

h)"

;

= (Q

—Qx—(2— 2/nO

+ 2m

A)"

;

Du = (Si +

—2m

+ 2m

A)"

(9.4.15)

Постоянная m

в (9.4.15) может быть выбрана из дополнительных

соображений. Например, можно считать, что т

= 2/2, тогда R

= —T

j и при Q

имеем #^ = 7^ = 0, т. е. схема превращается

в RD-схеыу.

Теперь получим В£7И-схему для сингулярных ячеек. Для простоты

предположим, что шаги Ах, Ау столь малы, что каждую ячейку пере-

а б в

Рис.

9.2. Элементарная ячейка с сингулярной характеристикой

200

секает лишь одна сингулярная характеристика. Обозначим (xl

у\

р)

)

и (х\

\ у\

п)

) координаты точек пересечения сингулярной характе-

ристики с границей ячейки n

тогда х\

р)

< х\

\ у]

р)

< у\

п)

(рис. 9.2).

Для нахождения решения в точке

разделим ячейку я

на подъ-

ячейки, зависящие от конкретной конфигурации (см. пунктирную ли-

нию на рис. 9.2), и применим следующую двушаговую процедуру.

Сначала, интегрируя уравнение (9.3.1) вдоль сингулярной харак-

теристики, получаем разностную схему порядка

О (h

)

(см. § 9.6), ис-

пользуя которую, находим

в точке (х[

\ у)

п)

) через значение

в точ-

ке (х\

у]

р)

Затем находим решение

точке

С по

#Г-схеме, включая

точки (xi

, у)

п)

) и (/, / + 1) в случае конфигураций а и б или точки

{х\

\ у)

п)

) и (i + 1, /) для конфигураций виг, Для этого случая RT-

схемы более удобны по сравнению с другими схемами. Порядок ап-

проксимации будет при этом

(К).

§ 9.5. УСТОЙЧИВЫЕ СИСТЕМЫ КОНЕЧНО-РАЗНОСТНЫХ

УРАВНЕНИЙ И НЕКОТОРЫЕ ИХ СВОЙСТВА

В § 9.2 были получены системы разностных уравнений для геомет-

рий

(г,

\л)

(г,

i|),

Y).

Эти

системы хотя

обладают достаточной степенью

локальной аппроксимации [порядка

(Дг

+ А[л

)1 на гладких реше-

ниях, но при их численном решении часто возникает неустойчивость

счета. Ю. А. Кузнецовым

А.

В. Протасовым [109] построены системы

разностных уравнений для сферической и цилиндрической геометрий

(для плоской геометрии такие системы были известны и раньше), сво-

бодные от упомянутого недостатка. Их локальная аппроксимация на

гладких решениях является величиной порядка

(Дг + Д[х). Рас-

смотрим эти системы, предполагая, что функции 2, 2

, gn / непрерыв-

ны во всей области.

Плоская (х, ^)-геометрия

В этом случае кинетическое уравнение и граничные условия имеют

вид

^|L

2cp=.-|L J гЫчф'+Л (9.5.1)

1*о=|*|*';

—1<|*<1;

о<*<#;

Ф|*

о.

ц > о

U-я.

= 0.

Введем сеточную область

= {(x

VJ) : x

= ih, i = 0, N + 1, |*j = — 1 + (/ - 1/2) A|*,-

/= 1, 2m},

где N,m — целые числа; h = H/(N + 1) и Дц = 1/m, и выпишем си-

стему уравнений, аппроксимирующую в узлах сетки задачу (9.5.1):

201

Ujj—Uj+if

■2,и,

AfxS

WN+I/ = 0; i = l,N\ j = l,m;

ij —

i-lj i у „

kyZsj

2 8iJl

il+

fij>

/= 1

(9.5.2)

Фоу= 0;

i=^l

+ l,2m.

Здесь и далее для любой непрерывной функции

*ф,

зависящей, напри-

мер,

от переменных |,

т],...,

у, запись

гр^-

есть значение функции

яр

точке (l

v)j

.... v

Если из уравнений (9.5.2) исключить значения и

в граничных уз-

лах, специальным образом пронумеровать неизвестные, то приходим

к системе линейных алгебраических уравнений

с матрицей

= f

(9.5.3)

А = Ъ

—

S,.

(9.5.4)

порядка

Пх

= 2mN, где

L =

"Л 0 "

0 Л

+2 и S = 2

(9.5.5)

Весьма удобным средством записи в явном виде матриц, которые в

дальнейшем будем рассматривать, являются тензорные произведения

и суммы. Тензорным произведением двух матриц В =

(Ь

}]

)

и С = (с^)

порядка k и / соответственно называется матрица

А=В®С=

~СцВ ... с

В~

_с

В ... с

В\

порядка tn = kl, а их тензорной суммой

—

матрица

Л=ВфС=

В 0

0 С

порядка п = k + /.

Теперь матрицы Л, 2, 2

, С

и С

могут быть записаны в следую-

щем виде:

= (l/A)ji®ft; 2 = (7

02)®/

; 2

= (/

® 2

)®/

Ci=>A(iIgA©...egfl/2; С

2в

=Дц[£*Ф...0£Г1/2,

(9.5.6)

202

где

[А

= (HagU^I, ..

м | i*

|>;

i; = diag{S

, ...

2i>;

K-

= diag{S

, ...

y \ L

-1 1]

[gni

•

gilm gillm

• • •

gilm+1

gt =

: gr = \ I

_eiml

•

• •

gimm _

LSim2m

•

• •

gimm + lj

(9.5.7)

(9.5.8)

i= I, N и I

обозначается единичная матрица порядка k ^ 1.

Сферическая (г, |л,)-геометрия

Рассмотрим уравнение

№

= №';

—

1<ц<1, 0<г<#

с граничными условиями (9.2.19). Если это уравнение умножить на

г\ то после несложных преобразований придем

двум

эквивалентным

формам его записи:

—

(9.5.9)

j г(1*о)чФ'+'*/.

(9.5.10)

J_„p

+ r

(l-n«) ^~цгф + 2^

= ^

—

Введем сеточную область D

= {(r

\i/) : r

= ih, i = 0,

+ 1,

fiy = —

+ /

|x, / = 1,2

m},

где

m — целые положительные чис-

ла, h = R/(N X 1) и А(я= 2/(2m + 1). Аппроксимируя в узлах об-

ласти D

при \i < 0 уравнение (9.5.9), а при |х > 0

—

уравнение

(9.5.10),

приходим к системе уравнений

Ып

nuij—rj+jUj+ij

+ r

(i-^)^-(i-^-i)^y-i

/= i

UN+U

= 0, t = I, N, j= 1, m;

(9.5.11)

203

+ rlf

l,N

l~m

1,2m,

где

\i=*

I/tn.

Так как

1 —

О, то в уравнении (9.5.11) величины {u

}iL\

и {u

}

jZm+\ входят

нулевыми коэффициентами. Поэтому, нумеруя

неизвестные

располагая уравнения точно так же, как это было

сделано в п. 1, приходим к системе (9.5.3) с матрицей

А = L

-S

(9.5.12)

порядка и

2mN

где

-Л

Л^

и S

= 2,

С 2

W J

(9.5.13)

Матрицы С

и С

здесь полностью совпадают с аналогичными мат-

рицами п. 1 (с точностью до замены переменной

на

и Ац на Aiji),

а элементы матрицы 2 и 2

отличаются от элементов аналогичных мат-

риц п.1 только соответствующими множителями п

/ =

1, N. Выпи-

шем матрицы Л и Л:

(1 /Л)

|х® (?Кг)+ /Са®г—jx<g>г;

го о

• ° 1

•

_00

•

т_

®г.

(9.5.14)

Здесь Д и /С — матрицы, введенные в предыдущем разделе,

diag {r

, ..., r

};

а =

diag К,..., <x

} (9.5.15)

а#

(1 — fi

)/A[x,

/ =

1т.

Цилиндрическая (г, ф, у)-геометрия

Рассмотрим кинетическое уравнение для бесконечного цилиндра

в переменных (г, ф, у):

A_-!lA.)

+Sf

=^[^'}^ЙМ^^

(

)

6 о

0<г<Я, 0<ij3<n, 0<Y<1

204

с граничным условием

(Я,

<Р,

У)

0; я/2

<р

я; 0

1, (9.5.17)

где

TY'

+ V[l-Y

][l-(YT]cos(i|)-^'). (9.5.18)

Умножая (9.5.16) на г и проводя некоторые несложные преобразо-

вания, приходим

уравнениям

Ег

= ^-Г^'|^'^ЫФ +

г/;

(9.5.19)

КГ^^/гсоБф —г

sinaJ)

+ -^M-b

дг дф 2 /

. +2г

= *^-

j<fy'|Й7^ЫФ+Г/.

(9.5.20)

Построим сеточную область

{{г

7ft)

Л, t

= 0, ЛГ+1,1р,

я—/Д<р,

/ =

0,2

m+1,

Yfc

(*—1/2)

А?,

к=

l,t).

Здесь Nytriyt — некоторые положительные целые числа: h

7?/(iV +1),

Д-ф

я/(2т

1) и

\lt.

Тогда если для значений я/2 <<p

<я

рассматривать уравнение (9.5.20),

для 0<^< я/2—(9.5.19),

то со-

ответствующие разностные уравнения имеют вид

[

Г|

1/2

Ов

^ | rl

Uijh-~r}l\ui+

^ sln^tfiifc-sln^-iMfZ-u

[

fijh,

«N+i/*=o, i=i, iv, /=i,m, ^ =777;

Л А\|?

_ .

2/71 ^

■^-«.л

/= ! s=l

+ 'i/u*»

i=UN,

/=l,m,

A=l,f,

(9.5.21)

где со

Kl— vl,

1~и А^

=_я/2т. Как

и в

предыдущем раз-

деле, величины

ioh

OJk

/V;

/ = т +

2т;

ft =

1, /) в

уравнениях (9.5.21)

не

участвуют, поскольку они входят

нулевыми

коэффициентами.

205

Располагая теперь уравнения и нумеруя неизвестные в порядке

нумерации узлов сетки, опять приходим к системе (9.5.3) с матрицей

вида

A =

-S

(9.5.22)

порядка п =

2mtN.

Основные свойства матриц построенных систем уравнений

Из пп.1—3 настоящего параграфа следует, что во всех рассмотрен-

ных случаях мы приходим к системе линейных алгебраических урав-

нений (9.5.3) с матрицей

Л = L — S (9.5.23)

порядка n

где

Изучим основные свойства указанных выше матриц и некоторых их

комбинаций. Для этого наложим естественные дополнительные ограни-

чения на функции 2, 2

и g (|i

), а также введем необходимые для

дальнейшего понятия.

Предположим, что функция g(ix

) положительна и удовлетворяет

1 1

условию нормировки: — J g (|л

) d\jJ =

в плоскопараллельной и сфе-

2 —1

1 я 1

рической геометриях и— J dty' \dy'g (\i

) =

в цилиндрической гео-

я о о

метрии, а функции 2, 2

неотрицательны и 2 — 2

^ 0 во всей об-

ласти. Кроме того, предположим, что (2 — S) е — вектор с неотрица-

тельными компонентами и Se—ненулевой вектор, где е — вектор, все

компоненты которого равны единице. Удовлетворения последнего

требования всегда можно добиться путем соответствующей аппрокси-

мации интегрального оператора.

Обозначим Е

пространство вещественных векторов с обычным ска-

лярным произведением (,). Тогда вещественную матрицу В порядка п

будем называть положительно определенной (полуопределенной) в

, если (Bv, v) >

(>0) для любого ненулевого вектора v £ Е

Очевидно, что матрица В положительно определена (полуопределена)

вЕ

тогда и только тогда, когда все собственные значения симметрич-

ной матрицы В + В

положительны (неотрицательны).

Матрица В = (Ьц) порядка п называется неразложимо диагональ-

но преобладающей, если она неразложима [341] с диагональным преоб-

ладанием, т. е, Ьц > 2 | b

|, t == 1, /г, и существует такое i

(1 <

i+i

A О

Л А

+ 2 и S = 2

206

<: t

)> что b

ioJO

> 2

j |. Далее, матрица В называется

Л1-матрицей, если она неособенна и все элементы матрицы В"

неот-

рицательны. Известно (см., например, [341]), что неразложимо диаго-

нально преобладающая матрица В, все внедиагональные элементы ко-

торой неположительны, является Л1-матрицей, Симметричная матри-

ца В с диагональным преобладанием

(ЬЦ

> 2 | b

|, i = 1, я) и

}*/

неположительными внедиагональными элементами является Л1-мат-

рицей (симметричные М-матрицы называются также матрицами Стил-

тьеса) в том и только том случае, когда она положительно опреде-

лена.

Перечислим сначала свойства, которые вытекают непосредствен-

но из вида выписанных в предыдущих пунктах матриц, сделанных вы-

ше предположений и приведенных выше фактов из теории матриц.

Л — нижняя треугольная матрица с положительными диаго-

нальными и неположительными внедиагональными элементами.

Л, 2 и 2

— диагональные матрицы с неотрицательными диаго-

нальными элементами.

и С

— симметричные матрицы с неотрицательными элемен-

тами.

Матрицы L и (L + L

) одновременно являются положитель-

но определенными Af-матрицами.

5. (2 — S) — симметричная положительно полуопределенная мат-

рица.

6. А и (А + А') — неразложимо диагонально преобладающие мат-

рицы с неположительными внедиагональными элементами, следова-

тельно, они одновременно являются положительно определенными

М-матрицами [более того, все элементы матриц Л"

и (А + А

по-

ложительны].

7. Матрица L^S является неотрицательной неразложимой прими-

тивной матрицей. Поэтому, согласно теории Перрона —Фробениуса

(см.,

например, [33, 341]), спектральная задача

Хи = L^Su (9.5.25)

имеет простое положительное собственное число Х

1у

которое больше

по модулю всех остальных собственных чисел (отсюда, в частности,

следует, что К

равно спектральному радиусу матрицы L^S) и ему со-

ответствует собственный вектор с положительными компонентами.

Более того, так как разложение (9.5.23) матрицы А является регуляр-

ным,

то

< 1.

Введем в рассмотрение матрицу (в рассматриваемых нами случаях

п всегда четно)

HI 'г]

207

Тогда так как

ГS

ARS

+ SiR

ARS

~[s

ЛЯ

Si flSj

ARS

+ S

(9.5.27)

—

симметричная матрица,

где

"=[A *°4

28)

то выполняется следующее свойство:

8. Если матрица GS положительно полуопределена,

то

для задачи

(9.5.25) справедлива теория Гильберта — Шмидта

(в

частности,

все

собственные числа матрицы L~

S вещественны,

она обладает полной

системой собственных векторов), поскольку она сводится к спектраль-

ной задаче

WSu = GSL^Su

(9.5.29)

с симметричной

положительно полуопределенной матрицей GS и сим-

метричной матрицей GSL^S, которые имеют общее ядро.

Сделаем дополнительное предположение,

что

матрица

— С

неотрицательна. Для исходной задачи

это

соответствует требованию

Но

I) ^ g

(Ро)

(преобладающее рассеяние вперед),

что

почти

всегда выполняется

на

практике.

Теперь представим матрицу

А в

виде

= £-§, (9.5.30)

где

^-4V°

J»

'[cl

с!]-

Для спектральной задачи

VLU=L~

SU (9.5.32)

справедливы следующие свойства:

Матрицы

L и (L + L

)

одновременно являются положительно

определенными Л1-матрицами.

10.

Матрица Х*

^ является неотрицательной неразложимой при-

митивной матрицей. Поэтому аналогично задаче (9.5.25) спектральная

задача (9.5.32) имеет простое положительное собственное число

равное спектральному радиусу матрицы

Z^S,

которое больше

по

модулю всех остальных собственных чисел,

и ему

соответствует соб-

ственный вектор

положительными матрицами. Кроме того,

так как

разложение (9.5.30) матрицы

Л,

аналогичное разложению (9.5.3),

208

является регулярным,

то

{3411 (ы^Ях

< 1. В

заключение заметим,

что задача (9.5.32) сводится также

спектральной задаче

\iSL~

= SL^L'

Su (9.5.33)

с симметричной положительно полуопределенной матрицей

SL^S

и симметричной матрицей

SI"

S=[s

(^^+«+^s

(«

2^+/?+j?)sJ®

(9.5.34)

§ 9.6. РАЗНОСТНЫЕ СХЕМЫ МЕТОДА ХАРАКТЕРИСТИК

ДЛЯ ОБЩЕГО СЛУЧАЯ

В методе характеристик численное интегрирование проводится

вдоль траекторий нейтронов (характеристик)

пространстве коор-

динат

х = (х

л:

принцип построения разностных уравнений

одинаков

для

всех геометрий. Для определенности будем проводить

дальнейшие рассуждения

для

уравнения, записанного

форме

(5.1.2),

переобозначив

на/. Специфика рассматриваемых задач

за-

ключается

том,

что

если функции

/, с, /,

входящие

уравнение,

кусочно-постоянные

(или

кусочно-непрерывные)

области

D (что

соответствует многозонности реактора),

то

решение ср имеет особен-

ности

первых производных

на

некоторых поверхностях

в Q х D\

более того, если граница

или поверхности разрывов функций

/, с,

содержат отрезки прямых линий,- то и само решение <р (х, Q) не будет,

вообще говоря, непрерывным

по Q. В

работах

В. С.

Владимирова

[41,

43] исследованы вопросы гладкости решений

показано, что даже

при сколь угодно гладких

/, с, / и Г

может случиться^ что решение

не

будет непрерывным

замкнутой области

fi X D.

Поскольку

разрывы

решениях распространяются вдоль характеристик,

а ло-

кальная аппроксимация конечно-разностных уравнений должна

рассматриваться

на

классе решений кинетических задач [таким клас-

сом является, например, класс 0

> построенный B.C. Владимировым

[41] (см. §5.1)1, конечно-разностные методы, неэквивалентные методу

характеристик, вообще говоря,

не

аппроксимируют локально

на

классе Do кинетические задачи. Для простоты выкладок предположим,

что коэффициенты уравнения

источник

(х,

являются кусочно-

постоянными функциями

от х. Это

предположение упростит вывод

разностных схем, а распространение метода вывода их на общий случай

очевидно. Чтобы

не

затруднять дальнейшее изложение громоздкими

выкладками, последовательно введем еще ряд упрощений. При наших

предположениях областьD разобьется

на

подобластиDj,

j =

1,2,..,

и,

в каждой

из

которых все коэффициенты уравнения (5.1.2)

источник

/ (х, ft)

постоянны по

х.

Характеристиками для дифференциального

Оператора этого уравнения будут

все

прямые линии

в R

Возьмем

любую из них, совпадающую по направлению

вектором Я, обозначим

ее L. На прямой

определим числовой параметр

характеризующий

209