UBASICを試みる


1. 素因数分解


整数の素因数を求めるには通常小さい数から順に割って割り切れるかどうかを調べる.
"ubasic"(注1)には,prmdiv() という組み込み関数が準備されているのでこれを使うと非常に有用である.
これは 2^34 - 1 = 17179869183 以下の数について,最小素因数を返し,最小素因数が見つからなかった場合は 0 を返す.
たとえば,次のubasic programを試みる.

10 input "N=";A
20 print=print+"soinsuu.dat" ' print命令でdata file "soinsuu.dat"に書き込む準備
30 print "N=";A;"=";
40 B=prmdiv(A):if B=0 then 70  
45 print "*";B;
50 A=A//B:if A=1 then print "*";B:goto 10
60 goto 40
70 print "*";A;"(?)":goto 10

<実行結果の例>("soinsuu.dat"にoutputされる)
N= 111111 =* 3 * 7 * 11 * 13 * 37 * 37            
N= 77777777 =* 7 * 11 * 73 * 101 * 137 * 137
N= 102456543 =* 3 * 7 * 199 * 24517 * 24517

等と得られるが,

N= 98765432123456789 =* 449 * 219967554840661(?)

では,prmdiv()が0を返す場合である.

219967554840661が素因子であるか否かは疑わしい.


木田祐司氏による"RHO.UB"(注2)というプログラムを実行すると,

98765432123456789 = 449 * 494927 * 444444443

が得られ,219967554840661は494927と444444443の素因子に分けられる.




注1: http: //www.rkmath.rikkyo.ac.jp/~kida/ubasic.htm

注2: 上のページに,「一般的な応用プログラム」として提供されている


2. 逆行列

2整数の間の四則演算のうち,割り算に関しては閉じていない.たとえば,7を3で割る場合に我々はa=7/3のように分数で表記するが,通常,計算を進めて,a'=2.333333を得る.その後,6を乗ずると,6×a'=13.999998を得る.ところが分数表記のままだと7/3x6=14のように正確な値が得られる.逆行列を求める計算において桁落ちを起こすのも同じような事情による.中学生の時よくやったように途中の計算は通分や約分をし,分数のままにしておき,最後に小数表示にすれば桁落ちを起こすことはない.Mathematicaやubasicは分数計算が得意である.
逆行列を掃き出し法(注3)で求めてみた. ubasicにおいては分数表現 (整数/整数)は(整数//整数)のように表される.例えば,

print 1//2+3//4+5//6, 7//8*9//10, 1//2//(3//4)

を実行すると,結果は次のように分数で与えられる.

25//12         63//80         2//3

"//"を"/"の1つ にすると計算は実数で遂行される.6x6 行列の逆行列を掃出し法で求める実行例を以下に示す.逆行列が分数を含むLotkin行列(注4)となる場合を扱っている(注5).右に付け加えられた6個の列は単位行列である.


Original Matrix

-6 630 -6720 22680 -30240 13860 1 0 0 0 0 0
105 -7350 88200 -317520 441000 -207900 0 1 0 0 0 0
-560 29400 -376320 1411200 -2016000 970200 0 0 1 0 0 0
1260 -52920 705600 -2721600 3969000 -1940400 0 0 0 1 0 0
-1260 44100 -604800 2381400 -3528000 1746360 0 0 0 0 1 0
462 -13860 194040 -776160 1164240 -582120 0 0 0 0 0 1

掃き出しの各ステップの結果が以下のように,分数の表現で与えられる.

Step (1)


1 -105 1120 -3780 5040 -2310 -1//6 0 0 0 0 0
0 3675 -29400 79380 -88200 34650 35//2 1 0 0 0 0
0 -29400 250880 -705600 806400 -323400 -280//3 0 1 0 0 0
0 79380 -705600 2041200 -2381400 970200 210 0 0 1 0 0
0 -88200 806400 -2381400 2822400 -1164240 -210 0 0 0 1 0
0 34650 -323400 970200 -1164240 485100 77 0 0 0 0 1

Step (2)

1 0 280 -1512 2520 -1320 1//3 1//35 0 0 0 0
0 1 -8 108//5 -24 66//7 1//210 1//3675 0 0 0 0
0 0 15680 -70560 100800 -46200 140//3 8 1 0 0 0
0 0 -70560 326592 -476280 221760 -168 -108//5 0 1 0 0
0 0 100800 -476280 705600 -332640 210 24 0 0 1 0
0 0 -46200 221760 -332640 158400 -88 -66//7 0 0 0 1


Step (3)

1 0 0 -252 720 -495 -1//2 -4//35 -1//56 0 0 0
0 1 0 -72//5 192//7 -99//7 1//35 16//3675 1//1960 0 0 0
0 0 1 -9//2 45//7 -165//56 1//336 1//1960 1//15680 0 0 0
0 0 0 9072 -22680 13860 42 72//5 9//2 1 0 0
0 0 0 -22680 57600 -35640 -90 -192//7 -45//7 0 1 0
0 0 0 13860 -35640 22275 99//2 99//7 165//56 0 0 1


Step (4)

1 0 0 0 90 -110 2//3 2//7 3//28 1//36 0 0
0 1 0 0 -60//7 55//7 2//21 4//147 3//392 1//630 0 0
0 0 1 0 -135//28 55//14 1//42 3//392 9//3920 1//2016 0 0
0 0 0 1 -5//2 55//36 1//216 1//630 1//2016 1//9072 0 0
0 0 0 0 900 -990 15 60//7 135//28 5//2 1 0
0 0 0 0 -990 1100 -44//3 -55//7 -55//14 -55//36 0 1


Step (5)

1 0 0 0 0 -11 -5//6 -4//7 -3//8 -2//9 -1//10 0
0 1 0 0 0 -11//7 5//21 16//147 3//56 8//315 1//105 0
0 0 1 0 0 -11//8 5//48 3//56 9//320 1//72 3//560 0
0 0 0 1 0 -11//9 5//108 8//315 1//72 4//567 1//360 0
0 0 0 0 1 -11//10 1//60 1//105 3//560 1//360 1//900 0
0 0 0 0 0 11 11//6 11//7 11//8 11//9 11//10 1


最後に付け加えられた6個の列に逆行列が分数の形で得られる.

Step (6)


1 0 0 0 0 0 1 1 1 1 1 1
0 1 0 0 0 0 1//2 1//3 1//4 1//5 1//6 1//7
0 0 1 0 0 0 1//3 1//4 1//5 1//6 1//7 1//8
0 0 0 1 0 0 1//4 1//5 1//6 1//7 1//8 1//9
0 0 0 0 1 0 1//5 1//6 1//7 1//8 1//9 1//10
0 0 0 0 0 1 1//6 1//7 1//8 1//9 1//10 1//11

行列式 =-31052236723200000 (注6)

注3: 「数値計算法 宇野利雄著」第2章 参照
注4n x n のLotkin行列An は次のように与えられる.この行列は数値計算プログラムの計算推移のチェックによく使用される.


最初に与えたテスト行列は A6-1である.

注5: プログラムリストは → こちら

注6: この行列式は「掃き出し」の各段階に使われる「枢軸」(イタリック体で表されている数字)を各段階において乗じて得られる.すなわち,
  (-6)×3675×15680×9072×900×11=-31052236723200000