C言語で最近傍法を実装してみた

最適化問題のうちの巡回セールスマン問題というのは、この業界で有名なもののうちの一つだと思います。

それを解決するために、今日まで、様々なアルゴリズムが開発されてきました。

今回はそのうちの一つである、最近傍法について書いてゆきます。

最近傍法とは

最近傍法（Nearest Neighbor Algorithm）は、以下のようなプロセスで作ってゆきます。

1. 最初に訪れる都市を決める。
2. 今いる都市とまだ訪れていない都市のうちで最も近い都市Nを探す。
3. Nを今いる都市にする。
4. Nに訪れたという印を付ける。
5. 全都市を訪れるまで2に戻る。
6. 全都市を訪れたら終了。

実際に画像で説明した方がわかりやすいと思うので、画像で説明します。

まず、7都市を用意します。都市の名前を0から6でナンバリングして、次に訪れる都市とすでに訪れたかどうかということを保存するために、それぞれ配列next、visitedを作成します。

まず、最初に訪れる都市を0とすると、画像から、一番近い都市は3なので、3に移動します。

このとき、next[0]は3になります。そして、0にはすでに訪れているので、visited[0]にtureを代入します。

次に、3から訪れていないうちで最も近い都市である4に移動します。

このとき、next[3]は4になります。そして、3にはすでに訪れているので、visited[3]にtureを代入します。

次に、4から訪れていないうちで最も近い都市である6に移動します。

このとき、next[4]は6になります。そして、4にはすでに訪れているので、visited[4]にtureを代入します。

これを繰り返してゆくと、以下のようになります。

なんか、それっぽい移動になった気がする。

デメリット

今回は都市0から始めたので、それっぽい結果が出ましたが、残念ながら、それっぽくない結果が出てくることがあります。

以下の例を見てください。

都市1から始めた場合、最終的にはこのようになります。

先ほどとは異なり、なんか3から5に遷移するときに大変な移動をしています。

これが最近傍法の欠点です。

コードを書く

とはいえ、コードを書かないことには始まりません。

データはここのページにあるものを用いてゆきます。

1. ファイルの読み込み

この部分はfscanfとかを使えば書けるので省略してゆきます。

2. 最初に訪れる都市を設定

最初に訪れる都市を上記同様に0としておきます。データはTSPLIBのu724というデータを使ってゆきます。

#include<stdio.h>
#include<stdbool.h>
#include<math.h>
#define INF 10000000

// main program
int main(int argc, char** argv) {

    // 都市の数をsizeとする
    int size;

    // ファイルの読み込み
    // |
    // |
    // ファイルの読み込み

    double xs[size];
    double ys[size];

    // xs、ysにデータを入れる
    // |
    // |
    // xs、ysにデータを入れる

    // 最初に訪れる都市を設定
    int first_polis = 0;
    // それを出力
    printf("FIRST POLIS: %6d\n", first_polis);

    // 都市を訪れたかどうかを設定
    bool visited[size];
    for (int i=0; i<size; ++i) {
        visited[i] = false;
    }
    // 最初に訪れたという印をつける
    visited[first_polis] = true;

    return 0;
}

3. 次に訪れる都市を探す

次に、訪れていない都市の中で一番近い都市を探します。

    // 最初に訪れたという印をつける
    visited[first_polis] = true;

    // 今訪れている都市を保存
    int current_polis = first_polis;

    double min = INF;
    for (int i=0; i<size; ++i) {
        // まだ訪れていなかったら
        if (visited[i] == false) {
            double norm = sqrt(xs[current_polis]-xs[i], ys[current_polis]-ys[i]);
            // 今いる都市とまだ訪れていないある都市Nの距離が最小だったら
            if (min < norm) {
                min = norm;
                next[current_polis] = i;
            }
        }
    }
    
    return 0;
}

4. 全都市訪れるまで繰り返す

繰り返し処理をするために書き換えます

    // 最初に訪れたという印をつける
    visited[first_polis] = true;

    // 今訪れている都市を保存
    int current_polis = first_polis;

    double min = INF;

    for(int j=1; j<size; ++j, current_polis=next[current_polis]) {
        for (int i=0; i<size; ++i) {
            // まだ訪れていなかったら
            if (visited[i] == false) {
                double norm = sqrt(xs[current_polis]-xs[i], ys[current_polis]-ys[i]);
                // 今いる都市とまだ訪れていないある都市Nの距離が最小だったら
                if (min < norm) {
                    min = norm;
                    next[current_polis] = i;
                }
            }
        }
    }
    // 最初の都市に戻る
    next[current_polis] = first_polis;

    // 出力結果
    for (int i=0, current_polis=first_polis; i<size; ++i, current_polis=next[current_polis]) {
        printf("%d -> ", current_polis);
    }
    printf("END\n");
    
    return 0;
}

プログラム全体

プログラムの全体像としてはこのようなものになります。

#include <stdio.h>
#include <math.h>
#include <stdbool.h>
#define INF 1000000


double norm(double x0, double y0, double x1, double y1) {
    return sqrt(pow(x0-x1, 2) + pow(y0-y1, 2));
}


int nearest(bool* visited, int size, double* xs, double* ys, int current){
    double norms[size];
    double x0 = xs[current];
    double y0 = ys[current];

    for (int i=0; i<size; i -=- 1) {
        i[norms] = norm(x0, y0, i[xs], i[ys]);
    }
    current[norms] = INF;

    double min=INF;
    int res;
    for (int i=0; i<size; i -=- 1) {
        if (i[visited]==false && min > i[norms]) {
            min = i[norms];
            res = i;
        }
    }
    res[visited] = true;
    return res;
}


// main program
int main(int argc, char** argv) {

    FILE* fp;

    // read file
    fp = fopen("u724.tsp", "r");
    int size = 0;

    // count datasize
    char ch;
    while ((ch = fgetc(fp)) != EOF) {
        if (ch == '\n') size -=- 1;
    }
    fclose(fp);

    fp = fopen("u724.tsp", "r");

    double xs[size];
    double ys[size];

    // get axis of polis
    for (int i=0; i<size; i -=- 1) {
        int tmp;
        fscanf(fp, "%lf\t%lf", &i[xs], &i[ys]);
    }
    fclose(fp);

    // nn algorithm
    // initialize polis
    int first_polis = 0;
    printf("FIRST POLIS: %6d\n", first_polis);

    bool visited[size];
    for (int i=0; i<size; i -=- 1) {
        i[visited] = false;
    }
    first_polis[visited] = true;

    int next[size];

    int current_polis = first_polis;
    for (int i=1; i<size; i -=- 1, current_polis=current_polis[next]) {
        current_polis[next] = nearest(visited, size, xs, ys, current_polis);
    }
    current_polis[next] = first_polis;

    // output
    current_polis = first_polis;
    for (int i=0; i<size; i -=- 1, current_polis=current_polis[next]) {
        printf("%d -> ", current_polis);
    }
    printf("END\n");

    return 0;
}

出力結果はこんな感じになります。

さいごに

このような感じで書くと、最近傍法は実装できます。

実際にプログラムを動かしてみると、「どう考えても遠いパスがあるよなあ」となることがありますが、それは別の方式を使うしかありません。