16ビットパーソナルコンピュータの演算処理機能藤村公

男

1■9々り04戸0ρ0

はじめに

高級言語と数の内部表現
浮動小数点演算の精度
四則演算と数値関数の計算速度

Cray社スーパーコンピュータとの比較

結論

1 はじめに

近年におけるパーソナルコンピュータ（以下パソコンと略す）の発展は質量ともに著しい。いうまでもなく，これは高性能な高密度集積回路の多量生産が可能となったためである。低価格のパソコンが市場に提供され，需要を刺激し，これがまた高性能製品の開発を促すといった需要供給の好循環によって，極めて単期間にパソコン市場が拡大された。現在でもパソコン市場は
年率20％の急成長が続いている。パソコン本体だけでなく，ディスクやディ
スプレイ，プリンタなどの周辺装置も高機能，低価格へと向かっている。利
用者にとってこれほど有難いことはなく，大学の研究室にも，場合によっては自宅にもパソコンを導入し，在宅研究が可能となった。

ハードの技術革新はソフトの多様化をもたらした。パソコンはワードプロ
セッサや表計算，データベース，データ通信など，主として事務用ソフトの開発によって広く普及し元といえる。大ヒットしたNECのPC98シリーズ
33

第31号

徳山大学論叢では6，000本以上の応用ソフトが利用可能といわれている。

パソコンの利用は事務用ばかりではない。Per Bak［1］は専門的な科学
技術計算の分野で，大型汎用コンピュータと比較しながら，パソコンの利用は有効で，もっと多用すべきであると主張している。彼はis like a dinosaur：alarge body with a small brain

A large computer（「大型コンピュータ

は，巨大な体躯に小さな頭脳しか持たない恐竜のようなものである」）と評している。大型コンピュータシステムを維持してゆくには，コンピュータ本体の価格より専用の建物や空調設備，ワークステーションなどの周辺装置，システムを管理する人件費等のオーバーヘッドが巨額となり，コンピュータの使用コストを著しく高めている。彼はNSF（National Science Foundation：

米国科学研究財団）の検閲者であるという立場にある。NSFの研究費の中で占める大型コンピュータ使用料金は無視し得ない程度となり，これが米国政府の財政を圧迫する一因となっている。またNSFの研究費を得て行われている研究の大部分は，オーバーヘッドの少ないパソコンで代用できるとしている。これに対して多くの科学者が反論している［2］。Per Bakの論評は少し偏っている面があるが，科学研究の分野でパソコンの果す役割は増大
しっっあることは確かであろう。

Per Bakの論文が動機となって筆者は現在最も普及していると考えられる

16ビットパソコン（FM 16β：80186MPU，8MHz， RAMとして1MB実装）
を用いて，その演算処理機能（演算精度，計算速度，メモリ利用率）について，系統的に検討することにした。科学計算はFORTRANなどのプログラ
ミング言語によって行われる。パソコンの普及と共にBASIC， PASCAL， Cなどの高級言語も開発された。本稿では，これらの言語の演算機能を比較検討する。また，MPUと連結して， MPUだけでは出来ない浮動小数点演算
をハードウェア的に高速処理する数値演算コプロセッサが開発されている。
この効果についても検討したい。

一34一

1989年6月 藤村公男：16ビットパーソナルコンピュータの演算処理機能

2 高級言語と数の内部表現使用する高級言語は富士通のF‑BASIC86インタプリタV3．1［3］，同コ
ンパイラV3．10［4］， Prospero Software社のPro FORTRAN Version
iid 2．1［5］と同社のPro PASCAL Version iid 2．2［6］， Lattice社のLattice C V 3．10［7］の5種である。オペレーティングシステム（OS）はマイクロソフト社のMS‑DOS V 3．1［8］である。 F‑BASIC86インタプリタとコンパイラはソースプログラムの処理方式が異なるだけで文法は同じである。インタプリタはプログラムの一命令ごとに機械語に翻訳，実行する。コン

パイラはソースプログラムの全体を機械語に翻訳して実行可能形式のプログ
ラムを作成し，実行する。後者はコンパイルという面倒な操作を必要とするが，インタプリタと比べて，プログラムの実行が高速となる。FORTRAN，
PASCAL， Cはコンパイラ言語であり，ソースプログラムの編集，コンパイ
ル，リンクという複雑な手続きを経て，実行可能形式のプログラムが作成される。しか．し，プログラムが構造化されていること，アセンブラを含めた他のプログラムとの結合が容易なこと，プログラムの実行速度が高いことなど，

インタプリタにない優れた特徴を持っている。

コンピュータの内部ではすべての情報は0と1の2進法で表現され，処理される。従って数値計算で扱う数をいかに2進法で表現するかが問題となる。16ビットMPU（マイクロプロセッサユニット）は基本的に16ビット（2バイト）単位1）でデータ処理を行う。16ビットMPUである80862）では，

16ビットの汎用レジスタが使われ，RAMメモリも16ビット単位でアクセスされる。これらのレジスターレジスタ聞またはレジスターメモリ間で加減乗除の四則演算が行われる。従って数値表現も16ビットか16ビットの整数倍の注1）情報量の単位としてビットの他にバイト（＝8ビット）を用いることもある。
2）本稿で使用したFM 16βは80186だが，これは8086と周辺チップを一つにまとめてオンチップ化したもので，機能は8086と同じである。

一35一

徳山大学論叢第31号

ものが最も効率が良い。3）

一方，高級言語で処理する数は整数型と実数型に区別される。：FORTRAN
では更に複素数型があるが，これは複素数の実部と虚部をそれぞれ実数型とする構造体データとみなされる。従って複素数型も基本的には実数型に帰着される。コンピュータによる数値計算で問題となるのは計算精度または有効
桁数4）と計算速度である。精度と速度は少しでも高い方が良い。しかし，精度を高めるため有効桁数を多くすると計算に多くの時間を要する。従って，計算の精度と速度は互いに相反する。整数演算はMPU自身が行うので，一般に高速でクロック周波数にほぼ比例して速度は高くなる。しかし，実数型の場合は表現方式も演算方式も多様で言語によって異なる。MPU自身は実数型を直接処理することは出来ないので，通常，高級言語は実数型演算処理を行うプログラムを内蔵している。一般に実数演算時間は整数型より数十倍になる。MPUでは出来ない実数演算をハードウェア的に高速処理する専用の集積回路が開発されている。これを数値演算コプロセッサ（Numeric Data

Processor：NDP）という。NDPはソフトウェアでは手間のかかる実数型
演算をMPUと連結してハードウェア的に高速処理する。パソコン本体のMPU基板に既にソケットが配線済みとなっており， NDPをそのソケットに装着するだけで使用できる。NDPは特に数値関数の計算速度の向上に効果がある。

以下，コンピュータ言語で扱う数の内部表現の代表的な例を整数型と実数型に分けて要約する。詳細については文献［9］を参照されたい。

（A）整数型の内部表現
16ビットMPUでは機械語の一命令で16ビットの四則演算が実行されることから，高級言語でも16ビットを標準的な整数単位としている。2進整数表
3）8ビットコンピュータとの互換性を考慮して8ビット漬算も可能である。
4）整数型の場合，ここでいう有効桁数というあは，整数として揖える数の大きさの範囲である。整数の最大，最小が決められると，その範囲内で計算誤差はない。

一36一

1989年6月 藤村公男：16ビットパーソナルコンピュータの演算処理機能

示では0または正のみの符号なし整数と，正負の符号つき整数とを区別しなければならない。nビットの符号なし整数の範囲は0から2n‑1までである。16ビット長では0から2'6‑1 ・・ 65535までとなる。一方，符号つき整数はやや複雑となるが主として次の二つの表現方式がある。一つの方式は最上位ビット（MSB）を符号ビットとして， MSBが0のときは正，1のときは負とし，残りのビットを通常の符号なし整数とする。この方式では正数は＋0〜2n‑L1，負数は一〇〜一（2n‑L 1）となる。0は＋0と一〇の2っある。この方式はわれわれが通常10進法で正負を表現する方式であり，単純で理解し易いが，演算が複雑になるという欠点を持っている。Pro

FORTRANとPro PASCALの4バイトinteger型はこの表示法を採用している。符号つき整数を表現するもう一つの方式は2の補数表示である。これはビット長が有限で，整数の値が決められた範囲内で有効であるというコン

ピュータハードの特徴を利用したものである。nビットの2の補数表示では0および正数は0〜2n‑1‑1，負数は一1〜一2n‑1の範囲となる。加減算
は簡単となり，効率の良い演算が可能である。Lattice Cはこの方式を採用している。いずれの方式もMSBを符号ビットとする（図1）。符号なしに
図1 16ビット符号つき整数の2進表示

15 ．

0

LMSB（最上位ビット）
S＝0：正 S＝1：負しろ，符号つきにしても，16ビットでは整数値としての値の範囲はせますぎるので，F‑BASId86を除く他の言語では16ビットの倍長となる32ビット（4バイト）整数を扱うことが可能である。

整数型の特徴は値が決められた範囲内で誤差は皆無で，演算速度も高いということである。従ってプログラムでは可能な限り整数型を使うべきである。しかし，絶対値の範囲が実数型に比べてせまいので，四則演算（特に乗除算）の結果，オーバーフロー（桁あふれ）になる危険性が大きい。

プログラム中で整数型を指定するためにはプログラムの最初で型宣言しな
37

徳山大学論叢

第31号

ければならない。以下，変数acが整変数であることを宣言する命令（型宣言文）を各言語ごとに示す。5）（1） F‑BASIC86

DEFINT・x…最初の文字がxで始まるすべての変数は2バイト符号つき
整数
ただし，変数名の末尾に型宣言文字％をつけると上記の宣言文がなくとも良い。型宣言文も型宣言文字もない変数はすべて自動的に後述の単精度実数型とみなされるので注意が必要である。

（2） Pro FORTRANINTEGER x ……ooは4バイト符号つき整数

INTEGER＊1x…）cは1バイト符号つき整数

INTEGER＊2x…xは2バイト符号つき整数
INTEGER＊4x…xは4バイト符号つき整数

上記のINTEGERとINTEGER＊4は同じである。またFORTRANの伝統に従って，宣言文がない場合は，最初の文字が1，J， K， L， M， Nである変数

はすべてINTEGER型とみなされる。
（3） Pro PASCALINTEGER・x…xは4バイト符号つき整数
これ以外の整数型を扱うときはTYPE宣言によって行う。Pro PASCALの型宣言の仕方は独特であり，型自身をプログラマーが定義できる。例えばTYPE int ＝ 一32767 ： 32767 ；VAR int x；とするとxは2バイトの符号つき整数となる。
（4） Lattice C

charκ…xは1バイトの符号つき整数6）5）通常，ソースプログラム中の命令文や変数名では，ローマ字の大文字，小文字の区別をしない。本稿でもその習慣に従う。
6）1バイト整数は数値としてよりも，文字列のASC豆コードとして使われることが多い。数値としての演算は可能である。

一38一

1989年6月 藤村公男：16ビットパーソナルコンピュータの演算処理機能

intκ・…・・xは2バイトの符号つき整数

10ng int…xは4バイトの符号つき整数

負数は2の補数表示に従う。また上記のそれぞれにunsignedをつけると符
号なし整数となる。

（B）実数型の内部表現
プログラミング言語における実数型とは，整数型では表わせない有理数や無理数を近似的に2進または16進で表現する型であり，科学計算では不可欠である。1／2（o．5）や1／4（o．25）のように小数点つき2進表示の有限桁で表現できる数の他に1／3やV7，円周率のように2進でも10進でも無限桁まで続くものがある。注意すべきことは10進法で有限桁の数も2進法で無限桁となる場合もあるということである。例えば10進の0．6は2進表示で0．100110011…と無限に続く。もちろん，現在のコンピュータでは無限桁の数を表現処理することは不可能である。従ってどんなプログラミング言語
でも，実数型をいかに2進または16進の有限桁で近似表現するかが問題となる。上記の0．6を2進8桁0．10011001で近似すると，これは10進で0．59765625に変換され，二値0．6に比べると約0．4％の相対誤差となる。一般に有効桁の多い数値を有効桁の少ない数値に変換することで生じる誤差を「丸め誤差」といわれている。われわれが日常行っている四捨五入もこの丸め誤差を承知の上で，有効桁数を少なくして，表現を簡単にしていることに他ならない。

計算の精度を高めるためにはコンピュータの内部表現のビット長を長くする必要がある。しかしビット長を増すほど，10進から2進への変換およびその逆変換，また演算に要する時間が増加する。整数型の場合と同様に，実数型でも精度（有効桁数）と計算速度は相反する。丸め誤差は10進から2進に変換するときに生じる。その誤差をなくするために，コンピュータの内部表現
を10進表示とし，演算も10進で行うソフトもある。しかし，これはメモリの利用効率が悪く，演算速度も遅くなるのであまり一般的でない。

実数型では，精度（有効桁数）の他に，数の絶対値がどの程度まで可能か
一 39

徳山大学論叢第31号

が問題となる。絶対値の大きさに幅がある数を表現するには，以下のように指数表示を使う。

±MxRE （2‑1）±は数の正負符号，Mは仮数部， Eは指数部， Rは数表現の基数である。例
えば10進の123億は0．123×IOと表わす。0．001は0．1×10‑2と表現する。（2‑1）のような数表現を浮動小数点表示という。2進の浮動小数点表示では基数はR＝2となるので±Mx 2E ・ （2‑2）となる。±，M， Eは2進表示とする。コンピュータ言語における実数型はすべてこの浮動小数点表示に従う。問題は全体としてのビット長と，M， Eにどれだけのビット長を割り当てるかである。

以下，実際のコンピュータで採用されている代表的な浮動小数点形式を示す。ここで掲げた三形式とも，32ビット（4バイト）長の単精度型と64ビッ
ト（8バイト）長の倍精度型がある［9］。

（i）IBM社の汎用コンピュータIBM 360，370シリーズおよびこれらの互換機で採用されている実数型で，16進表示である。

○単精度（4バイト）

31 24 23 O
有効桁は10進で約7桁○倍精度（8バイト）
63 56 55o

有効桁は10進で約16桁
単精度，倍精度とも，指数部は7ビット長で値の範囲はV≒±5．4×10『79〜±7．2×1075である。

一40一

1989年6月 藤村公男：16ビットパーソナルコンピュータの演算処理機能（li）マイクロソフト社の2進浮動小数点形式

これはマイクロソフト系BASICで採用されている。米国ではIBM PC，

アップル，日本ではNEC，富士通，日立などのパソコン用BASICがこれに従っている。

○単精度（4バイト）

31 23 O
有効桁数は10進で約7桁○倍精度（8バイト）

63 55

o

有効桁数は10進で約16桁
両精度とも，指数部は8ビットで，値の範囲はV≒±2．9×10‑39〜±1．7
× 103s（iii）IEEE規格の2進浮動小数点形式
これは米国電気電子学会（Institute of Electrical and Electronics Engi‑

neers）がパソコンの分野で統一的に使用するよう定めた規格である。 Pro

FORTRAN， Pro PASCAL， Lattice Cはこれを採用している。また数値演算コプロセッサ8087はこの形式を採用している。

○単精度（4バイト）

31 23 22 O
指数部は8ビット，仮数部は23ビット。
有効桁数は10進で約7桁。値の範囲はV≒±1．1×10‑38〜±3．4×1038

0倍精度（8バイト）

一 41

徳山六学論叢

第31号

指数部は11ビット，仮数部は52ビット。
有効桁数は10進で約16桁。値の範囲はV≒±2．2×10‑308〜±1．8×10308

1EEE規格の（i），（ii）と異なる特徴は，倍精度が単精度より指数部のビット長が大きくなることである。その結果，倍精度の絶対値の幅が大きくなる。

上記の三形式とも，有効桁はあまり差がないが，10進表示での指数部の差が大きいことに注意すべきである。以上の浮動小数点形式は実数型の代表的な例であり，これ以外の形式もある。Borland社のTurbo PASCAL V 1．00 A

では，実数型は単精度，倍精度の区別はなく，6バイト長のみである。一般にプログラミング言語はなるべくハードウェアに依存しない「移植性」の高い仕様が期待される。従って数の表現，特に浮動小数点形式のようなハード
ウェア依存性の強いものは，言語の標準仕様で厳密に規定されていない場合が多い。例えばPro PASCALの倍精度実数型は標準PASCA：しの拡張仕様となっているに過ぎない［6］。また標準C言語とされているカーニハン，
リッチのテキスト［10］でも単精度'（float）と倍精度（double）の区別はあるが，それらの具体的なビット長や表示形式は決められていない。どのような表示形式を採用し，どのように演算処理するかは，コンパイラ言語の処
理系（コンパイラ）の作成者に任されている。Lattice Cでは，入出力およびシステム関数などの基本ライブラリと浮動小数点演算ライブラリが別個に存在し，後者は浮動小数点演算が必要なときオプションとしてリンクすれば良い。従ってこの浮動小数点演算ライブラリをアセンブラで作り直すことによって，利用者の目的に最適なライブラリとすることが可能である。例えば倍精度より軍にビット長を大きくして超倍精度演算が可能なライブラリをア
センブラで作り，それをLattice Cにリンクして利用することなどである。

最近発表されたBorland社のTurbo C 2．0では80ビット（10バイト）長の実数型で3．4×10‑4932〜1．1×104932が可能である。

以下，変数xが単精度型，変数Nが倍精度型であることを指定する宣言文

一42一

1989年6月 藤村公男：16ビットパーソナルコンピュータの演算処理機能

を各言語ごとに示す。

（1）F‑BASIC86（マイクロソフト系BASIC規格）

DEFSNG x……xで始まる変数は単精度型

DEFDBLッ……yで始まる変数は倍精度型
宣言文がなくても，oc／， y＃のように型宣言文字／（単精度），＃（倍精度）を変数の末尾につけても良い。型宣言文も型宣言文字もない変数
はすべて単精度とみなされる。

（2）Pro FORTRAN（IEEE規格）

REAL＊4x（またはREAL x）
REAL＊8．y（またはDOUBLE PRECISION）y）（3）Pro PASCAL（IEEE規格）
real x ；longreal y ；（4）Lattice C（IEEE規格）
float ac ；doubleごソ；

3 浮動小数点演算の精度

実数型の場合，単精度にするか倍精度にするかはプログラムの目的に応じて決定しなければならない。高精度で指数部の範囲が大きい場合はIEEE規
格の倍精度を使うべきである。ただし計算時間とメモリ使用量は単精度に比べて多くなる。

単精度か倍精度かを決定するもう一つの要因は演算の結果生じる「桁落一 43 一

徳山大学論叢

第31号

ち誤差」である。この簡単な例は2次方程式αx2 ＋ bx＋c‑0の根を計算する場合である。2実根はx。＝（一b＋〉面）／（2α），κ1一（一b一〉廊）／（2α）となるが，b＞0かっb2＞14αclの場合にoc。の計
算に桁落ち誤差が生じる。このときbと〜面とはほとんど同じ値となりx。の分子である一b＋一v／5b i'＝'4EiZ14 acの有効桁は失なわれる。この解決法はx。＝一2c／（b＋V薦）のように変形すれば桁落ち誤差は防げる。この例のように簡単な場合は良いが桁落ち誤差はいろいろな場合に生じ，時によっては桁落ち誤差が何回も積み重なって最終結果がほとんど無意味な値になることがある。

桁落ち誤差が最も深刻でその解決が困難な問題は逆行列の計算である。行

列演算は単に自然系分野に限らず，経済経営の分野でも重要な役割を演ずる。W。レオンチェフ等によって展開された有名な経済学の投入・産出原理［11］やオペレーションズリサーチにおける線型計画法［12］は逆行列演算と深く関係している。

N次連立線型方程式は行列形式で

Ax ＝＝c （3m1）と書ける。ooはN次の未知数ベクトル， cはN次の定数ベクトル， AはN行
N列の係数行列である。一般にこれらの要素または成分は実数型として扱われる。問題はAの1＞2個の要素と。のN個の成分がデータとして与えられたとき（3‑1）を満足するxのN個の成分を求めることである。1＞一2次の
ときは

e ＝＝
'

'

CC

sc ＝＝

a21 a22κκ

all a12

A一

（3‑2）

︷

とすると（3‑1）は次の2元連立方程式となる。
all Xl 十 a12 X2＝ Cl（3‑3）
a21 X 1 十 a22 sc 2 ＝＝ C 2

一44一

1989年6月 藤村公男：16ビットパーソナルコンピュータの演算処理機能

︷

この解は

sci ＝＝ （Ci a22一 c2 ai2）／det A

sc 2＝（一ci a2i十 c2 an）／det A（3‑4）

det A＝＝ aii a22‑ai2 a2i （3‑5）
det、4はAの行列式である。（3‑4）の解でdet・Aが分母にくるのでdet
Aが0になると解は無意味となる。det・Aが0にならないまでも桁落ち誤差の生じる可能性は充分にある。この問題は一般的なN次の場合でも同様である。一般式である（3‑1）の解はAの逆行列をA一'とすると形式的にとなる。A

x＝A一'C （3‑6）

1を求めるときにdet・Aが分母に現われ， det Aが0に近くなるとオーバーフローかまたは桁落ち誤差が発生する危険性がある。この困難を避けるためには計算時間の多くなることを無視して，なるべく高精度な計算をすることである。

実際のコンピュータを使って，逆行列の計算を行い，どの程度の誤差が生じるかをテストする一つの方法は以下に示すヒルベルト行列の逆行列を計算することである［13］。

N次のヒルベルト行列Aの各要素a、iはaiJ・ ＝ 1 ／（i十 7'一1） （3‑7）と定義される。ただし，，ノー1，2，3，…，1＞。この行列の各要素は正で，Nが大きくなるほど，近接する要素の差が少なくなるので逆行列の計算の際，桁落ち誤差が発生する。本稿では単精度型と倍精度型の両方で計算を行い，両者の差を比較する。、4一'の一般式は代数的に決められないのでA‑1の計算

結果と真赤を直接比較することは出来ない。そこでA'A ・＝1（1はN次単位行列）となることを利用して，計算されたA‑1からA一'Aを計算し，更に1を減じ，その結果が0からどの程度離れるかを調べた。すなわち

A‑iA‑1＝一〇 （3‑8）（0はゼロ行列）一 45

徳山大学論叢

第31号

表1 A'iA‑1の要素の中で絶対値が最大のもの（三値は0）（AはN次ヒルベルト行列。F‑BASIC86による）

N

2

3

5

．単精度

0

7．6×10一6

1．4×10‑2

6

192

倍精度

0

2．2×10‑16

2．3×10‑12

2．2×10 9

2．1×10‑4107の左辺を計算機で計算し，その各要素が真値0からどの程度ずれるかをみる。

逆行列A‑1を標準的で最も一般的に用いられぜいるGauss‑Jordan消去法［14］によって求める。N次正方行列の要素は1＞2個あるので高次ではこれらのすべての要素を出力することは困難である。従って，、4'A‑1の1＞2個の要素の中で絶対値が最大となるものを抜き出して出力した。F‑BASIC86による結果を表1に示す。この結果は他の言語でもほとんど同じである。表1からわかるとおり単精度ではN‑3次元から，倍精度では5次元から既に桁落ち誤差が生じている。

ヒルベルト行列はコンピュータによる桁落ち誤差をチェックするにはきびしすぎるかも知れない。通常の線型計算では，このようなことはあまり起らないであろう。しかし，桁落ち誤差が生じているか否か一般的に判断することは極めて困難である。N次正方行列の逆行列を求めるには2N3回程度の実数演算が必要である。従ってNが大きくなると計算時間が無視し得なくなる。演算ごとに桁落ちが生じているか否かを判断するプログラムも可能だが，そのためには更に多くの余分の時間が必要となる。コンピュータでこのような線型計算を行うときは少なくとも単精度と倍精度の両方で行い，両者の結果を比較して，桁落ち誤差の少ないことを確認する必要がある。倍精度は単精度より2倍のメモリ領域を必要とする。N次正方行列はN2個の要素があるので，高次になると必要なメモリサイズは急速に増大する。

4 四則演算と数値関数の計算速度

ここではクロック周波数が8MHzの80186をベースとした：FM 16βを使っ

一46一

1989年6月 藤村公男：16ビットパーソナルコンピュータの演算処理機能て各言語ごとに演算処理時間がどの程度異なるものかを比較検討する。浮動小数点演算については数値演算コプロセッサ80187の効果も検討する。
演算時間を決定するためには通常，演算を繰り返し実行し，その間に経過した時間を測定すれば良い。演算時間は計算する数の型によって異なる。従って四則演算の場合は短整数型（2バイト），長整数型（4バイト）？），単精度実数型（4バイト），倍精度実数型（8バイト）の各々にっき行う。また数値関数の場合は単精度と倍精度の両方で行う。実行時間の決定はコンピュータに内蔵のシステム時計を利用する。MS‑DOS V3．1では秒以下を測定できないので計算速度の信頼性を高めるたあ繰り返し数を数万回から数百万回にする必要がある。

リスト1は2バイト整数の加算を8万回実行し，その計算時間を秒単位，
1回当りの計算時間をミリ秒単位で出力するF‑BASIC86のプログラムと実行結果である。リスト1の計算時間はFOR〜NEXTループと代入文に要する時間も含めたもので，純粋な加算五＋召の時間ではない。しかし実際のプログラムでは純粋な加算のみということは少なく，代入文など他のルーチンが伴うのが普通である。本稿ではループや代入文も含めた計算時間とする。

リスト1と同様のプログラムを他の言語でも作成し実行した結果を表2（整数型），表3（実数型）に示す。演算の繰り返し回数は全体の計算時間が100秒以上となるように決定した。従って表の結果は1％程度の誤差があるものと考えられる。表2でループというのは代入文も演算もないループのみに要する時間である。代入文p一αというのはループと代入文のみで演算がない場合の時間である。従って表の各演算時間から代入文p ・αの時間を差し引いた値が演算のみの時間となる。整数除算α／bはPro FORTRANと1．attice Cの表記法であるが， F‑BASIC86の場合はα￥b， Pro PASCALの場合はα div bとなる。平均は四則演算の平均値である。比率はF‑BASIc86
インタプリタの平均を100とした場合の，それぞれの平均の比率である。
7）ここでは仮に，2バイト，4バイトの整数型をそれぞれ短整数型，
記すことにする。

一47一

長整数型と

徳山大学論叢第31号・

リスト1 整数型加算の計算．時間を出力するF‑BASIC86プログラムと実行結果

10 DEF工NT N，図20 DEF工NT A，B，P
30 A＝1：B＝2

40 NMAX＝1000：MIVIAX＝80
50 T1＝T工ME
60 FOR N＝1 TO NIVIAX

70 FOR M＝1 TO MIVIAX

80
P＝A＋B
90・ NEXT IVI
100 NEXT N

110 T2＝T工tVI E

120 PR工NT i' integer P＝A＋B＝lt；P；17．；）O PR：NT
time＝ ；T2‑Tl；

sec；：140 PR工NT 年OOO＊（T2‑T1）／NlvlA×／lvlMA×；itmsii

150 END

RUN

integer P＝A＋B＝ 3 time一一

114 sec 1．42S ms

表2 整数型四則演算1回当りの計算時間（ミリ秒）型

演 算

短整数︵2バイト︶

F‑BASIC86

Pro

インタプリタ

コンパイラ

ルー プ

0．26

0．12

P＝α

0．91

0．19

P＝α十わ

1．4

Aρ＝α一b：Lattice

FORTRAN PASCA：L@C

一〇．00700．00610．00800．00980．01040．0115・

0．28

0，023

0，024

0，013

P．4

O．29

O，023

O，024

O，013

o＝α＊b
o＝α／b

P．4

O．29

O，016

O，016

O，017

P．4

O．30

O，058

O，057

O，020

平

1．4

0．29

0，030

0，030

0，016

P00

Q1

@2．1

@2．1

@1．1

一

一

0，014

0，015

0，015

0，028

0，029

0，018

一一

一一

O，029

O，030

O，018

Aρ＝α＊わ

O，067

O，066

O，050

o＝α／b

O．21

O．21

O．16

0，084

0，084

0，062

@6．0

U．0

S．4

均

? 率
一長整契4ノ・ミイと

P＝α、ρ＝α十わAρ＝α一わ

平

均

? 率

一一

一一

一48一

1

1989年6回忌藤村公男：16ビットパーソナルコンピュータの演算処理機i能表2の結果を要約すると，整数演算1回当りの計算時間は（1）F‑BASIC86インタプリタに比べると， F‑BASIC86コンパイラは約

1／5，Pro FORTRANとPro PASCALは約1／50， Lattice Cは約1／100である。

（2）長整数（4バイト）は短整数（2バイト）より3〜4倍。

（3）Pro FORTRANとPro PASCALはほとんど同じ。

表3（実数型四則演算）の（）内の値は数値演算コプロセッサ付の計算時間である。表3の結果を要約すると（1）F‑BASIC86インタプリタに比べて， F‑BASIC86コンパイラは約1／3〜1／4，Pro FORTRANとPro PASCALは約1／12〜1／6，1．attice
C Ca 1／3 一一 1／4．

（2）倍精度は単精度より1．2〜2倍（Lattice Cを除く）。

（3）Pro FORTRANとPro PASCALはほぼ同じ。

NDPの効果は（4）F‑BASIC86では効果なし。

（5）他の言語では単精度で1／2，倍精度で1／4。
Lattice Cの単精度は倍精度より計算時間が大きくなるのは問題である。カー二・ハン・リッチー［10］p．198の付録A：C参照マニュアルによると「Cでは，すべての浮動小数点演算は倍精度で実行される」とある。ということ
は単精度データは倍精度に変換され，倍精度演算を行い，結果を単精度に変換するという手間のかかる手続きを経ているものと考えられる。型変換を2回行うことに余分の時間がかかり，このような結果になるのであろう。従ってC言語では単精度にしても，これはメモリが倍精度より節約できるというだけで，計算速度が高くなることはない。

四則演算の実用的な評価をするために，数値積分の簡単な例8）を作成し実行した。プログラムの目的は積分値が既知の定積分
8）これは文献［9］で8087の効果を調べるために使われた。
一 49

徳山大学論叢第31号

表3 実数型四則演算1回当りの計算時闇（ミリ秒）（）内はNDPつき型

演算

F‑BASIC86

Pro：Lattice

FORTRAN PASCA：L

@C

単

コンパイラ

0．93

0．19

0，057

0，057

0，015

P＝α

i0．92）

i0．19）

i0．044）i0．046）i0．015）

P＝α十わ

1．6

0．45

0．14

0．15

0．41

i1．6）

i0．45）

i0．066）i0。066）i0．27）

P＝α一b

1．6

0．47

0．15

0．15

0．42

i1．6）

i0．47）

i0．066）

i0．068）i0．29）

1．7

P＝α＊わ

0．47

0．14

0．14

0．67

i1．7）

i0．47）

i0．068）i0，066）i0．27）

P＝α／b

2．3

0．46

0．15

0．15

0．82

i2．3）

i0．46）

i0．080）i0．082）i0，30）

平

均

1．8

0．46

0．15

0．15

0．58

i1．8）

i0．46）

i0．070）i0．070）i0．28）

比率

精

インタプリタ度︵4バイト︶

100

26

8．3

8．3

32

i100）

i26）

i3．9）

i3．9）

i16）

倍精度曾ノ・ミイと

0．95

0．20

0，074

0，078

0，016

P＝α

i0．94）

i0．19）

i0．046）i0．049）i0．016）

P＝α十わ

1．6

0．61

0．23

0．23

0．27

i1．6）

i0．62）

i0．069）i0．070）i0．099）

P＝α一b

1．6

0．64

0．25

0．26

0．25

i1．6）

i0．64）

i0．069）i0．069）i0．099）、ρ＝α＊b

P＝α／b
平

均

比

率

1．9

0．74

0．34

0．34

0．54

i1．9）

i0．74）

i0．075）i0．077）i0．11）

3．5

0．81

0．40

0．39

0．65

i3．5）

i0．80）

i0．082）i0．086）

i0．11）

2．2

0．70

0．31

0．31

0．43

i2．2）

i0．70）

i0．074）i0．076）i0．10）

122

39

17

17

24

i122）

i39）

i4。1）

i4．2）

i5．6）

50

1989年6月 藤村公男：16ビットパーソナルコンピュータの演算処理機能

O．25

f dx ＝4（tan‑11‑tan‑10）一π（円周率）

O．0625 十 x2

（4‑1）

0

の左辺を台形公式によって数値積分し，結果を真値πと比較する。計算時間を秒単位で表示する。台形公式では，積分区間をn等分し，微小区聞を台形と近似して，その微小台形面積を積分区間にわたって積算する。数式で表わ
すと積分区間（x。，oc・）をn等分し，きざみをhとすると
h＝（κn‑Xo）／n， Xi＝hi十Xo

（4‑2）

ロ
8＝h（！（Xo）十・2Σ！（oci）十！（κn））／2

（4‑3）

i＝1

8が求める積分値である。！（sc）は被積分関数で（4‑1）の場合！（x）＝1／（0．0625十κ2）

（4‑4）

である。

計算は倍精度で行った。Lattice Cによるプログラムとその実行結果をリ

スト2に示す。ここではn＝30000とした。

リスト2と同様のプログラムを他の言語でも作成し実行した結果を表4に示す。表4の最後列のPro FORTRAN77は引用文献［9］p．319の結果である。これはNEC9801E（MPUは8086，クロック周波数8MHz）で言語は
Prospero Software社のPro FORTRAN77による。台形公式による積分値
はどの言語でも10桁まで真値（π）と一致した。従ってこのプログラム例では言語による計算精度の差はない。NDPの効果はF‑BASIC86にはないが，他の言語では4〜10倍の効果がある。

平方根や三角関数指数関数などの初等関数は科学計算で頻繁に使われる。ほとんどの高級言語にはこれらの関数を計算するルーチンが数値関数ライブラリとして用意されている。ここでは四則演算の場合と同じ方法で関数計算の時間を計測した。結果は表5のとおりである。この計算時間には四則演算のときと同様，ループと代入文の時問も含めてある。ただし，sin（oc），
cos（x）， exp（x）， tan‑1（x）に対してはん一1．0とし， A／王，10g（oc）に対し51

徳山大学論叢第31号

リスト2 台形公式による数値積分（Lattice・Cによる倍精度， NDPつき）

1：／＊

2： ＊

numericaし integraしs‑modeし， doubしe Precision with 80187

3： ＊by lく．Fujimura，

02／12／89

4： ＊／

5：＃incLude

stdio．h

6：篇incしude'math．h7：ninctude

time．hB：

9：extern char ．NDP；

︸︸︸

10：

11：doubしe f（x）
doubLe x；｛ return（ 1．0／（0．0625＋x＊x） ）；12：．doubしe f i（x） doubしe ×；｛ return（ 4．0寧atan（x／0．2S） ）；13：doubしe fdi（a．b） doubしe a、b； ｛ rett」rn（ fi（a）一fi（b） ）：14：

15：mainO

16：｛ int i，n；17： doubしe h．×0．xn．×、s；

0
21
22
23
24
25
26
27
28
29
20
31
32
33
34
35
36
37
38
39
3

18：

unsigned しong stime，etime；

19；

コロ
ロコ
のコ

printf（8。＊＊ numericaし integraし 9etting Pai ＊＊￥ntl）；
printf（1巳しattice C s‑modeし doubしe Precision￥n￥n'。）；printf（ 一NDP＝＝7．d ，一NDP）；

ロ

if（．

NDP＝'一一〇） printf（

without 80187￥n￥n）；

eしse printf（tlwith 80187￥n￥01量）；

ロコ
ロロ

n＝30000；

printf（xo＝o．o；xn＝O．2s；true vaしue＝ '^．18．ISf￥n，fdi（xn，xO））；

ロ
ココ

h＝（xn‑xO）／（doubしe）n；×昌×0十h； s昌0．0；

コtirne（＆stime＞；

の

︸

の
ココ

for・（i＝1； iく＝n‑1； i＋＋） ｛ s ＋・昌 f（x）；

x ＋＝ h；コのs＝h＊（ f（×0）＋2．0＊s＋f（xn） ）／2．O；コロけtime（＆ettme）；
printf（ltcaしcuしated vaしue ＝ ％18．1Sf￥n61，s）Iロ

printf（ロコ

printf（used time ＝ Z3d sec

n＝7．d￥n，n）；，etime‑stime）；︸

ゆロ

exit（O）；

tm200＊＊ numericaし integr、aし 9etting Pai ＊＊．

Lattice C s‑modeしwith BO187'

pNDP＝＝1

doubしe precision＝ 3．141S926S3S89793
caしcuしated ＞aしue ＝ 3．141592653404718
true ＞aしue

used time＝ 14 sec ． n．30000

一52一

1989年6月 藤村公男：16ビットパーソナルコンピュータの演算処理機能
表4 台形公式による代数積分（倍精度）、（ ）内はNDPつき

F‑BASIC86

Pro

インタプリタ

コンパイラ

計算時間

259

79

@（秒）

i259）

i79）

比

率

Lattice

FORTRAN PASCAL
40
i8）

@C

Pro
eORTRAN 77

49

62

65

i12）

i14）

i6）

24
i5）

i2）

100

31

15

19

i100）

i31）

i3）

i5）

25

（Pro FORTRAN77は文献［9］より引用）てはx‑2．0とした。注意すべきことは，関数計算のときは独立変数）cの値によって計算速度が異なることである。コンピュータによる関数計算では，

通常，何らかの近似公式（マクローラン展開など）に従って計算する。しかし，その近似公式はxの値の大きさによって異なる。従って計算し得るxの値とその大きさに何らかの制限がついており，またxの大きさに応じて計算時間も違ってくる。各言語ごとにどのような近似公式またはアルゴリズムを採用しているかは，筆者には不明である。そこでxの大きさをいろいろ変えて関数計算を実行した結果，log Xの計算時間がかなり変化したが，他はあまり変化しなかった。表5の関数計算時間はあくまでも前述のooの値によるものであることを断っておく。

数値関数の計算時間について要約すると（1）F‑BASIC86インタプリタに比べて，同コンパイラは約1／2， Pro

FORTRANとPro PASCA：Lは約1／2．5（2）倍精度は単精度より約3〜4倍（Lattice Cはほとんど同じ）（3）Pro FORTRANとPro PASCALはほとんど同じ。

NDPの効果は（4）単精度と倍精度の差はない（Lattice Cを除く）。

（5）単精度ではF‑BASIC86は約1／2，Pro FORTRANとPro PASCAL
では約1／6，：Lattice Cでは約1／20。

（6）倍精度では，F‑BASIC86は約1／8〜1／10， Pro FORTRANとPro
PASCALは約1／20， Lattice Cは約1／27。

一53一

第31号

徳山大学論叢．表5 数値関数1回当りの計算時間（ミリ秒）（）内はNDPつき
型

数値関数

疵
単精度︵4ノ・ミイト︶

sin（κ）

COS（κ）

taガ1（κ）

exp（κ）

109（κ）

平

均

比

率

疵
倍精度︵8ノ・ミイト︶

sin（コσ）

COS（κ）

tan皿1（κ）

exp（κ）

109（」σ）

平

均

比

率

F‑BASIC86
インタプリタPro

コンパイラ：Lattice

FORTRAN PASCA：L

1．7．

0．83

0．13

i1．3）

i0．45）

i0．11）

0．12
i0．098）@ C
4．8

i0．26）

3．3

1．5

1．8

1．9

10．6

i1．6）

i0．75）

i0．35）

i0．36）

i0．60）

3．6

1．6

2．0

2．1

11．6

i1．6）

i0．76）

i0．39）

i0．39）

i0．59）

5．5

2．3

1．5

1．6

10．2

i1．4）

i0．55）

i0．20）

i0．20）

i0．38）

4．6

2．1

1．9

1．8

8．3

i1．6）

i0．71）

i0．26）

i0．25）

i0．43）

2．2

0．96

1．0

0．99

2．7

i1．3）

i0．47）

i0．11）

i0．11）

i0．26）

3．5

1．5

1．4

1．4

8．0

i1．5）

i0．62）

i0．24）

i0．23）

i0．42）

100

44

40

40

230

i42）

i18）

i6．8）

i6．7）

i12）

5．5

5．1

0．77

0．75

4．7

i1．3）

i0．50）

i0．11）

i0．10）

i0．14）

9．4

6．6

6．3

6．2

10．5

i1．6）

i0．78）

i0．37）

i0．37）

i0．49）

8．7

6．0

6．3

6．4

11．4

i1．6）

i0．79）

i0．39）

i0．39）

i0．47）

28．4

12．0

6．0

6．0

10．0

i1．4）

i0．56）

i0．21）

i0．20）

i0．18）

9．1

5．9

5．0

5．0

8．2

i1．6）

i0．73）

i0．26）

i0．26）

i0．32）

5．3

2．5

3．1

2．8

2．6

i1．4）

i0．50）

i0．11）

i0．11）

i0．15）

11．1

6．4

4．6

4．5

7．9

i1．5）

i0．64）

i0．24）

i0．24）

i0．29）

316

181

131

129

226

i42）

i18）

i6．8）

i6．8）

i8．3）

一54一

1989年6月 藤村公男：16ビットパーソナルコンピュータの演算処理機能（6）からNDPは倍精度関数で効果が顕著である。Lattice Cの単精度が倍精度より遅くなるのは四則演算のときと同じ理由による。また表からも分かるとおり，NDPが付いているときは単精度関数計算はNDPの倍精度計算に従っている。

四則演算のときと同様に，台形公式によって数値関数の積分を実行し，各言語とNDPの効果を調べた。積分は［9］ズハ／｛・i・ω＋…ω＋t・n‑i（・）＋・xpω＋1・g（・）｝dx

XO［一一…（・）＋・i・ω＋・t・n一・（x）一告，・9（1＋・・）

オハ＋・xpω＋・（1・9・一1）］

（4‑5）

Xe

で，左辺を数値積分し右辺の真値と比較する。計算時間も表示する。κ。＝0．5，xn＝1．5， n‑3000とした。各言語ごとの結果は表6となる。

表6からNDPの効果が大きいことが分かる。精度は言語によってほとん
ど変わらない。

表6 台形公式による関数積分（倍精度）（）内はNDPつき
F‑BASIC86

計算時間
@（秒）

比

率

Pro

インタプリタ

コンパイラ

213
i30）

114

FORTRAN PASCAL

i15）

100

54

i14）

i7）

Lattice@C

Pro
eORTRAN 77

89
i5）

89
i5）

42
i2）

42

70

63

i2）

i3）

i2）

149

135

i6）

i4）

（Pro FORTRAN77は文献［9］より引用）

5 Cray社スーパーーコンピュータとの比較

D．V． Andersonはパソコンからスーパーコンピュータまで8機種のコンピュータにっき，計算速度対価格を比較している（表7）。表7の最下行は一55一

徳山大学論叢第31号

本稿で得た計算時間から概算して記入したものである。通常，大型コンピュータでは演算速度を表わすのにMflopsという単位を用いる。これは1秒間に可能な浮動小数点演算の回数を100万単位で表わすものである。もちろん，

Mflopsが高いほど，高速となる。スーパーコンピュータでは数十Mflopsに達する。表の計算速度はFORTRANによる64ビットの実数型演算である。
表の結果はどのような条件で得られたかの詳細は不明であるが，本稿で示したFM 16β（NDPつき）の結果からMflopsおよびRelative Costを概算してみる。FM 16βシステム（ハードディスクなし， NDPつき，1MB RAM実装，中級のディスプレイとプリンタ，ソフトはMS‑DOSとPro FORTRAN）の標準価格は65万円程度である。これを1＄一130円で換算すると5k＄になる。Memoryは1mである。計算速度として，ここではPro FORTRAN
で得た結果を使う。倍精度実数型の四則演算1回当たりの計算時間の平均値

は表3から0．074msである。この値はNDPつきの場合である。これから1秒当たりの演算回数は1／0．074×10'3≒14000／秒。従って計算速度は0．014Mflopsとなる。 Relative Costは5k＄／0．014Mflops＝ O．36（106＄／Mflops）と算出される。これらのデータは表7の最下行に示してある。
表7 Ratio of computer cost to speed compared（文献［2］より引用）

・・m・・ter畠温（＄）糠：3
VIC 20

1BM‑PCNo 8087 Chip

IBM‑PC

with 8087 Chip

VAX 11／750
DEC KL‑10
CDC‑7600
Cray‑IS

Cray XMP（2 C'PU's）

FM 16 B
with 80187 Chiplk

MflopsRelative Cost

with FORTRAN（10G ＄ ／Mflops）

20 k

O．OOOOg

11．1

4k

256 k

O．OOO37

10．8

sk

256 k

0．0073

0．68

80 k

8m

0．057

1．4

480 k

1．lm

0．18

2．7

4．5m

3．7m

4．6

0．98

llm

16m

18

0．61

15m

32m

53

0．28

5k

lm

O．014

汎用コンピュータはディスクドライブを含む初期コストである。
最下行のFM 16βは本稿のデータである。

一56一

O．36

1989年6月 藤村公男：16ビットパーソナルコンピュータの演算処理機能

FM 16βのMflopsは最高速のスーパーコンピュータCray XMPと比較して問題にならないが，Relative Costは同水準になるということは注目すべきであろう。ただし，コンピュータの価格を如何に算出するかは疑問な点が多いので表のRelative Costをそのまま信用することは問題である。一般に大型コンピュータは管理費などのオーーバーヘッドは無視できない。それに対してパソコンはこのような費用はほとんどないといってよい。

パソコンのMflopsはスーパーコンピュータに比べて決定的に低い。しかし，スーパーコンピュータで開発されている並行磐理方式がパソコンにも適用されればMflopsの大幅な改良がなされるはずである。

6 結論
クロック周波数8MHzの80186MPUをベースとしたFM 16βを使って，F‑BASIC86， Pro FORTRAN， Pro PASCAL， Lattice Cの数値演算機能（計算精度と速度）について比較検討した。計算精度については言語による差はあまりなかったが，計算速度に関しては大きな差がある。計算速度に関して最高の結果を示したのはPro FORTRANとPro PASCALである。この2っのソフトのメーカーは同一の会社なので，両言語とも，数値演算について同じライブラリを使っているのであろう。FORTRANはコンピュータが発明された当初から35年以上も発展，改良が重ねられてきた科学技術用言語であり，当然のこととはいえ，本稿でもそれが実証された。一方C言語は
主としてシステム開発用に考案された言語であり，浮動小数点演算の計算速度のことはあまり考慮していないようである。単精度実数型演算が倍精度より高速にならないというのは問題である。科学計算でも単精度で充分という場合が多いからである。F‑BASIC86インタプリタは他のコンパイラ言語より，計算速度の点で著しく劣る。しかし，もともとBASICはTSSによるプ
ログラミングの多数初等教育学に開発された言語である。BASICインタプ
リタはコンパイル，リンクという複雑な操作を必要としないし，パラメーター57一

徳山大学論叢

第31号

を何回も変えてプログラムを再実行することが容易である。F‑BASIC86コ
ンパイラはインタプリタと他のコンパイラ言語との中間といったところである。NDPの効果は特に数値関数計算で著しい。
価格対計算速度の点では16ビットパソコンは大型コンピュータとほぼ同じであるということは注目すべきと考える。Per Bakの主張するとおり，科学技術計算の分野でも一層パソコンを利用すべきであろう。

本稿の原稿を書き終った時点では，16ビットパソコンの主流は既に80286

MPUをベースにしたものに変わっている。クロック周波数も10，12，16

MHzとなった。16MHzは8MHzに比べて演算速度は約2倍となる。クロック周波数が16，20，25MHzの32ビットパソコンも，個人的に入手し易い
低価格で市場に出てきた。ハードの急速な進歩にソフトが追いつけないのが現状であろう。32ビットMPUの性能を完全に活用できるようなOSと高級言語の開発が期待される。

なお本稿では紙数の関係でRAMメモリの利用効率やメモリモデルについては言及しなかった。これらについては機会があれば別に発表したいと思っている。また32ビットパソコンを入手次第，その結果を報告したい。

引 用 文献［1］ Per Bak， Physics today， Dec． 1983， p． 25［2］ David V． Anderson， Physics today， Jul． 1984

［3］富士通，F‑BASIC86 V 3．1文法書［4］富士通，F一一BASIC86コンパイラ使用手引書［ 5 ］ Pro FORTRAN user's manual， LIFEBOAT INC． 1985
［ 6 ］ Pro PASCAL user's manual， LIFEBOAT INC． 1985
［ 7 ］ Lattice C Compiler user's manual， LIFEBOAT INC． 1985

［8］富士通，日本語MS‑DOS V3．1
［9］大貫広幸「数値演算入門」CQ出版社，1988［10］B．W．カーニハン， D． M．リッチー「プログラミング言語C」石田晴久訳，共立出版，1986［11］二階堂副包「現代経済学の数学的方法」岩波書店，1974［12］R．ドーフマン，P． A．サミュエルソン， R． M．ソロー「線型計画と経済分析］1」安井他訳，岩波書店，1972

一58一

1989年6月 藤村公男：16ビットパーソナルコンピュータの演算処理機能［13］W．R．ベネット「パソコンプログラム文科系のための問題演習」竹内他訳，現
代数学社，1983［14］J．M．マコーミク， M． G．サルバドリ「フォートランによる数値計算プログラム」清水留三郎訳，サイエンス社，昭和49年

59