CC BY
134
Известия Тульского государственного университета Естественные науки. 2013. Вып. 1. С. 132-141
= ИНФОРМАТИКА =
УДК 519.632.4
О табулировании с высокой точностью нулей функций Бесселя
С. Д. Алгазин
Аннотация. Приводится алгоритм вычисления нулей функций Бесселя, который позволяет вычислить несколько первых десятков нулей с 15-30 знаками после запятой.
Ключевые слова: уравнение Бесселя, численные алгоритмы без насыщения.
Введение. С появлением современных ЭВМ задача табулирования изменилась. Примером нового подхода к задаче табулирования служит справочник Ю. Люка (Yudel L. Luke) [1]. В предисловии К.И. Бабенко к этой книге сказано: «... Раньше, во времена ручных вычислений, мы были, как правило, заинтересованы в том, чтобы иметь обширную таблицу, которая позволяла бы вычислять значение функции с помощью простейшей интерполяции. Теперь, используя при вычислениях ЭВМ, мы очень часто заинтересованы в том, чтобы иметь таблицу минимального объёма, даже если алгоритм восстановления функции достаточно сложен. Таким образом, возникает задача построения оптимальных по объёму таблиц». Далее показано, что возможно вычислить несколько десятков первых нулей функций Бесселя с 30 знаками после запятой.
1. О табулировании нулей функции Бесселя J0
Рассматривается задача на собственные значения для уравнения Бесселя нулевого индекса:
(xy1)1 + Xxy = 0, x € (0,1), (1)
У(1) = 0, (2)
|y(0)l < те. (3)
Краевая задача (1)—(3) эквивалентна интегральному уравнению
( x + 1) X f ( x + 1 ^ + 1) ^ + 1 ( ^ + 1)
4—) = 2j Ч—-^г) — Ч ~)
-1
G(x,0 = - ln[(x + с + |£ - x|)/2].
Применим для функции [(£+1)/2]у[(£+1)/2] интерполяционную формулу вида
T (т)
(Pnf )(x) = V f (Tk )lk (т) + Rn(x; f), lk (т) = 7------—Г, k = 1,2,...,n,
k= (т - Tk )Tn(xk)
(4)
Tn(x) = cos(n arccos т), Tk = cos[(2k — 1)^/2n], где Rn(x; f) — погрешность интерполяции. Проводя вычисления, получаем
n
Vj = ^Bjk Vk + rn(Tj; y), Bjk = Bk (Tj), yk = y(Tk), k,j = 1,2,...,n, k=i
(5)
где
Bk(т) = k+i j, £±i) lk(c)*,
-1
n.(T,y) = I ja ( ^, Ц2-1) Rn (; (i±±)V) >k.
-1
Отбрасывая в (5) погрешность дискретизации, получаем приближенную задачу на собственные значения y = \By. Здесь у — вектор приближенных значений искомой собственной функции у(х) в узлах сетки, Л — приближенное собственное значение. Легко видеть (см., например [2, с. 189]), что max |rn(x; y)| ^ c |Л| (1 + wn)En(y), где с — абсолютная постоянная,
|ж|<1
(ши = 0(ln(n)) — константа Лебега интерполяции, а En(y) — наилучшее приближение функции у многочленом степени не выше (п — 1) в норме С). Заметим далее, что собственные функции задачи (1)—(3) целые, а поэтому lim ПEn(y) = 0 (см. [3, с. 254]), т.е. величина погрешности дискретизации
n—
очень быстро стремится к нулю. Возмущение, вносимое в собственные значения отбрасыванием погрешности дискретизации, будет оценено ниже. Здесь обсудим результаты численных расчетов. При п = 5 получим первое собственное значение с четырьмя знаками после запятой, а третье собственное значение - с одним знаком после запятой. При n = 20 первое собственное значение вычисляется с 22 знаками после запятой, а 14-е собственное значение вычисляется с одним знаком после запятой. Последний расчет проводился на ЭВМ БЭСМ-6 с использованием двойной точности (длина мантиссы 80 бит). Время расчета 4 мин. 40 с. вместе с вычислением матрицы. В [4] опубликована программа BESSEL, по которой проводились эти расчеты, а также результаты расчета собственных функций краевой задачи (1)—(3).
Эти результаты получены автором в соавторстве с К.И. Бабенко давно, целью настоящей работы является их уточнение с использованием совре-
менных вычислительных средств. Сейчас появился транслятор с Фортрана (Intel Visual Fortran 9.1), который позволяет вести расчёты с учетверенной точностью REAL*16. С использованием этого транслятора были проведены расчёты с числом точек N = 3-23, N = 110. Эти расчёты сравнивались с таблицей VI работы [5]. Ниже приведены значения \/Л^, i = 1, 2,..., 13 на сетке из N = 23 узлов:
N = 23; EPS= 0.33E-28
2.40482555769577276862163187932315 5.52007811028631064959995531048766 8.65372791291101229166741865684604 11.7915344390140843686267083669051 14.9309177084599940892731158995010 18.0710639711203998901712950375023 21.2116367066011374633859922408520 24.3524696616047016958823129522958 27.4934601396705964187069258148688 30.6347976897899335897734472259312 33.7759497154294753471218322997414 36.9231498885250203330212676349967 39.9510270188728602686403023637232
Курсивом выделены знаки, совпавшие с таблицей VI работы [5]. Таким образом, сетки из N = 23 узлов достаточно, чтобы получить первое собственное значение со всеми знаками из указанной выше таблицы (29 знаков после запятой). Расчёты на сетке N = 110 дают 32 собственных значения (всю первую колонку таблицы VI) со всеми знаками таблицы (28-29 знаков после запятой).
Далее проводились расчёты на сетках из N = 3-23 узлов. Ниже приведена погрешность определения первого собственного значения: 3) 0.11;
4) 0.53e-02; 5) 0.37e-03; 6) 0.54e-04; 7) 0.2e-05; 8) 0.23e-06; 9) 0.64e-08; 10)
0.68e-09; 11) 0.15e-10; 12) 0.14 e-11; 13) 0.27e-13; 14) 0.23e-14; 15) 0.36e-16; 16) 0.28e-17; 17)0.39e-19; 18) 0.28e-20; 19) 0.35e-22; 20) 0.23e-23; 21) 0.26e-25 ; 22) 0.15e-26; 23) 0.33e-28. По этой таблице подбиралась аналитическая зависимость е = e(N). Получена зависимость lne=a+bN+ cN2+dN3+eN4+fN5; a=0.39307047, b=0.58615539, c=0.40323914, d=-0.599229377, e=0.12421718, f=-0.0076898089. Таким образом, скорость убывания погрешности экспоненциальная. Заметим, что метод конечных разностей и метод конечных элементов дают только степенное убывание погрешности.
2. Вычисление функций Бесселя целого индекса
Уравнение Бесселя имеет вид:
1 /к\2
—(v"(r) +—v'(r)) + ( -j v(r) = Лу(т), v(1) =0, |v(0)| < то.
Поэтому легко видеть, что матрица дискретной задачи (матрица, у которой необходимо вычислить собственные значения) имеет вид:
Ло = В-1,
где В — матрица дискретного оператора (обратного к оператору Бесселя), построенная в [4];
Лд = Л0 + к2 К, к = 1,2,... ,т;
К — диагональная матрица с числами (1/г^)2ти = (1 + еоз((2и — 1)п/2т))/2,
V = 1, 2,... ,т на диагонали.
Рассмотрим результаты численных расчётов. На сетке из 23 узлов получено:
МАСНЕР = 1.92Е-4; т = 23, к=1; уХ собственные значения:
1) 3.83170597020751231561443588605205
2) 7.01558666981561875353120866886788
3) 10.1734681350627220491724007823358
4) 13.3236919363142570433728670520907
5) 16.4706300508824154776364636843869
6) 19.6158585101973123876830625220157
7) 22.7600843701865157377758772901100
8) 25.9036722639552867175841920011564
9) 29.0468258315682780761054280922193
10) 32.1898292345730302050572370405816
На сетке из 140 и 280 узлов получено: т = 140 к = 1; у/\Ї собственные значения
1) 3.83170597020751231561443588630900
2) 7.01558666981561875353704998147703
3) 10.1734681350627220771857117767758
4) 13.3236919363142230323936841269486
5) 16.4706300508776328125524604709891
6) 19.6158585104682420211250658841377
7) 22.7600843805927718980530051521825
8) 25.9036720876183826254958554459800
9) 29.0468285349168550666478198835323
10) 32.1896799109744036266229841044603
11) 35.3323075500838651026344790225194
12) 38.4747662347716151120521975577165
13) 41.6170942128144508858635168050605
14) 44.7593189976528217327793527132122
15) 47.9014608871854471212740087225073
16) 51.0435351835715094687330346332243
17) 54.1855536410613205320999662145337
18) 57.3275254379010107450905042437506
19) 60.4694578453474915593987498083826
20) 63.6113566984812326310397624178737
21) 66.7532267340984934153052597500424
22) 69.8950718374957739697305364354993
23) 73.0368952255738348265061175690918
24) 76.1786995846414575728526146235348
25) 79.3204871754762993911844848724880
26) 82.4622599143735564539866106487815
27) 85.6040194363502309659494254933808
28) 88.7457671449263069037359164348545
29) 91.8875042516949852805536222144906
30) 95.0292318080446952680509981871741
31) 98.1709507307907819735377591608520
ш = 280 к = 1; собственные значения:
1) 3.83170597020751231561443588630041
2) 7.01558666981561875353704998148102
3) 10.1734681350627220771857117767752
4) 13.3236919363142230323936841269496
5) 16.4706300508776328125524604709880
6) 19.6158585104682420211250658841390
7) 22.7600843805927718980530051521805
8) 25.9036720876183826254958554459816
9) 29.0468285349168550666478198835296
10) 32.1896799109744036266229841044623
11) 35.3323075500838651026344790225187
12) 38.4747662347716151120521975577171
13) 41.6170942128144508858635168050603
14) 44.7593189976528217327793527132127
15) 47.9014608871854471212740087225063
16) 51.0435351835715094687330346332256
17) 54.1855536410613205320999662145325
18) 57.3275254379010107450905042437520
19) 60.4694578453474915593987498083827
20) 63.6113566984812326310397624178735
21) 66.7532267340984934153052597500420
22) 69.8950718374957739697305364354993
23) 73.0368952255738348265061175690917
24) 76.1786995846414575728526146235340
25) 79.3204871754762993911844848724880
26) 82.4622599143735564539866106487806
27) 85.6040194363502309659494254933803
28) 88.7457671449263069037359164348535
29) 91.8875042516949852805536222144899
30) 95.0292318080446952680509981871736
31) 98.1709507307907819735377591608515
Таким образом, имеется совпадение с 30 знаками после запятой. Отметим, что с таблицей VI работы [5] расхождений нет, но приведённые результаты немного точней.
3. Пример применения функций Бесселя целого индекса
В произвольной области Г € Я2 с достаточно гладкой границей дГ рассмотрим задачу (7), (8):
Ди(г) + / (г) = 0, г € Г, (7)
и\дГ = °- (8)
Здесь функция /(г) либо задана, либо /(г) = [д(г) + Ар(г)]и(г), где д(г) и р(г) — заданные функции, и в этом случае имеем задачу на собственные значения для оператора Лапласа. В дальнейшем будем считать, что /, д и р
— гладкие функции.
Пусть г = ф((), \С\ ^ 1 — конформное отображение единичного круга на область Г; тогда в плоскости ( формально получаем те же соотношения (7)-(8), где, однако, вместо и(г) и /(г) следует писать и(() = и(г(()) и
\Ф'«)\2/(г(С)).
Дискретизация нулевого уравнения Бесселя рассмотрена в [4]. Оказывается, что этого достаточно, чтобы построить дискретизацию двумерной
задачи для уравнения Лапласа [6]. В круге радиуса 1 дискретный оператор Лапласа имеет вид:
2 п
Н = N^21 Лк ® Нк, (9)
к=0
где штрих у знака суммы означает, что слагаемое при к = 0 берется с коэффициентом 1/2, Лк, к = 0,1,... ,п — матрица размера т х т (матрица к-го дискретного уравнение Бесселя); Нк, к = 0,1,... ,п — матрица размера N х N: Нк^ = ео8[к2п(г — ] )/N)], г] = 1, 2,... N, через ® обозначено кро-некерово произведение матриц. Здесь по в выбирается равномерная сетка: вг = 2пl/N, N = 2п + 1, I = 0,1,..., 2п. По г сетка может быть произвольна,
в данном случае выбираем ти = (1 + еов(^ — 1)п/2т))/2, V = 1, 2,..., т. То есть в единичном круге выбираем т окружностей, и на каждой окружности равномерную сетку по в:
Ло = В-1,
где В — матрица дискретного оператора (обратного к оператору Бесселя), построенная в [6];
Лд = Л0 + 2 К, к = 1,2,...,п;
К — диагональная матрица с числами (1/г^ )2, V = 1, 2,... ,т на диагонали.
В произвольной области матрица дискретного оператора Лапласа имеет вид: Z-1H, где Z — диагональная матрица с числами zvl= \ф'(£^)|2, (иі = = ти ехр(ів1), V = 1, 2,... ,т; I = 0,1,..., 2п [6].
Примеры численных расчётов. Для аналитически заданного конформного отображения вычисления проводились для эпитрохоиды, т.е. области получающейся из круга конформным отображением г = ^ (1 + є$Пр), \^\ ^ 1, є < 1/(пр + 1) (в программе обозначаются ЕРБ1 и №). На сетке т=50, N = 61 получено 73 первых собственных значения (см. ниже) для є = 1/6, пр = 4. Жирным шрифтом выделены знаки, совпавшие с расчетными на сетке 30 х 41:
М = 50 N = 61; ЕРБ1 = 0.1(6); NP = 4; время счета около 12 часов;
собственные значения:
1) 2.38444650949604970517317500817529 27
2) 3.73481160264336270303800592942815 28
3) 3.73481160264336270303800592946108 29
4) 4.60299170636341650906486673149750 30
5) 5.21305408447636594415167361552140 31
6) 5.40987176983340399646015402762260 32
7) 5.95284020254957074266667664726953 33
8) 5.95284020254957074266667664735062 34
9) 6.85046245675845252009714841665951 35
10) 6.85046245675845252009714841671300 36
11) 7.13745065059079443266158516815009 37
12) 7.25964235762084739781928108470845 38
13) 7.43089698228522450380843999247177 39
14) 8.20641593737488872010002074487221 40
15) 8.46703033293640948409921693554772 41
16) 8.46703033293640948409921696650293 42
17) 8.64626933595839253986265920247716 43
18) 8.78022273361123217092732024702811 44
19) 8.78022273361123217092732073102321 45
20) 9.43239467155670143601071573004845 46
21) 9.88298891926391094289884123527280 47
22) 9.88298891926391094290885516854269 48
23) 9.96025412843765989976019019427798 49
24) 9.96067883497458249296432391828820 50
25) 10.0677169120140574151369945503258 51
26) 10.3810329681107550983274521904108 52
10.9363318197115286486061588596269
10.9363318197115286490566933688893
11.1612248721630717492558347679048
11.2703687840812923758623562592516
11.2703687840812924379441899497323
11.7170723532432755105961602064307
11.8116162760450746824244865820387
11.8116162760450750359238289782441
11.8373790709381441799296422710614
12.2479455414825741525876370390611
12.4578735789763446630870262555313
12.4840161602227866181362114410024
12.7174707310530584401839103924169
12.7174707310530591099976162945176
13.1733648658278880535709423991943
13.2265396248205802726012953547913
13.2576192483420038095986918397910
13.3763092234868143064268852369909
13.3763092234873412793268306810793
13.7284082241398722372390888657113
13.7284082241405295639937794836824
14.0845624387645828041014741776930
14.2896890432957022566899381892198
14.3468569599344682193714374693753
14.3660088868535119687103057007520
14.3660088870628754963926158794226
53 14.5260032444333450006408510849440 64 16.0388925568050805309490904588994
54 14.7823390227650328822139182427564 65 16.1932121816405151627558326565322
55 14.7823390227701605589790191475342 66 16.3265204999443634036312371287739
56 15.0607377174290784301671781394397 67 16.3265205609221678113712393690573
57 15.1087606415694495735270739799535 68 16.3429023570436304015647080529037
58 15.3618732658649948819304045893759 69 16.6567315100655668236466221185369
59 15.4412325849722392654167903907534 70 16.8040692371885309000559415659651
60 15.4412325852883705119617736204168 71 16.8338073184134952489995022549795
61 15.6236132210873097648998985968366 72 16.8914598929759515336473804892064
62 15.7224016836042896639936517952767 73 16.8914602236018460187330833459450
63 15.7224017024130703061339497646196
Обсуждение полученных результатов. Расчёты производились на ПЭВМ Pentium IV с тактовой частотой 3,00 ГГц и объемом оперативной памяти 1 Ггб. Время последнего расчёта на сетке 50 x б1 около 12 часов. Предыдущий расчёт на сетке 30 x 41 занимает около 30 минут. Как видно из сравнения результатов на двух последних сетках, надежно с б-7 знаками после запятой определяется 30 собственных значений. Первое собственное значение определяется с 21 знаком после запятой. Расчёты проводились с учетверенной точностью REAL*^ с использованием транслятора с Фортрана Intel Visual Fortran 9.1. В приведенной выше таблице оставлены все знаки, которые даёт учетверенная точность. Знакам, выделенным жирным шрифтом, можно определенно верить. Первые б собственных значений получены на сетке 30 x 41 на ЭВМ БЭСМ-б [б, таблица 3.2, стр. 48] . Все знаки, приведённые в этой таблице, совпали (кроме последней цифры, по которой проводилось округление).
4. Нули функций Бесселя полуцелого индекса Jn+1/2
Вначале приведём результаты вычислений в Maple 11 для n = 0,..., 10: eval^BesselJZeros^^,! ..10));
3.141Б92бБ4, б.2В31ВБ307, 9.4247779б1, 12.Ббб370б1, 1Б.7079б327, 1В.В49БББ92, 21.99114ВБВ, 2Б.13274123, 2В.274333ВВ, 31.41Б92бБ4
evalf (BesselJZeros^^, 1 ..10));
4.4934094БВ, 7.72Б2Б1В37, 10.904121бб, 14.0бб19391, 17.2207ББ27, 20.3713029б, 23.Б194Б2Б0, 2б.ббб0Б42б, 29.В11Б9В79, 32.9Бб3В904
evalf (BesselJZeros^^, 1 ..10));
Б.7б34Б9197, 9.09Б011330, 12.32294097, 1Б.Б14б0301, 1В.бВ903б3б, 21.ВБ3В7422, 2Б.012В0320, 2В.1б7В2971, 31.32014171, 34.47048833
evalf (BesselJZeros(3.Б, 1 ..10));
б.9В7932000, 10.41711ВББ, 13.б9В0231Б, 1б.923б2129, 20.121В0б17, 23.30424б99, 2б.47б7б3бб, 29.б42б04Б4, 32.В0373239, 3Б.9б140БВ0
evalf (BesselJZeros(4.Б, 1 ..10));
В.1В2Бб14Б3, 11.7049071Б, 1Б.039бб471, 1В.3012ББ9б, 21.Б2Б41773, 24.727БбБББ,
27.91ББ7б20, 31.09393321, 34.2бБ39009, 37.43173б77
eval^BesselJZeros^^, 1 ..10)); 9.3ББВ12111, 12.9ббБ3017, 1б.3Б4709б4
29.332Бб2БВ, 32.Б24бб129
evalf (BesselJZeros^^, 1 10.Б12В3Б41, 14.2073924б 30.7303В073, 33.93710В30
evalf (BesselJZeros^^, 1
11.бБ703219, 1Б.4312В921 32.11119б24, 3Б.3331941В
evalf (BesselJZeros^^,! . 12.79078171, 1б.б41002ВВ 33.47бВ00В2, 3б.714Б2913
evalf (BesselJZeros(9.Б, 1 13.91Б822б1, 17.838б4320 34.828б9бБ4, 38.08247909
evalf (BesselJZeros(1С.Б, 1 1Б.0334б930, 19.02Б8Б3Б4 3б.1б81Б714, 39.43821448
3Б.707Б7б9Б
.10));
17.б4797487
37.13233172
.10));
18.92299920
38.Б413б48Б
10));
20.1824707б
39.93б12781
.10));
21.42848б97
41.3178б4б9
..10));
22.бб2720бб
42.б87бБ128
Б34б3411, 32
Ниже приводятся расчёты, аналогичные пункту 3: MACHEP =
19.бБ31Б210,
38.883б309б
22.904ББ0бБ, 2б.1277Б014,
20.9834б307, 24.2б27б804, 27.Б078б83б, 40.318892Б1
22.29Б34802, 2Б.б028ББ9Б, 28.8703733Б, 41.7390Б287
23.Б9127482, 2б.92704078, 30.2172б271, 43.14Б42Б02
24.87321392, 28.2371343б, 31.ББ018838, 44.Б39144б3
2б.1427б7б4, 29. 4Б.9212017б
870Б34б0,
1.9E-34.
m = 320, п=0.Б, л/ХЇ, i = 1,..., 10:
1) 3.14159265358968148253314596484760
2) 6.28318530717980998847751005721168
3) 9.42477796076904444882659617826803
4) 12.5663706143596199745004153513290
5) 15.7079632679484074188026742442073
6) 18.8495559215394299556118042407448
7) 21.9911485751277703949203122405416
8) 25.1327412287192399293501443109994
9) 28.2743338823071333796442247984913
10) 31.4159265358990498932469544081448
m = б40, п=0.Б, л/ЛЇ, i = 1,..., 10:
1) 3.14159265358978974602450302670173
2) 6.28318530717959346180037664523438
3) 9.42477796076936923807827248853378
4) 12.5663706143591869235912261908801
5) 15.7079632679489487301463333079038
6) 18.8495559215387803853630192704717
7) 21.9911485751285282222382168356239
8) 25.1327412287183738471062219818291
9) 28.2743338823081077143634600934317
10) 31.4159265358979673088112936148603
Таким образом, получается, что квадратные корни из собственных значений кратны п. Для сравнения приводим значение п с 4Б знаками, посчитанное в Maple 11: 3.141Б92бБ3Б897932384б2б4338 32 79Б028841971б940. Следовательно, получить больше 1Б верных знаков после запятой не получается. Задача с полуцелым индексом более сложная для вычислений. Дальнейшие вычисления приведены ниже. Отметим, что для n = 1, 5 проводились расчёты с m = 2240, в первом нуле получено 23 знака после запятой, а в 31 — 20 знаков после запятой.
5. Пример применения функций Бесселя полуцелого индекса
Рассмотрим задачу о собственных колебаниях сферы радиуса г0 с нулевыми граничными условиями первого рода. Эта задача сводится к отысканию собственных значений и собственных функций уравнения Аь + Хь = = 0 с граничным условием на поверхности сферы V = 0, [7]. Обозначим
(га) (га) (га) т ґ \ г\
VI , V2 ,...^т корни трансцендентного уравнения ^га+1/2^) = 0, нахо-
дим собственные значения Хт,п = . Каждому собственному значе-
нию Хп,т соответствует 2п+1 собственная функция. Введём обозначение
Вычисления проводились с помощью подпрограмм: QBESSEL вариант с учетверенной точностью программы BESSEL [4,6], QTELM фортранный вариант с учетверенной точностью алгольной программ elmhes [8], N — размер матрицы, H — исходная матрица заменяется на выходе матрицей в верхней форме Хессенберга; QHQR фортранный вариант с учетверенной точностью алгольной программ hqr [8], N — размер матрицы, H — исходная матрица в верхней форме Хессенберга (не сохраняется), MACHEP — машинная е, равная наименьшему положительному числу для которого 1 + е > 1, WR, WI выходные массивы длины N, которые содержат на выходе действительную и мнимую часть вычисленных собственных значений, CNT — целый массив длины N, который содержит на выходе число итераций для вычисления каждого собственного значения. QMINV — вариант с учетверенной точностью подпрограммы MINV для обращения матрицы [9].
Заключение. По поводу получения полных версий программ обращайтесь по электронному адресу: [email protected] или на адрес Института проблем механики РАН, 119526, Москва, проспект Вернадского д.101, к.1.
1. Люк Ю. Специальные математические функции и их аппроксимации. М.: Мир, 1980. 608 с.
2. Бабенко К.И. Основы численного анализа. М.: Наука, 1986. 744 с.
фп(х) = \[П^п+1 /2(х). Тогда собственные функции, рассматриваемой задачи, можно представить в виде:
(п = 0,1,m = 1, 2,j = -n,..., -1, 0,1,..., n), уП0)(в,ф) = Pn(cos в), уП-^(е,ф) = pnj)(cos в) cos jф,
Уп](в,Ф) = Pij)(cos в) sin jф, (j = 1,2, ...,n).
Pnj — присоединенная функция Лежандра j-го порядка [7]. Примечание.
Список литературы
3. Гончаров В.Л. Теория интерполирования и приближения функций. М.: Гостех-теориздат, 1954. 328 с.
4. Алгазин С.Д., Бабенко К.И., Косоруков А.Л. О численном решении задачи на собственные значения. Препринт. М.: ИПМ, 1975. №108. 57 с.
5. Таблицы нулей функции Бесселя. Библиотека математических таблиц. М.: ВЦ АН СССР, 1967. Вып.44. 95 с.
6. Алгазин С.Д. Численные алгоритмы классической математической физики. М.: Диалог-МИФИ, 2010. 240 с.
7. Тихонов А.Н., Самарский А.А. Уравнения математической физики. М.: Наука, 1977. 742 с.
8. Уилкинсон Дж, Райнш Ц. Справочник алгоритмов на языке АЛГОЛ. Линейная алгебра. М.: Машиностроение, 1976. 391 с.
9. Сборник научных программ на фортране. Матричная алгебра и линейная алгебра. М.: Статистика, 1974. Вып.2.
Алгазин Сергей Дмитриевич ([email protected]), д.ф.-м.н., ведущий научный сотрудник, Институт проблем механики им. А.Ю. Ишлинского РАН, Москва.
About tabulation with split-hair accuracy of zero of cylindrical
functions
S. D. Algazin
Abstract. The algorithm of calculation of zero of cylindrical functions which allows calculating some first tens zero from 15-30 signs after a comma is reduced.
Keywords: Bessel equation, numerical algorithms without saturation.
Algazin Sergey ([email protected]), doctor of physical and mathematical sciences, leading researcher, Institute for Problems in Mechanics of RAS, Moscow.
Поступила 17.01.2013
