自作TTLコンピュータで円周率を計算する（第1話）

なぜ、円周率なのか？

円周率の計算自体が目的ではないのです

大昔（1940-1950年代）は、コンピュータとは「特定の計算を行う機械」でした。その後もしばらくの間は「科学計算」と「事務計算」を同一のマシンで行うことは不可能もしくは非効率と信じられていました。しかしIBM360の登場(1964)を契機に「ソフトウェアの変更のみで何でもできるのがコンピュータ」という考え方が広まり、その考え方は現在では「当り前」とまでなりました。

さて、私がTTLコンピュータの製作を決意した時、最初に悩んだのは「コンピュータとは何か？」という問題でした。「特定の処理に特化した論理回路を作るのは比較的簡単だが、それをコンピュータと呼んでよいのか？」という自問でもありました。（詳細は設計編参照）

結局、「円周率の計算」と「オセロゲームの相手」という様な、全く異なる処理をプログラムの変更のみで実行可能な論理回路を作ることができれば、「それがコンピュータだ」という結論に落ち着きました。

ですから、本当は全く異なるプログラムであれば、「素数の算出」と「五目並べ」でも、「電子オルガン」と「家計簿」でも、何でも良かったのです。

オセロゲームのプログラムの方は、まだ何も着手していません。円周率の算出プログラムを優先して考えています（2011-7-30)

甘く見ていた円周率

実は円周率の計算アルゴリズムはその辺にいくらでも転がっているだろうと楽観視していました。
しかし実際に検索してみると「プログラム化を全く考慮していない数学の公式」と「高速CPUと大容量メモリのマシンを対象としたアルゴリズム」ばかりです。

実効メモリが4096バイト、算術演算は8bitの加減算のみ、という超貧弱自作TTLコンピュータに適用できるプログラムなどどこにもありません。せめて本機にサブルーチンコール命令があれば乗除算サブルーチンを作り何とかするのですが、当然、そんな高級命令があろうはずもありません。

「やはり自作TTLコンピュータではまともな処理は無理」と白旗をあげることも考えましたが、Twitterで『円周率、頑張ってください』と先に言われてしまい、後に引けなくなりました。
それに、1940年代のエンジニアはもっと貧弱なハードウエアで様々な数学計算を実現しています。ここは一念奮起して、ロートル技術屋の意地を見せることとしました。

まずは、C言語で組んでみる

結局、いきなりRETROF-8のアセンブラで円周率算出プログラムを組むのは逆立ちしても無理と悟りました。

ここからが、半世紀ぶりの猛勉強。まずはエクセルで「級数」なるものを理解し、とりあえず、自作TTLコンピュータは置いといて、32bitのWindowsマシンとＣ言語（雑誌のオマケのC++コンパイラ）で、円周率を求めるプログラムを幾つか書いてみました。

その中で、ソースコードを工夫しながら変形することを繰り返すことによって、なんとかRETROF-8の貧弱命令セットに変換できそうなのが、下記のコードです。
　VC++Express2008用のソースですが、コピペ後にヘッダーファイルを適切なものに変更することにより、他の処理系でも動くと思います。

円周率を求めるプログラム　Ｃ＋＋版

ソースプログラム

#include "stdafx.h" // ヘッダは処理系に依存する
const int ARRAYMAX = 255 ;

int main(void) {

    unsigned short nom[ARRAYMAX+1];  // 分子
    unsigned short dig[ARRAYMAX+1];  // 小数各桁の値
    unsigned short n ;
    
    for(unsigned short i=0; i <= ARRAYMAX; ++i) nom[i] = 2; 
  
  
    for(unsigned short i = ARRAYMAX ; i > 0 ; --i) {
        n=0;
            for(unsigned short j = i-1; j > 0; --j) {
            n      = n * j + nom[j] * 10;
            nom[j] = n % (2*j-1);
            n      = n / (2*j-1);
            }
        dig[ARRAYMAX-i] = n;
    }
  
    //以下は結果の表示用のコード。
    //TTLコンピュータでは結果を配列に格納するので以下のコードは不要
        for(unsigned char  i=ARRAYMAX; i>0; --i) 
            if (dig[i]>9) dig[i]-=10,dig[i-1]++ ;
        for(unsigned char  i=0; i < 80; ++i) {
            if (!(i%10)) printf(" ") ;
            if (!(i%50)) printf("\n") ;
            printf("%d",dig[i]) ; 
        }
    getchar() ; //一時停止用。
}   

実行結果

やはり精度がイマイチで75桁目あたりから誤差がでていますが、50桁算出できればヨシとして、次のステップ（乗除算の展開）に進むこととしました。

前言撤回、やはり100桁を

個人の趣味としてのモノ作りであっても、目標策定後に目標を下げてはいけない

前項で、50桁算出できればヨシと書きましたが、「『百桁求める』と言ったやないけ！」とのお叱りのメイルを頂戴しました。私にも「意地」はありますので、100桁算出できるようにプログラムを書き直しました。（どうだ、まいったか！　^^）

個人の趣味としてのモノ作りは、当初の目標を下げるのが実に簡単です。でも、それに甘えてしまうと陳腐な作品しか残らず、腕も上がらないのも事実です。メイルを頂いたＡ氏に感謝いたします。目が醒めました。

このプログラムは配列を倍の大きさにしたので、150行までは正確に算出できます。ARRAYMAXの値を大きくすれば更に多くの桁を求められるのですが、RETOROF8の主メモリは4KBしかありませんので、short int（16bit)の配列宣言は要素数1000程度が限界です。以下のプログラムでは、実行時間や他の変数領域、プログラム領域（命令領域）も考慮して、配列のサイズは512としています。

#include "stdafx.h" // ヘッダは処理系に依存する
const int ARRAYMAX = 512 ;
unsigned char dig[ARRAYMAX];  // 小数各桁の値
void debug_print(int from,int to) {
    for(int p=ARRAYMAX-1; p>0; --p)  //(注）
        if (dig[p]>9) dig[p]-=10,dig[p-1]++ ;
    for(int p=from; p < to; ++p) {
        if (!(p%10)) printf(" ") ;
        if (!(p%50)) printf("\n") ;
        printf("%d",dig[p]) ; 
    }
    getchar(); //一時停止(VC++コンソールアプリ)
}

void main(void) {
    //*******************************************
    unsigned short i,j,n,nom[ARRAYMAX];  
    for(i=0; i < ARRAYMAX ; ++i) nom[i] = 2; 
    for( i = ARRAYMAX-1 ; i > 0 ; --i) {
        n=0;
        for( j = i-1; j > 0; --j) {
            n      = n * j + nom[j] * 10;
            nom[j] = n % (2*j-1);
            n      = n / (2*j-1);
        }
        dig[ARRAYMAX-i-1] = n;
    }
    //*******************************************
    debug_print(0,ARRAYMAX) ;
}   

関数分割

RETOROF8での実行は、//****で挟まれた部分のみをRETOROF8の機械語に変換します。
debug_print（）はC++で、結果配列の内容を10進値で表示するためルーチンです。RETOROF8単体で実行する際には、dig[ARRAYMAX+1]に対応するメモリ領域に各桁が格納されます。ロータリーSWでアドレスを指定し、LEDの点灯で各桁の値を確認することが可能ですので、debug_print()に対応する機械語プログラムは不要です。

プログラム中の //（注）と書かれた行のforループは小数点以下の10進表現の正規化処理です。
このループ文を削除しても正しく計算はされますが、dig[]に格納される小数各桁が10を超える場合があります。例えば3.14は、本来ならdig[]に「３」「１」「４」と格納されるのですが、この正規化処理を省くと、「３」「１」「４」ではなく「３」「0」「1４」と格納されます。一瞬、戸惑いますが計算結果自体は正しいので、RETOROF8で実行する際にはこの処理は割愛することとします。

実行結果

152桁目までは正しい様なので、当初公言した（してしまった）100桁はクリアです。
でも、short int や乗除算を含むプログラムを、加減算しかできない8bitコンピュータの機械語に変換（コンパイル）することができるのでしょうか？？？　
何となく、大法螺を吹いてしまった気もします。（挫折しなければ、第２話に続く...）