アルゴリズム教室
Ｐｒｏｇｒａｍｉｎｇ ａｌｌ ｎｉｇｈｔ
コンピュータの世界は日進月歩で、今日覚えたソフトの操作方法が明日には役に立たなくなってしまいします。
ところが、変わらないものもあります。それが、アルゴリズムです。
アルゴリズムとは、問題を解く手順のようなものです。現在のコンピュータの仕組み（ノイマン型コンピュータで、物事を順番に解いていく）が変わらない限り（並列型や量子コンピュータ）永遠です。
たとえば、あなたが先生で生徒の順位をつけたい場合、成績の良い順に生徒を並び替えなければなりません。ここでは、ソートのアルゴリズムを必要とします。
このようにコンピュータを用いて何か問題を解こうとする際に、大きな役割を果たすのがアルゴリズムです。
アルゴリズムの中で一番多く使われるのが、「ソート（整列）」と「検索」です。
「ソート」には実にいろいろな方法があり、ざっと挙げてみると、選択法、挿入法、バブルソート、クイックソート、シェルソート、ヒープソートなどがあります。
選択法と挿入法やバブルソートは、アルゴリズムは単純ですが、速度が個数（ｎ）の２乗に比例して遅くなります。クイックソートやヒープソートは、アルゴリズムは複雑ですが、速度がｎ＊ｌｏｇ（ｎ）に比例します。このように、問題のサイズによってソートの方法を選ぶ必要があります。
「検索」にはソートほどバラエティにとんではいませんが、最初から順番に調べていく順次検索法と、あらかじめソートされたデータを調べる二分探索法（バイナリサーチ）が有名で、そのほかにハッシュ法（データ自体から検索の位置を作り出す方法）や二分木探索法（データを二分木に構成して探索する）などがあります。
ここでは、これらのアルゴリズムの解説は市販の書籍や雑誌に任せて（大きな本屋に行けばいっぱいある）書籍にあまり載っていない方法を紹介します。
コムソート
数あるソートの中で最速と言われるのがクイックソートですが、それに迫る速度でしかもどんなデータの条件でも安定（クイックソートは最悪の場合がある）しています。速度は、ｎ＊ｌｏｇ（ｎ）に比例します。出典は、「日経バイト」の記事です。
procedure ComSort();
var
i,j,gap,sw: integer;
n,buf: integer;
data array of integer;
begin
gap:=n; //データの個数
repeat
gap := Trunc(gap/1.3);
if gap<1 then gap := 1;
sw := 0;
for i:=1 to n-gap do
begin
j := i + gap;
if data[i]>data[j] then //データ交換
begin
buf := data[i];
data[i] := data[j];
data[j] := buf;
sw := sw + 1;
end;
end;
until (sw=0) and (gap=1);
end;
ラディックス（基数）ソート
クイックソートが最速と言われていますが、このソートはさらに速く、速度は個数に比例します。ただし、データが整数である条件が付きます。次に、１０進を基数とした２桁の整数での例を示します。
Number 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
Data 19 01 26 43 92 87 21 38 11 55 21 64 54 70 36
まず、0〜9までの10個の配列a[ ]を用意して、１桁目の数字に該当する場所に入れる。同じ場所に入るときには、リンクでつなぐ。
a[ ] 0 1 2 3 4 5 6 7 8 9
Data 70 01 92 43 64 55 26 87 38 19
21 54 36
11
21
続いて、同様に配列b[ ]を用意して、２桁目の数字に該当する場所に入れる。
b[ ] 0 1 2 3 4 5 6 7 8 9
Data 01 11 21 36 43 54 64 70 87 92
19 21 38 55
26
このアルゴリズムは、整数という制限やデータの個数と同程度の領域を必要とすることからあまり使われていません。整数という条件を越えて実数や文字列までできるように応用したＬＳＯＲＴ（リニアソート）が、「The Lifeboat PERSPECTIVE 1989.4」という小冊子に載っていましたが、詳しいことはわかりません。
乱数を順次抽出する方法
ある範囲の乱数の中からランダムに数値を取り出すとき、同じ数値を取り出さないように乱数を選ぶにはどうしたらよいのでしょうか。たとえば、壺の中に１から１００までの番号が付いた玉が入っていて、それをランダムに取り出して再び壺の中に戻さない場合を考えてみてください。単に１から１００までの乱数を発生させて、それに基づいて玉を取り出していたのでは、同じ値の乱数が発生したときには無駄になってしまいます。それに、残りの玉が少なくなるにつれて、壺の中に残っている玉の番号に該当する乱数が発生する確率は低くなり、時間がかかりすぎてしまいます。
アルゴリズムの考え方としては、一度発生した乱数はまだ発生していない他の乱数を示すようにして、同じ乱数が再び発生しても重複しないようにします。そして、乱数が１個発生するたびに乱数の範囲を１だけ減らしていきます。
procedure GetRnd();
var
i,r,n,guide,x: integer;
ptr array of integer;
begin
for i:=1 to n do ptr[i] := -1;
guide := 1;
randomize();
while (guide<=n) do
begin
r := random(n-guide+1) + guide;
if ptr[r]=-1 then x := r
else x := ptr[r];
//xが取り出された乱数
if ptr[guide]=-1 then ptr[r] := guide
else ptr[r] := ptr[guide];
guide := guide + 1;
end;
end.
このアルゴリズム乱数順次抽出法は、私がパソコン雑誌「Ｉ／Ｏ 1991.11」に発表したものです。
データの空を詰める方法
配列に入っているデータのいくつかを削除した場合に、空になった配列を１パスで詰める方法です。
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
a[ ] 19 01 26 43 92 87 21 38 11 55 21 64 54 70 36
x x x x x
state[ ] 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
state[ ] 0 0 0 0 -1 -1 0 0 0 -1 -1 0 0 -1 0
たとえば、配列 a[ ]の 5、6、10、11、14番目のデータを削除するとき、あらかじめ削除データの有無を判断する配列 state[ ]を用意して 0 を入れておきます。削除するデータに -1 を入れます。
次に、変数 counter=0 を用意して、state[1] から state[15] まで順に調べて行きます。state[ ]=0 のときには、counter を +1 にして、state[ ]=-1 のときには、counter はそのままにします。このとき state[ ]=0 のデータは counter の示す位置に移動させます。ただし、現在の配列の番号と counter の値が等しいときには移動はせずにそのままにします。
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
a[ ] 19 01 26 43 x x 21 38 11 x x 64 54 x 36
state[ ] 0 0 0 0 -1 -1 0 0 0 -1 -1 0 0 -1 0
counter 1 2 3 4 4 4 5
↑ ↓
←←←
これは7番目まで調べたときの状態を表したものですが、state[7]=0 で、しかも配列の番号７とcounterの値が等しくないので、データ a[7] を a[counter]=a[5] の位置に移動させます。
procedure PackData();
var
i,n,counter: integer;
begin
counter := 0;
for i:=1 to n do
begin
if state[i]=0 then
begin
counter := counter + 1;
if i<>counter then data[counter] := data[i];
end;
end;
end.
このアルゴリズムは自分で考えました。
小数を分数で近似する
有理数は分数で表すことができますが、無理数（√2など）や超越数（πやｅなど）は分数で表せず無理小数となります。そこで、このような小数を分数で近似することを考えました。
小数を分数で表すというのは、たとえば π＝３．１４１５９２６５３５８９…を 3＋10/71 = 3.1418511… （誤差0.008%）とか3＋16/113 = 3.1415929… （誤差1/200万） などのように近似することをいいます。
では、もとの小数の誤差が最小となる分数を発見する方法について考えてみます。考え方のポイントとしては、もとの小数（以後真値 Ｔ ）と近似分数 Ｂ／Ａ との誤差、Δ＝Ｔ−Ｂ／Ａ が小さくなる方向に探索を行います。
まず、最初の近似分数を真値より小さくなるように設定（ Ｔ＞Ｂ／Ａ ）します。このとき、Δ＝Ｔ−Ｂ／Ａ＞０ となります。次に、近似分数の分子を固定して分母から−１すると、近似分数は大きくなり真値との誤差Δは小さくなります。
そして分母をどんどん小さくしていけば、誤差Δもどんどん小さくなり、ついには Δ＜０ にないます。このときの分数の値もしくはこの寸前の分数の値が Δ＝０ に最も近く、誤差が最小となるので、これを現時点での誤差最小分数として記録しておきます。
次に、これ以上分母から−１してもΔは増える一方なので、分子から−１して近似分数の値を小さくして再び Δ＞ ０にします。そして、この操作を繰り返していけば、最後にもっともΔが小さな近似分数が記録されることになります。
以上の方法を「最小誤差探索法」と呼ぶことにします。この方法を具体例でもっと詳しく説明します。例としてπを取り扱うことにします。
まず、分数の分母を何桁にするかを決めます。ここでは３桁としてπの小数点以下４桁以降を切り捨てて、0.1415926…を0.141（整数部分は簡単にするため省略します）とします。よって出発分数は 141/1000 となり、これが最小分数となります。そして、真値に最も近い分数は 141/1000 〜 15/100 の間にあることになります。
さて、出発分数の分母から−１して分数の値を大きくします。そして、真値より大きくなりまで−１します。このときのひとつ前の分数またはこのときの分数が、真値との差が現時点での最小になるので、これを記録しておきます。この例では、141/1000 から出発して 141/999、1141/998 … となり、141/995 が真値より大きくなる分数です。この場合は、ひとつ前の 141/996 が真値との誤差が最小の分数です。
次に、現在の分数の分子から−１して分数の値を真値より小さくします。ここでは、140/996 となります。そして、これをまた出発分数としてこれらの操作を繰り返し、15/100になったら終了します。そうするうと、誤差が最小の分数が記録されていることになります。この例では、16/113 です。
もっと探索の数を減らす方法を考えました。前述の方法により、真値Tとの誤差が負となるひとつ前の分数（初期分数）を求めます。これを B/A とします。ここで真値との誤差は、Δ＝T − B／A＞０ となります。
次に、B/A の分子から−ｍ して、分母から−ｎ した （Ｂ−ｍ）／（Ａ−ｎ） について、Δ（ｍ，ｎ）＝Ｔ−（Ｂ−ｍ）／（Ａ−ｎ）＞０ とすると、Δ（ｍ，ｎ）−Δ（０，０）＝Ｂ／Ａ−（Ｂ−ｍ）／（Ａ−ｎ）＝（ｍＡ−ｎＢ）／Ａ（Ａ−ｎ） です。ここで、Δ（ｍ，ｎ）が初期誤差Δ（０，０）より小さくなるには、Δ（ｍ，ｎ）＜Δ（０，０）で、変形してΔ（ｍ，ｎ）−Δ（ｍ，ｎ）＜０ だから ｍＡ＜ｎＢとなります。
ｍを１，２，３…と変化させて、条件式ｍＡ＜ｎＢに合う ｎ を求めます。そして、このときの（Ｂ−ｍ）／（Ａ−ｎ）を計算して記録しておきます。この例では、初期分数は 141/996 で ｍ＝１，２，３、４…にたいして ｎ＝８，１５，２２，２９…となります。
さらに計算を続けてゆき、15/100 になったら終了します。そうすると、誤差が最小の分数が記録されていることになります。この例では、前述の方法と同じ 16/113 です。この方法を「最小誤差計算法」と呼ぶことにします。
πの近似値の例を示します。 1/7→3.1128571428571428 、14/99→3.1414141414141414 、16/113→3.1415929203539823 、144/1017→3.1415929203539823 、14093/99532→3.1415926536189366 、140914/995207→3.141592635886504 …となります。色の付いた部分がπと一致する部分です。
この小数を分数で近似するアルゴリズムは、私がパソコン雑誌「Ｉ／Ｏ 1991.9」に発表したものです。
点の集合を三角形で分割する
三角形で分割するすぐれた方法としては、点の集合をボロノイ図に分割してそこからドロネー三角形を作成する方法がありますが、難しいので他の方法を自分で考えてみました。結構面倒なのでここでは概要だけを示します。
まず点の集合の外周をＪａｒｖｉｓの凸包アルゴリズムで求めます。内部に残った点の集合に対して再びその凸包を求めます。そうすると外周と内周が求まります。外周の線分（隣り合う２点）から最短距離の内周の１点を求めて三角形を作ります。全ての外周の線分に対してこの操作を行い、三角形を表に登録していきます。できた表をスキャンして同じ三角形を取り除きます。次に内周を外周にしてさらにその内周を求めて、三角形を表に登録する操作を繰り返します。そうすると最終的に内部に点が、なし・１点・２点 の場合になりそれぞれに対応した処理をすると三角形の分割が完了します。この方法を「外周内周による三角形分割」と呼ぶことにします。（2004.3作成）
さらに別の方法で「逐次生成法による三角形分割」を考えましたので、その概要を示します。点集合の任意の１点とそこから最も近いもう１点を選び線分を作ります。この線分から最も近い１点を探して三角形をひとつ作り、辺の情報を登録します。登録した各辺からすべての点への距離を計算してソートします。短い点から始めて、作られる三角形が正しいかチェックしながら新しい三角形を作成します。この操作を繰り返して調べる辺がなくなったら終了します。最後にこの方法では作られる三角形がいびつな形になり易いので、「フリッピング」という手法で最適な形の三角形にします。（2004.7作成）
さらに詳しい説明(PDF) フリッピングの説明(PDF) サンプルプログラム
サンプリング探索法
ソートされたデータの探索法としては、ほとんど二分探索法が用いられますが、他にないか自分で考えてみました。結構面倒なのでここでは概要だけを示します。
昇順にソートされたデータを対象として、ｘ軸にデータの順番（ｘ＝１〜ｎ）をとり、ｙ軸にデータの値をとった平面を考えます。このデータの集合（ｘｉ、ｙｉ）（ｉ＝１〜ｎ）の中から適当にｍ個をサンプリングして、最小二乗法により３次近似式 ｙ＝ａ ｘ^3＋ｂ ｘ^2＋ｃ ｘ＋ｄ の係数を求めます。
ここで探索したいデータ ｙ0 を３次近似式に代入して根を求めます。この根が求めるデータの順番になりデータが発見されたことになります。ただし近似式なので、根も近似解となり根の周辺をさらに探索して真の解を求めます。データの分布状況により近似解の精度が変わってきます。おおむねランダムや正規分布のデータでは精度は比較的良くなります。また、サンプリングの個数を増やせば精度は良くなります。（2005.4作成)
三角形分割によるコンター図（等高線図）
三角形に分割されたされた領域にコンター図を描く方法を自分で考えてみました。ここでは概要だけを示します。
三角形の頂点 A, B, C にデータのレベル a<b<c が与えられているとき、求めたいレベル x が a<x<b とすると当然 a<x<c となるから、x は辺 AB および辺 AC 上に交点があり２点を結ぶ辺が引ける。このように三角形に分割した場合求めたいコンターの交点は、三角形の辺上に２点あるか全くないかのどちらかなのでコンターの線分が一意に定まる。
線形予測ソート
新しいソート（外部ソート）の方法を自分で考えてみました。ここでは概要だけを示します。
ランダムなデータがソートされた後、その値がほぼ直線上に並ぶと仮定します。昇順にソートされたデータについて、ｘ軸にデータの順番（ｘ＝1〜n ）をとり、ｙ軸にデータの値をとった平面を考えます。このデータの集合（ｘｉ、ｙｉ）（ｉ＝1〜n ）が一次式ｙ＝ａ ｘ ＋ｂ を満たすとします。ここで元のランダムなデータ ｙ０を １次式に代入して根を求めます。この根が求めるデータの順番になります。このソートでは１次式で計算して求めた値を記録しておくために別の領域を必要とします。いわゆる外部ソートと呼ばれるものです。新しく記録したデータの集合は空きやリンクがあるので、これを修正して元の集合に戻します。（2007.6作成)
データに順位をつける
データに順位をつける方法を自分で考えてみました。ただしデータの配列はそのままです。ここでは概要だけを示します。
データ（n個）をサーチして最大値と最小値を求めてそれぞれ 1 と n という順位をつけます。ふたたびデータをサーチして最大値と最小値を求めて今度は 2 と n-1 という順位をつけます。 これをデータがなくなるまで繰り返します。 （2017.7作成)
データ（n個）から2個ずつ取り出し2つのグループどうしで順位をつけなおし1つのグループにする。これをデータがなくなるまで繰り返します。（2016.8作成)
また順位の情報に基づいてデータを並び替えればソートすることができます。
鉄道の路線間どうしの最短接続数を求める
鉄道の路線データ（駅のつながりのデータ）があるときにすべての路線間の最短接続数を求める。ここでは概要だけを示します。
隣り合う路線の距離を１と定義する。路線どうしの距離の対応表を作成する。あらかじめ距離１の路線がわかっているとする。路線対応表を用いて、距離１のすべての路線とさらに距離１でつながる路線の距離を２に書き換える。次に距離２の路線と距離１でつながる路線の距離を３に書き換える。これを繰り返して対応表がすべてうまれば終了する。 （2017.7作成)
鉄道の路線間どうしの最短接続数を求めるアルゴリズムを用いたいて駅間の路線を探索する。
最長接続数を求める
路線対応表を用いて、距離１でつながる路線を調べる。路線を見つけたらさらに距離１でつながる路線を探す。これを繰り返して対応表がなくなれば終了する。 （2018.4作成)