И про NUMERIC DIGITS
Автор: Виктор Смирнов
Дата: 23.05.2015
Источник: vasm.livejournal.com
Прочитал на Хабре интересную статью Рекуррентное соотношение Мюллера: проблемы с округлением чисел с плавающей точкой и тут же вспомнил про оператор NUMERIC DIGITS, задающий точность выполнения арифметических операций. В документации-то явно отмечено: "Не существует ограничений на значение DIGITS (кроме объёма доступной памяти)".
---
И ведь никогда особо и не задумывался над смыслом этой фразы, а тут как по темечку тюкнуло. И вся моя физическая сущность :) возопила: "В REXX, в старом добром REXX можно выполнять вычисления с ЛЮБОЙ заданной точностью!!!" Только мозгов добавляй в случае чего... Да физики на него молиться должны!!!
А рекуррентное соотношение Мюллера - лучший способ проверки возможностей REXX!
Тест № 1.
/* muller.cmd - Jean-Michel Muller's Recurrence */
/* see https://habrahabr.ru/post/258483/ */
/* Usage: muller dimension [precision] */
parse version _v
say _v
parse arg _n _p .
if (\ownCheck(_n, _p)) then
exit;
if (_p \== '') then
numeric digits _p
say 'Dimension = '||_n
say 'Precision = '||Digits()
say ''
signal on halt
signal on syntax
_t = Time('E');
do i=0 to _n
select
when (i == 0) then
_xn = 4.0;
when (i == 1) then do
_xp = _xn;
_xn = 4.25;
end; otherwise do
_xpp = _xp;
_xp = _xn;
_xn = 108.0-(815.0-1500.0/_xpp)/_xp;
end;
end;
say Right(i, 5, ' ')||' '||ownTime(Time('R'))||' '||ownNum(_xn)
end;
exit;
halt:
say 'Halted...'
exit;
syntax:
say 'Interrupted...'
exit;
ownUsage: procedure
say ''
say 'Usage: muller dimension [precision]'
return 0;
ownCheck: procedure
_n = Arg(1);
_p = Arg(2);
if \DataType(_n, 'W') then do
say 'Error - dimension is non-whole number'
return ownUsage();
end;
if (_n < 0) | (_n > 99999) then do
say 'Error - dimension is out range (should be from 0 to 99999)'
return ownUsage();
end;
if (_p == '') then
_p = Digits();
if \DataType(_p, 'W') then do
say 'Error - precision is non-whole number'
return ownUsage();
end;
_f = Fuzz();
if (_p <= _f) then do
say 'Error - precision is out range (should be > '||_f||')'
return ownUsage();
end;
return 1;
ownTime: procedure
_t = Arg(1);
if (_t == 0) then
_t = '0.000000';
return Right(Left(_t, Length(_t)-4), 9, ' ');
ownNum: procedure
_n = Arg(1);
if (Length(_n) > 43) then
_n = Left(_n, 30)||'...'||Right(_n, 10);
return _n;
|
Тесты выполнялись под OS/2 Warp v4.50, установленной на виртуальной машине VirtualBox под OpenSUSE Linux v13.1.
20 итераций с точностью по умолчанию (9 знаков)
C:\>muller.cmd 20
OBJREXX 6.00 18 May 1999
Dimension = 20
Precision = 9
0 0.00 4.0
1 0.00 4.25
2 0.00 4.470588
3 0.00 4.644731
4 0.00 4.770412
5 0.01 4.853056
6 0.00 4.856383
7 0.01 3.824443
8 0.00 -24.340356
9 0.01 125.369755
10 0.00 101.007675
11 0.00 100.049759
12 0.00 100.002483
13 0.00 100.000124
14 0.00 100.000006
15 0.00 100.000000
16 0.00 100.000000
17 0.00 100.0
18 0.00 100.0
19 0.00 100.0
20 0.00 100.0
|
100 итераций с точность 65536 знаков
C:\>muller.cmd 100 65536
OBJREXX 6.00 18 May 1999
Dimension = 100
Precision = 65536
0 0.00 4.0
1 0.00 4.25
2 0.01 4.4705882352941176470588235294...5882352941
3 34.59 4.6447368421052631578947368421...7368421048
4 73.38 4.7705382436260623229461756373...4617563639
5 76.95 4.8557007125890736342042755344...7007123788
6 77.44 4.9108474990827932004402592637...8474947548
7 76.95 4.9455374041239167247733838031...4639023381
8 75.54 4.9669625817627005987119384872...6732428062
9 75.01 4.9800457013556311612686079942...4289882547
10 74.16 4.9879794484783922601401328939...3263612238
11 73.64 4.9927702880620680974895925483...7698473716
12 72.90 4.9956558915066340266240282615...3394321302
13 72.56 4.9973912683813441128938690216...9930742186
14 72.40 4.9984339439448169190138971781...2401474625
15 72.28 4.9990600719708938678168396062...9486544354
16 72.24 4.9994359371468391479978508864...7267981961
17 72.38 4.9996615241037675377868879857...2279237545
18 72.26 4.9997969007134179126629930746...4941997026
19 72.16 4.9998781354779312492317320300...8880024436
20 72.34 4.9999268795045999044664529810...4347400137
21 72.27 4.9999561270611577381190152822...1374158392
22 72.14 4.9999736760057124445790151493...5953103632
23 72.19 4.9999842055202727079241797881...1228303459
24 72.27 4.9999905232822276594072074700...5063049141
25 72.37 4.9999943139585595936498116895...1576223065
26 72.26 4.9999965883712560237063808793...3595927231
27 72.16 4.9999979530213569079884128681...4678182904
28 72.23 4.9999987718123113299993645145...1528744657
29 72.30 4.9999992630872057845553229449...8351857232
30 72.20 4.9999995578522583058676361926...6551923789
31 72.41 4.9999997347113315241634489886...4868002178
32 72.15 4.9999998408267904691283066991...4486641201
33 72.28 4.9999999044960712411436113485...8704788268
34 72.19 4.9999999426976416501660968977...5031419126
35 72.32 4.9999999656185845960724209285...5966073504
36 72.38 4.9999999793711506157936445604...1232930199
37 72.33 4.9999999876226903184102552956...9329327332
38 72.28 4.9999999925736141726624177373...9624765115
39 72.23 4.9999999955441684969793058578...7147474753
40 72.20 4.9999999973265010958050513865...0702375476
41 72.18 4.9999999983959006566253192645...1865476973
42 72.26 4.9999999990375403936664153942...1250623819
43 72.20 4.9999999994225242360886898172...1356679051
44 72.22 4.9999999996535145416131964994...8841051144
45 72.30 4.9999999997921087249535116388...3972876061
46 72.20 4.9999999998752652349669207294...4691499747
47 72.07 4.9999999999251591409782853862...7817492083
48 72.13 4.9999999999550954845862990932...5575579629
49 72.26 4.9999999999730572907515374861...0008980724
50 72.27 4.9999999999838343744508353825...4991641484
51 72.15 4.9999999999903006246704698702...0564482585
52 72.27 4.9999999999941803748022706327...2294009992
53 72.18 4.9999999999965082248813583155...6663132023
54 72.21 4.9999999999979049349288135262...7536204001
55 72.37 4.9999999999987429609572875890...9482591658
56 72.21 4.9999999999992457765743723637...0064627367
57 72.15 4.9999999999995474659446233500...8361959129
58 72.34 4.9999999999997284795667739854...4962400471
59 72.33 4.9999999999998370877400643824...0699283659
60 72.34 4.9999999999999022526440386262...8599098056
61 72.26 4.9999999999999413515864231746...9204422343
62 72.37 4.9999999999999648109518539043...2902916176
63 72.47 4.9999999999999788865711123424...8647087333
64 72.38 4.9999999999999873319426674054...0072394395
65 72.27 4.9999999999999923991656004432...8555107408
66 72.29 4.9999999999999954394993602659...0320107125
67 72.29 4.9999999999999972636996161595...6067790579
68 72.25 4.9999999999999983582197696957...9437269546
69 72.26 4.9999999999999990149318618174...0037104102
70 72.27 4.9999999999999994089591170904...1339667416
71 72.23 4.9999999999999996453754702542...6140165342
72 72.31 4.9999999999999997872252821525...2419642925
73 72.37 4.9999999999999998723351692915...4170596592
74 72.42 4.9999999999999999234011015749...6745151377
75 72.41 4.9999999999999999540406609449...3744711367
76 72.27 4.9999999999999999724243965669...4312954169
77 72.14 4.9999999999999999834546379401...4789610486
78 72.25 4.9999999999999999900727827641...0974086739
79 72.17 4.9999999999999999940436696584...4138295731
80 72.22 4.9999999999999999964262017950...8260590372
81 72.43 4.9999999999999999978557210770...0637363413
82 72.23 4.9999999999999999987134326462...9550980251
83 72.39 4.9999999999999999992280595877...4570752688
84 72.43 4.9999999999999999995368357526...1976807769
85 72.02 4.9999999999999999997221014515...9172410019
86 72.20 4.9999999999999999998332608709...3434742889
87 72.30 4.9999999999999999998999565225...4058714846
88 72.21 4.9999999999999999999399739135...4273356348
89 72.08 4.9999999999999999999639843481...6628737559
90 72.14 4.9999999999999999999783906088...3998065416
91 72.21 4.9999999999999999999870343653...6219062318
92 72.13 4.9999999999999999999922206191...3063095021
93 72.13 4.9999999999999999999953323715...1892394725
94 72.14 4.9999999999999999999971994229...9371948898
95 72.10 4.9999999999999999999983196537...4295034201
96 72.15 4.9999999999999999999989917922...7598393402
97 72.26 4.9999999999999999999993950753...2078349622
98 72.24 4.9999999999999999999996370452...4175332678
99 72.18 4.9999999999999999999997822271...5030611944
100 72.18 4.9999999999999999999998693362...9257595911
|
Результаты говорят сами за себя. Можно поиграться с точностью расчётов, чтобы увидеть насколько быстро накапливающиеся ошибки округления приводят к неправильному результату (ряд должен сходится к теоретически обоснованному значению 5).
Также тест прекрасно отрабатывает под Open Object REXX (проверялось в v4.2, как для Windows, так и для Linux).