ぱたへね

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

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の量子コンピュータ対応に期待したい。