ぱたへね

はてなダイアリーはrustの色分けができないのでこっちに来た

組合せ最適化アルゴリズムの最新手法―基礎から工学応用まで

組合せ最適化アルゴリズムの最新手法―基礎から工学応用までを読みました。

80年代後半、90年代の論文がよく出てきており、20年前の最先端をまとめた本です。題材に半導体関連のネタが多く、興味を持って最後まで読めました。シミュレーティド・アニーリング手法(SA)、遺伝的アルゴリズム(GAs)、タブー・サーチ手法(TS)、シミュレーティド・エボリューション手法(SimE)、確率的進化手法(StocE)について分かりやすく書いてありました。 ただ、この本を読んで実際の問題を解けるかどうかといわれると難しい気がします。

組み合わせ最適化について

この本にのっているアルゴリズムは、基本的にはすべてこのような手順になっています。

  • 上手い初期値を設定する
  • 上手く少しだけ良い方向へ変化する
  • 上手く局所解を抜けて、大域解へ向かう。

局所解を抜けるために、どこかで乱数を使います。各アルゴリズムの違いはどこを重視するかです。 遺伝的アルゴリズムであれば、初期値を大量に作り生き残りを作っていくことで、上手い初期値としています。 タブーサーチ方では上手く良い方向へへの変化を求めるときに、過去の失敗をタブーとして参考にします。

直行している概念も多く、最後の章ではそれぞれの組み合わせの考察もありました。

実際はこのアルゴリズムにかける前に、上手いモデルの定義があって、それはそれで見るからに面倒でした。上手く行った事例が結構載っているのですが、それ以外は上手く行かなかったんだろうなという感じです。

並列コンピューティングやSIMDへの期待が高いのも興味深く読めました。

当時のニューラルネットワークの評価

面白いことに半導体の配置問題をニューラルネットワークを使って解くという89年の取り組みが紹介されています。30年前のニューラルネットワークの評価として貴重な資料です。

本の中では、このようにまとめられています。

ホップフィールド型ネットワークを配置問題に適用したユーの結果は、有望なものではなかった。いくつかの問題点は、シミュレーション時間が長いこと、解の質が低いこと、および、ネットワークのパラメータに対して、解の感度が高いことである。現段階では結論として、ニューラルネットワークを困難な最適化問題にてきようできるかどうかについては、さらに多くの研究が必要であると、述べるにとどめる。

RustでSimulated Annealing

本棚を見たら組合せ最適化アルゴリズムの最新手法―基礎から工学応用までという本を見つけました。Simulated Annealingのアルゴリズムが載っていたので、せっかくだからRustで実装してみました。

問題

f:id:natsutan:20181117082906p:plain

この図の回路を2つのエリアにわける方法を考えます。簡単のために、10個あるノード(ゲート)を5個、5個に分けることにします。ノードには1~10の通し番号が入っていて、それぞれのネットには分割時のコストが設定されています。

分ける時のコストがこのように定義されています。

net node cost
N1 {1, 2, 4, 5} 1
N2 {2, 3, 5} 1
N3 {3, 6, 10, 4} 2
N4 {4, 8, 3, 7} 1
N5 {5, 7, 1, 6} 3
N6 {6, 4, 7, 2} 3
N7 {7, 9, 5} 2
N8 {8, 2} 3
N9 {9, 10, 5} 2
N10 {10, 5} 4

ノードを接続しているネットにもN1~N10の通し番号が打ってあります。例えばn1は、1, 2, 4, 5のノードが接続されています。 ノードを2つに分割した結果、1, 2, 4, 5が同じエリアに存在する場合はコストは0、2つのエリアに分割された場合は1のコストが発生します。分割した結果のコストの和を最小化することが目的です。

このようにオレンジと緑に分割したとします。

f:id:natsutan:20181117082937p:plain

2つのエリアをまたぐネットは赤線にしています。赤線のネットに設定されているコストを全て足すと10になり、この分割をしたときのコストになります。

こういう問題を「最小集合2分割問題」と呼び、意外にもNP困難です。そこで、Simulated Annealingの出番です。

Simulated Annealingアルゴリズム

そんなに難しくない。初期値はランダムに2つに分割します。

  • 現在のコストを計算する。
  • ランダムに1つノードを入れ替える。
  • コストが下がった場合は、一つ上に戻りノードを入れ替える。
  • コストが上がる場合は、以下の確率でその分割を採用する。二つ上に戻り、ノードを入れ替える。

コストが悪くなってもその分割を採用する確率はこのように決められています。

P = e^{-\Delta C /T}

ΔCが新しいコストと現在のコストの差、Tは現在の温度です。時間と共に温度が下がっていきます。温度という概念がでてきましたが、基本的にはループの回数です。

やってみた結果

本に載っているとおりにやってみました。

f:id:natsutan:20181117090309p:plain

この回路の場合分割後のコストが10が最小のようなのですが、何回かやると10になっていることが分かります。ちなみにここで求めた10が最小かどうかはSimulated Annealingでは分かりません。

ソースコード

extern crate rand;

use rand::thread_rng;
use rand::Rng;
use rand::seq::SliceRandom;


struct Net {
    nodes:Vec<u32>,
    weight:u32,
}

type Block =  [u32; 10];

//ネットリストの定義
fn make_netlist() -> Vec<Net> {
    let mut netlist :Vec<Net> = Vec::new();

    let n1 = Net{nodes:vec![1, 2, 4, 5], weight: 1};
    let n2 = Net{nodes:vec![2, 3, 5], weight: 1};
    let n3 = Net{nodes:vec![3, 6, 10, 4], weight: 2};
    let n4 = Net{nodes:vec![4, 8, 3, 7], weight: 1};
    let n5 = Net{nodes:vec![5, 7, 1, 6], weight: 3};
    let n6 = Net{nodes:vec![6, 4, 7, 2], weight: 3};
    let n7 = Net{nodes:vec![7, 9, 5], weight: 2};
    let n8 = Net{nodes:vec![8, 2], weight: 3};
    let n9 = Net{nodes:vec![9, 10,  5], weight: 2};
    let n10 = Net{nodes:vec![10, 5], weight: 4};

    netlist.push(n1);
    netlist.push(n2);
    netlist.push(n3);
    netlist.push(n4);
    netlist.push(n5);
    netlist.push(n6);
    netlist.push(n7);
    netlist.push(n8);
    netlist.push(n9);
    netlist.push(n10);

    netlist
}

//ランダムにブロックを振り分ける
fn initial_block() -> Block {
    let mut rng = thread_rng();
    let mut new_block: Block = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10];
    new_block.shuffle(&mut rng);
    new_block
}

//引数のnetが、ABどちらかのブロックに入っていればtrue
//AB両方のブロックにまたがっていればfalse
fn same_area(block:&Block, net:&Net) -> bool {
    let block_a = &block[0..5];
    let block_b = &block[5..10];
    //全てtrue、もしくは全てfalseの時trueを返す
    let contain_as: Vec<bool> = block_a.into_iter().map(|x| net.nodes.contains(x)).collect();
    let contain_bs: Vec<bool> = block_b.into_iter().map(|x| net.nodes.contains(x)).collect();

    let all_b = !contain_as.into_iter().fold(false, |all_a, x| all_a | x);
    let all_a = !contain_bs.into_iter().fold(false, |all_b, x| all_b | x);

    all_a || all_b
}

//分割後のcost計算
fn cost(block:&Block, netlist:&Vec<Net>) -> u32 {
    // blockの前半が、分割したエリアに全て入ればそのまま、そうで無ければそのネットに対応したweightを加算する。
    netlist.into_iter().fold(0,|cost, net| if same_area(&block, net) {cost} else {cost + net.weight} )
}

//ランダムに2つのノードを入れ替える
fn swap_node(block:&Block) -> Block {
    let mut new_block = block.clone();
    let mut rng = thread_rng();

    let choice_a = rng.gen_range(0, 6);
    let choice_b = rng.gen_range(6, 10);

    new_block.swap(choice_a, choice_b);
    new_block
}

#[test]
fn cost_test(){
    let netlist = make_netlist();
    let mut block:[u32; 10] = [2, 4, 6, 7, 8, 1, 3, 5, 9, 10];
    let cur_cost = cost(&block, &netlist);
    assert_eq!(cur_cost, 10)
}

//メトロポリス基準
// random < exp(- d_cost / t)
fn random(delta_cost:i32, t:f32) -> bool {
    let mut rng = thread_rng();
    let rand_value = rng.gen_range(0.0, 1.0);
    let d_cost:f64 = delta_cost as f64;

    rand_value < std::f64::consts::E.powf(- d_cost / t as f64)
}

//シミュレーテッドアニーリングの実行
fn simulated_annealing(mut t:f32, alpha:f32, beta:f32, max_time:f32, mut m:f32 , netlist:&Vec<Net> ) -> Block
{
    let mut cur_s :Block = initial_block();
    let mut time:f32 = 0.0;
    let mut best_s: Block = cur_s.clone();

    while {
        let result = metropolis(cur_s, best_s, t, m, netlist);
        cur_s = result[0];
        best_s = result[1];

        time = time + m * beta;
        t = alpha * t;
        m = m * beta;
        time < max_time
    } {}
    best_s
}

fn metropolis(mut cur_s: Block, mut best_s: Block, t:f32, m:f32, netlist:&Vec<Net>) ->[Block;2]
{
    let cur_cost = cost(&cur_s, netlist);
    let mut best_cost = cost(&best_s, netlist);


    for _x in 0 .. m as u32 {
        //ノードを一つランダムに交換
        let new_s = swap_node(&cur_s);
        let new_cost = cost(&new_s, netlist);
        let delta_cost = new_cost as i32 - cur_cost as i32;

        println!("new_cost = {} best_cost = {} block_a = {:?} block_b = {:?}", new_cost, best_cost,  &new_s[0..5], &new_s[5..10]);

        if delta_cost < 0 {
            cur_s = new_s;
            if new_cost < best_cost {
                best_s = new_s;
                best_cost = cost(&best_s, netlist);
            }
        } else {
            if random(delta_cost, t) {
                cur_s = new_s;
            }
        }
    }
    [cur_s, best_s]
}


fn main() {
    let netlist = make_netlist();

    let t0 = 10.0;
    let alpha = 0.9;
    let beta = 1.0;
    let m = 10.0;

    let best_block = simulated_annealing(t0, alpha, beta, 100.0, m, &netlist);
    let best_cost = cost(&best_block, &netlist);

    //最終結果
    println!("cost = {} block_a = {:?} block_b = {:?} ", best_cost,  &best_block[0..5], &best_block[5..10]);[f:id:natsutan:20181117082906p:plain]

}

まとめ

Vivadoの量子コンピュータ対応に期待したい。

rustでマクロ

rust standard library cookbookから。

マクロを使って、簡単な可変長引数を実現する例。

macro_rules! multiply {
    ( $last:expr ) => { $last };

    ( $head:expr, $($tail:expr), + ) => {
        $head * multiply!($($tail), +)
    };

}

fn main()
{
    let val1 = multiply!(2, 4, 8);
    println!("2 * 4 * 8 = {}", val1);

    let val2 = multiply!(2, 4, 8, 16);
    println!("2 * 4 * 8 * 16 = {}", val2);
}

マクロの中で再帰するのか。すごいな。

Rustでグラフ

ほとんどライブラリの使い方レベルだけど、なんとか動いた。

extern crate petgraph;

use std::fs;
use std::io::Write;

// https://docs.rs/petgraph/0.4.13/petgraph/
use petgraph::Graph;
use petgraph::algo;
use petgraph::dot::Dot;

fn main() {

    let mut graph = Graph::new_undirected();
    let na = graph.add_node("a");
    let nb = graph.add_node("b");
    let nc = graph.add_node("c");
    let nd = graph.add_node("d");
    let ne = graph.add_node("e");
    let nf = graph.add_node("f");

    graph.add_edge(na, nb, 20);
    graph.add_edge(nb, nc, 20);
    graph.add_edge(nc, nd, 10);
    graph.add_edge(na, ne, 10);
    graph.add_edge(nb, ne, 10);
    graph.add_edge(nb, nf, 30);
    graph.add_edge(nd, nf, 10);
    graph.add_edge(ne, nf, 20);

    let mut f = fs::File::create("graph.dot").unwrap();
    let dot_str = format!("{:?}", Dot::with_config(&graph, &[]));
    f.write(dot_str.as_bytes()).unwrap();

    let path = algo::dijkstra(
        &graph,
        na,
         Some(nd),
        |e| *e.weight()
    );
    println!("result: {:?}", path);
}

実行結果

f:id:natsutan:20181019131430p:plain

なんとなく動いているっぽい

result: {NodeIndex(1): 20, NodeIndex(5): 30, NodeIndex(0): 0, NodeIndex(3): 40, NodeIndex(4): 10, NodeIndex(2): 40}

久しぶりにKerasを動かしたときのメモ

import keras しただけでこういうエラーがでる。

ImportError: libcublas.so.9.0: cannot open shared object file: No such file or directory

libcublasが見つからないだけなので、.bashrcにLD_LIBRARY_PATHを追加する。

export LD_LIBRARY_PATH=$LD_LIBRARY_PATH:/usr/local/cuda-9.0/lib64/

ターミナルから動いてPyCharmで同じエラーが出る場合は、PyCharmのrun configurations(working directoryとかを設定する画面)から、Environmet variables に上のLD_LIBRARY_PATHを追加する。

TVMで使われている最適化の探索手法

TVMとは

ペーパーはこちらからダウンロードできます。 https://arxiv.org/abs/1802.04799

TVMは、ASICやAPGAの様な組込様CPU等、様々なバックエンドでDeep Learningを動かすことを目的としたコンパイラです。 グラフレベルや演算子レベルの様々な最適化を行います。

興味本位でペーパーを読んでみたら、僕の知らない最適化テクニックが使われていたので紹介します。 他にも面白いことは書いてあったのですが、最も気になった2つを紹介します。

TVM stack

f:id:natsutan:20181013084651p:plain

TVMは図のようなスタック構造になっています。上から順番に最適化をしていって、最後にLLVMやCUDAのソースコードを生成します。High level Graph Rewritingの話(Section 3)や、Tensorとスケジューリングの話(Section 4)で最適化した後、それらの組み合わせを考慮するときに、機械学習を使って最適化をしています。

TVMは、このOptimizerに来る前に可能な限りの最適化を考えます。その結果、10億以上の最適化の組み合わせ(configuration)から、実際のワークロードに適した組み合わせの検索をOptimizerが行います。つまり 、Optimizerは実際に最適化をするのではなく、最適化の組み合わせの探索に機械学習を使っています。

ML based Cost Model

最適化の探索には、様々なファクターがあります。 - メモリアクセスパターン - データの再利用 - パイプラインの依存性 - スレッドのパターン - etc 最適化のするに当たっては計算のコストモデルを作る必要があります。 従来のAuto-turningの手法では良い組み合わせを探すのに多くの実験が必要になります。定義済みのコストモデルを使う場合は、バックエンドのハードウェアへの依存性が高くなります。そこでTVMでは、機械学習を使ったアプローチを取っています。

f:id:natsutan:20181013084809p:plain

図のように、ループの構文ツリーをTreeRNNに入れ、特徴量を出しコストを予想しています。

Schedule Exploration

コストモデルが決まった後は、そのコストモデルに従って一番良い組み合わせを探します。理想的には全てのコンフィグレーションを試し上位のいくつかを選択すれば良いのですが、それには強大な探索空間が必要になります。その問題を解決するために、TVMは a parallel simulated annealing algorithmを使います。

まずはランダムに選んだコンフィグレーションから開始し、ランダムウォークをしながら、コストモデルが予測するコストが低い方へ移動していきます。ランダムウォークはコストの低いところを見つけ、移動しなくなるまで続けられます。

感想

コンパイラの最適化技術がここまで進んでいるとは思いませんでした。個別のターゲットに対して、個別の最適化を行う時代は過ぎているだと思います。AIによって世界中の誰でも、どのデバイスでも、最適化職人の技術の恩恵にあずかれると思うと素晴らしいですね。

コンパイラ勉強会が来月に開催されるようで、申し込んではみましたが抽選に通るかどうか難しいそうです。 https://connpass.com/event/103976/

Rustで文字列の連結

Rust Standard Library Cookbookを読み始めました。

Amazon CAPTCHA

この本ではRustの文字列連結として3つの方法が紹介されています。

moveによる連結

一番標準の方法。わかりやすく効率も良いが、連結後はmoveによりhelloは使えなくなる。

fn by_moving() {
    let hello = "hello ".to_string();
    let world = "world";

    let hello_world = hello + world;

    println!("{}", hello_world);
}

helloはStringで、worldが&str。この2つを連結するとStringになるというのもよくわかってない。(なんで to_string()が必要なの?)

cloneによる連結

clone()を使うことで、一時変数が作られ、それをmoveすることで文字列を連結している。連結後もhelloが使えることがポイント。

fn by_cloning() {
    let hello = "hello".to_string();
    let world = "world";

    let hello_world = hello.clone()  + world;
    println!("{} ", hello_world);
}

mutingによる連結

helloをmutableで宣言して、push_strを使う方法。エレガントじゃないらしい。

fn by_mutating() {
    let mut hello = "hello".to_string();
    let world = "world";

    hello.push_str(world);
    println!("{} ", hello);
}

今まで勉強したプログラミング言語の中で一番文字列がめんどくさい。

元のソースはこちらです。 github.com