Марчук Г.И., Лебедев В.И. Численные методы в теории переноса нейтронов
Подождите немного. Документ загружается.
А (Я), для которых метод (10.21.2) является оптимальным. Тема же
этой задачи связана с решением неполной проблемы моментов [210] и из-
ложена в § 10.22.
§ 10.22. КОМБИНИРОВАННЫЕ ИТЕРАЦИОННЫЕ МЕТОДЫ
Представляется заманчивым исследовать возможности следующего
итерационного алгоритма [140, 143]. Сначала решаем задачу (10.2.16)
(Sp (А) £ [m, М]) методом сопряженных градиентов (10.21.2), (10.21.40).
В процессе этих итераций вырабатывается информация о границах
спектра и распределении начальной ошибки. А затем, используя эту
информацию, переходим на линейный оптимальный с весом метод
типа (10.18.2), (10.21.2), вес которого аппроксимирует спектральную
функцию начального приближения.
Пусть для оператора А >0 мы сделаем / + 1 итераций (10.21.2),
(10.21.40). Метод сопряженных градиентов на каждом шаге оптималь-
но минимизирует ошибку в некоторой унитарной норме, зависящей от
оператора задачи и начального приближения: он порождает ортогональ-
ную,
удовлетворяющую рекуррентным формулам (10.21.6) систему
многочленов {p
k
(/)}o
+1
относительно некоторого неизвестного семей-
ства весов, характеризующих начальное распределение ошибки.
Пусть М
fe
, т
к
— соответственно наибольший и наименьший корень
Pk (0>Mfc
+1
>M
ft
, 0<m^+i<m
fe
. Выбирая / достаточно большим и
пользуясь асимптотическими формулами для нулей ортогональных мно-
гочленов
[210],
определяем границы спектра m, М. Перейдя в (10.21.21)
от переменной t к переменной х по формуле (10.21.20) и нормировав
затем многочлены р
к
(х) так, чтобы коэффициент при x
k
был равен
единице, получим многочлены p
h
(х)>
для которых известны коэффи-
циенты 6
П
, К в (10.21.21), 0 < п < / + 1.
Теперь сделаем допущение, что построенная система многочленов
{Рп (*)}о
+1
принадлежит бесконечному семейству многочленов, ор-
тогональных на [— 1, 1Гс весом вида
W {x)=VT=xVB
2l+1
(х),
(10.22.1)]
де B
2
i+i (х)—подлежащий определению неотрицательный на [—1,1
многочлен степени 21+ 1. Тогда, если в формулах (10.21.21) поло-
жить б
п
= 0, %
п
= 1/4 при п >1 + 1, получим систему многочленов
{рп (*)}?> ортогональную с весом (10.22.1) при
В
21
+
1
(х) = р?^ (х)-хр
т
(х)р
г
(х) + р! (х)/4. (10.22.2)
Этому выбору параметров при k > I + 1 соответствует линейный
метод (10.21.22), оптимальный с весом (10.18.28) при б = б
2
, Q
m
(х) =
= B
2
i+i (х). „
Многочлен
B
2l+1
(t) =
B
2l+1
(х (t)) будет сравнительно мал в ок-
рестности тех точек спектра, где величины | е,? | сравнительно велики.
Таким образом, по поведению многочлена
B
2l+1
(t) можно определить
расположение той части спектра оператора Л, которой соответствует:
а) существенная доля в начальной ошибке, разложенной по собствен-
ным функциям; б) спектральный состав ошибок округления.
320
Дальнейшие итерации можно продолжить и чебышевским с весом
методом типа (10.18.2). Численные эксперименты, проведенные
Ю.
А. Власовым
[143],
подтвердили высокую эффективность такого
рода методов.
§ 10.23. МЕТОД ПОСЛЕДОВАТЕЛЬНОГО ПРИБЛИЖЕНИЯ
ОБРАТНОГО ОПЕРАТОРА
Этот метод основан на построении последовательности операторов,
приближающихся к Л"
1
(см. [332, 287, 260]). Он является обобщени-
ем итерационной формулы для отыскания обратных чисел, получен-
ной по методу Ньютона. Метод содержит достаточно экономичную схе-
му вычисления А~
1
и предполагает возможность хранения в памяти
ЭВМ
оператора, приближающего Л"
1
.
Пусть для решения уравнения (10.13.3) имеем оператор С (Я),
такой, что С
(X)
А
(X)
~ /. Тогда если положить и
0
= С/, то можно
ожидать, что е° ~
Сг°
и что приближение и
1
= и
0
+
Сг°
= С
(21 —
—
AC)f лучше аппроксимирует и, чем и
0
. Таким образом, по началь-
ному оператору С найден новый оператор С
г
= С (2/ — АС), кото-
рый должен лучше приближать Л"
1
. Разовьем эту идею. Определим
последовательность операторов согласно формулам
В
р
= /
—
С
Р
А;
С
р
+
г
= (/ +
Вр)
Ср
у
где С
0
= С; В
0
= /
—
С
0
А. Тогда
В
р+1
= / -
С
Р+1
А
= /-(/ + В
р
) (I - Вр) =
В*р,
т. е. Вр = Я
0
Ч
Если Л
(X)
>0, то можно всегда, не уменьшая общности, считать,
что 0 < Л
(X)
< 1. Тогда 0 <
1 —
А < 1. Следовательно, можно най-
ти такие операторы С
0
, что | /
—
С
0
А |< 1, например С
0
= /,а тогда
\В
Р
\
= \ (I — С
0
А)
\
2
Р
достаточно быстро стремится к нулю при
р-> оо. Поскольку
(I
+ B
k
)(I + B
k
_
1
)...(I +
B
0
)
=
(I
+ B?)(I + B?-^
=
1
+
В
0
+
...
+ В?
+1
-1-+(1-В
0
у^А-1С
0
,
последовательность операторов (/ + B
h
)
(1
+
В
к
^
г
)....
(I + В
0
) С
0
стремится к
Л*"
1
.
Значит, для нахождения и можно предложить итера-
ционный процесс
wo=,C
0
/;
u
k
+
i=(I + Bf)u
k
, (10.23.1)
в котором оператор B\
k
образуется умножением BJ*""
1
на Bf~
x
.
В книге Хаусхолдера [236] предложено использовать этот метод при
Другой вычислительной схеме. Пусть
«•-«•-С/;
гО
= /-Ли» =
В/;
1
v
k+i
=
cr
k
; r
k
+
i
=r
k
~At^+
l
=Br
k
. |
321
Тогда
u*
=
xfl
+ v
1
+...+v
k
9
r
k
=
B
k
+
l
f.
(10.23.3)
Этот процесс состоит в последовательности определенных операций над
невязками r
k
и нахождении новых добавок к приближенному решению.
§ 10.24. ИТЕРАЦИОННЫЕ МЕТОДЫ, ОСНОВАННЫЕ
НА АППРОКСИМАЦИИ НЕСТАЦИОНАРНЫХ ЗАДАЧ
Пусть
Лг(М)=£ a
h
(X)t* (10.24.1)
есть многочлен по t с коэффициентами, зависящими от
X
£ Л, в кото-
ром а
0
(X)
= А (X). Пусть g
ft
(X)
(k = 1,..., п) — корни Р
п
(Я, t), отно-
сительно которых известно, что
6* (Х)< р < 0; Я 6 А. (10.24.2)
Рассмотрим следующее нестационарное уравнение:
УМ*)
d><l
^° =/(*)- (Ю.24.3)
Для уравнения (10.24.3) поставим задачу Коши: найти решение,
удовлетворяющее начальным условиям
dt
l
=
<Pi(b),
/ = 0 п—1. (10.24.4)
Решение этой задачи запишется в виде
и (К 0= £ c
h
(X)exp[t
h
(X)t] +
A-i(X)f
t
(10.24.5)
где коэффициенты c
h
(X)
определяются из системы уравнений
£ c
k
(X) + A-i(X)f(X) =
<p
0
(X);
п
2 M*)[MW=<p,(A)
f
/= 1,2 /i-i.
Л«1
(10.24.6)
В силу условия (10.24.2)
\и(Х, t)-A'
1
(X)f(X)\<0[exp(^t)]npn t-+oo. (10.24.7)
Видно, что решение (10.24.5) нестационарной задачи стремится при
t-+-
оо к решению стационарной задачи (10.13.3). Чем больше по аб-
солютной величине значение р, отделяющее в неравенстве (10.24.2)
l
k
(X)
от нуля, тем быстрее стремление и (X, t) к своему пределу. Если
заменить дифференциальное уравнение (10.24.3) разностным по /
уравнением, определенным на некоторой сетке t = tj, / = 0..., то
получатся рекуррентные соотношения, связывающие значения фуню-
322
ции
и (Я, t) в
некоторых точках сетки
tj.
Если решение разностной схе-
мы обладает
тем же
свойством (10.24.7),
то
полученную систему рекур-
рентных соотношений можно рассматривать
как
итерационный метод
решения уравнения (10.13.3).
Существует другой способ получения системы разностных уравне-
ний,
в
котором функции
exp [l
k
(X)t]
заменяются дробно-рациональ-
ными функциями. Разберем
его при п = 1 и Л > m > 0.
Тогда задача
(10.24.3), (10.24.4) запишется
в
виде
B(X)^-^-[A(X)u-f}> и|««о
=
Фо(Ь),
(Ю.24.8)
где
В
(X)
> b >0.
Решение задачи (10.24.8) дается формулой
и(К /)={1—ехр[—S(X)t]}A'
1
(X)f(X)
+
exp[-S(X)t]^
0
(X)
1
(10.24.9)
где
S(X)
=
В^(Х)А (X). (10.24.10)
Приблизим функцию
ехр (— г)
дробно-рациональными функциями
exp
(— z) « n
q
(z)/d
p
(z),
(10.24.11)
где
n
q
(z), d
p
(2)—многочлены степени
q и p
соответственно, имеющие
р
+ q+
1
подлежащих выбору коэффициентов, которые определяются
условием
dp
(z) exp (—z) — я
д
(z) =
О
(| z
\Р+Ч
+1
)
при z -> 0.
Дробно-
рациональные разложения Падэ были известны давно,
для
получения
их можно, например, пользоваться разложением
ехр (—z) в
цепную
дробь.
Вот
разложения
с
первыми номерами:
^\^
я
р ^^*\^
0
1
2
0
1 1-2
1—z
+
*V2
1
1
2—2
2+2
6
—
42 + 22
1
1
+
z
2—2
2+2
6+22
2
1
6
—
22
12—62
+
22
2
1+2 + 2
2
/2
6+42+22
12
+
62 + 22
Используем
их для
построения итерационных формул. Формулу
(10.24.9) перепишем
в
виде
и(Х
у
t+ At) = ехр (— SAt) и
(X,
t) + (1 — ехр (- SAt)) Л"
1
/-
(10.24.12)
Пусть
р = 0,
<7
= 1,
тогда
ехр (— SAt) « 1 — SAt.
Подставляя
это разложение
в
(10.24.12), получаем, считая
u
k
= и(К t
k
),
B(X){u
k
+
x
—u
k
)IAt=
—[Au
k
—f\.
(10.24.13)
323
Мы построили разностную схему для интегрирования уравнения
(10.24.8). Для нее известно ограничение на шаг: 0 < At < 2/max |S|,
к
при котором |1 — SAt\ < 1. Считая At
k
= h —
t
h
„
x
параметром, мы
получим итерационный метод типа (10.18.1).
Пусть р = 1, (7 = 0, тогда ехр
(—SAt)
« (1 + SAt)~
l
и
5(u*-H — а*-
1
)/Д/=— (Au
k
+
l
—
/).
Это неявная разностная схема (или итерационный метод) без ограни-
чений на шаг At >0, в которой для оператора перехода Т = (5 +
+ ДМ)"
1
В всегда выполнены неравенства 0 < Т < 1 при Д* > 0.
Пусть р = ? = 1, тогда ехр (— SAt) « (1 + SA</2)"
1
(1 — SA//2)
и
B(u
k
+i—u
k
)/At=
X
<-A(u
k
+
u
k
+*)
+ f.
Это схема Кранка — Николсона, для нее Т = (1 + SA//2)"
1
X
X(l— SA*/2) и |Г|<1 при Д*>0.
Легко видеть, что при соответствующем выборе Д* и В таким путем
можно получить методы, описанные в § 10.14—10.15, а при At = At
k
(1 < k < N) — циклические методы, или методы расщепления опера-
тора.
Рассмотрим еще раз некоторые аналоги методов расщепления опе-
ратора. Пусть р =
gr
= 1, А
==
А
г
+ А
2
и В =1. Тогда
ехр
(—AAt)
= ехр (— A
x
At) ехр
(—Л
2
Д*)
»
« (1 + Л
2
Д^/2)^ (1 — Л
х
Д//2) (1 + Л
1
Д//2)-
1
(1 — Л
2
Д//2).
(10.24.14)
Формула (10.24.14) дает нам оператор Т в методе Писмана — Рэкфор-
да
[320].
Метод изучения решений стационарной задачи путем сведения ее к
нестационарной был предложен еще М. В. Остроградским. В работе
Дж. Шелдона [334] исследован вопрос о наиболее целесообразном вы-
боре шагов Atk
9
размеры которых оказались связанными с корнями
полиномов Чебышева. Как всегда, чрезвычайно интересны обратные
задачи. Они заключаются в построении замыкания вычислительного
алгоритма (по определению С. Л. Соболева [219]) для конкретных ите-
рационных методов и исследований свойств решений получаемых не-
прерывных аналогов. В работе П. Гарабедяна [280] использован такой
подход для исследования метода последовательной релаксации, а в ра-
боте Е. Г. Дьяконова и В. И. Лебедева [80] исследованы непрерывные
аналоги методов дробных шагов [249] и расщепления оператора [77].
В работе Г. И. Марчука и К. Е. Сарбасова [172] исследован итера-
ционный метод,
в
котором осреднение (см. §10.20) определяется реше-
нием нестационарной задачи
-Au
= f\
dt
n
dr
324
= 0,
&
= 0, 1,...,
m —
1.
(10.24.15)
Пусть
V
T
==
S
T
U^~
{dtAdt*...
f udt,
т
и
о
о
тогда
Av
T
^f-^rU(KT).
(10.24.16)
Следовательно, если и (X,
/)
равномерно ограничена по
к
для
всех 0</<оо, то ф
г
стремится при Т->■ оо
к
решению и°°
=
= lim и (Я, /) стационарной задачи.
/-►оо
§ 10.25. СЛУЧАЙНЫЕ ИТЕРАЦИОННЫЕ ПРОЦЕССЫ
Последовательные приближения
в
итерационных методах можно
трактовать как осуществление некоторого специально подобранного
случайного процесса, приводящего к вычислению приближенного ре-
шения с требуемой точностью и значительной степенью вероятности.
Естественно, что в соответствии с таким подходом меняются и матема-
тический аппарат исследования, и понятие сходимости (сходимость по
вероятности, в среднем, в смысле Бернулли, почти достоверная сходи-
мость), и критерии качества итерационных процессов. Случайные ите-
рационные процессы были исследованы
в
работах 10. В. Воробьева
[45—47].
Рассмотрим итерационный метод (10.18.2) при условии, что R (А)
£
£ [т, М],
т,
М >0. Пусть а — случайная величина,
а
<х
1у
...,
..., a
h
— ее значения, получающиеся
в
результате независимых ис-
пытаний. Тогда после N итераций многочлен
PN(/)>
выраженный фор-
мулами (10.18.4), будет случайным многочленом. Пусть p
h
=
|1 —
п
— a
k
t\, тогда
In \P
N
\ = 2
lnp
k
.
k=* l
Случайные величины lnp
h
имеют одну и ту же функцию распределе-
ния, определяемую функцией распределения случайной величины
а,
и они получаются в результате независимых испытаний. Поэтому вы-
полнены условия центральной предельной теоремы, и, следовательно,
при возрастании N закон распределения In
\PN\
приближается
к
нор-
мальному с плотностью распределения
К(х)
-
[l/V2nND) exp [-(x-Na)V(2ND)l
где а
=
М (In 11 — a
h
t\)
и
D — математическое ожидание
и
дис-
персия случайных величин In p
k
.
Пользуясь формулами, связывающими функции распределения, на-
ходим, что функция распределения \P
N
(t)\ стремится
к
функции
при
х >
0 и S (х)
=
0 при х < 0.
325
Оценим степень отклонения от нуля случайного многочлена
P
N
(t),
вычислив
Р (\PN\ >е),
—вероятность того,
что
\PN{1)\>
е.
В нашем случае при больших
N
оо
Р(\Ры\>е)ъ
'
Гехр[—
(In*
—
Nafl2ND\—
=
~[/2nND
J *
=
1—
Ф[(1пе—NayVmD], (10.25.1)
л:
где Ф
= -L.
Г ехр
(—*
2
)
df.
Т/я
J
Из (10.25.1) следует,
что для
того чтобы вероятность отклонения
PN
(t) от
нуля стремилась
к
нулю
с
возрастанием
N,
необходимо
и до-
статочно, чтобы математическое ожидание
а=М (In |1 — at\)
было
меньше нуля.
Как выбрать плотность распределения
ф (х)
случайной величины
а,
чтобы
и
к
сходились
по
вероятности
к
искомому решению? Естест-
венно потребовать, чтобы
ср
(я)
была отлична
от
нуля лишь
на
отрезке
[ЛГ\
/и"
1
].
Тогда
1/ш
1/т
а=
\
\n\tx—l\y{x)dx
9
D= f
(\n\tx—
11—
a)
2
y(x)dx.
l/M
\JM
Поскольку предполагается,
что
ничего
не
известно
о
распределе-
нии
R (А)
внутри
[т, М]
у
естественно потребовать, чтобы
а,
харак-
теризующая скорость сходимости,
не
зависела
от t. Это
требование
приводит
к
однородному сингулярному интегральному уравнению
\/т
dt
J
tx—
1
l/Af
1/m
решением которого
при
условии
f w
(х)
dx = \
является функция
Ф(х) = [nxV^l
—тх)(Мх—
1)]~*
при М^^х^т'
1
;
Ф
(х)
= 0 при
Л:
§Ё
[Af""
1
, mr
l
].
Для этой функции
ф (х)
можно показать,
что
а=
—
\n[(VM + Vm)]/(VM
—
Vm), (10.25.2)
а дисперсия оценивается неравенством
С<д
2
+ 8УтШ1п2 + 0[ИМ)
3
/
2
].
Пусть d= max D (Л. Тогда из неравенства (10.18.9)
и
формул
326
- 2008 — 2025 «СтудМед»