Metalシェーダによるマルチエージェントシミュレーション

Metalとは

MetalはiOS8で登場したローレベルなグラフィックAPIで、画像処理や並列演算を得意とします。
このMetalによる並列演算を利用してマルチエージェントシミュレーションを試してみたいと思います。

マルチエージェントシステムとは

マルチエージェントシミュレーションとは、複数のエージェント（人や生物など）で集団を形成し、行動ルールのもとに相互作用をシミュレーションするシステムです。

例えば以下の図のように個人（エージェント）を見ると「進む」というルールに従って行動していますが、集団でみると、近くの人に「追従」するというルールで行動しています。
この相互作用が積み重なると、みんなが同じ方向に進むという結果になります。

このような相互作用は「鳥の群れ」「道路の渋滞」「流行」など身近にも存在しています。
マルチエージェントシミュレーションは、これらの現象についてその仕組みを解析するのに用いられます。

今回は "Dynamic Models of Segregation"（住み分けの動学的モデル）のシミュレーションを行いました。

たとえ白人と黒人が隣同士で暮らすことに抵抗がなくとも、いつの間にか白人が多く居住する地域と黒人が多く居住する地域に分かれてしまう理由を考察するためのモデルで、1971年にThomas Crombie Schellingによって発表されました。

ルールは以下の通りです。

2種類の人種が存在していて最初はランダムに配置されます。
各人は自分の周囲にいる人のうち同じ人種が50%を満たしていれば満足してその場にとどまり、50%に満たない場合はその場を移動します。

これを繰り返すことでどのような結果が得られるかをシミュレーションします。

実装

それではやっとコードの話に進みます。
先にサンプルプロジェクトのリンクを貼っておきます。サンプル
ところどころ省略されているので分かりづらかったらサンプルを動かしてみてください。

今回のマルチエージェントシミュレーションを実行するための流れは以下のようになります。

1. 各エージェントをオブジェクトと UIViewをエージェントの数分用意。
2. 生成したオブジェクトの配列をシェーダに渡し、ルールに従って座標を計算。
3. UIView に計算した座標を反映。
4. 以降繰り返し。

準備

Metalで並列演算を動かすための初期化処理は5ステップです。
以下のものを作成しています。

1. エージェントを表すオブジェクト
2. MTLDevice
3. Metalシェーダ
4. MTLComputePipelineState
5. MTLCommandQueue

1. エージェントを表すオブジェクト

以下のようなstructで表現することとします。
エージェントの数だけ生成しておきます。

struct Object {
    var group: UInt32 = 0
    var x: Float = 0
    var y: Float = 0
    var angle: Float = 0
}

angleはエージェントが移動する際の方向として使用します。

2. MTLDevice

まずは当たり前ですがMetalフレームワークをimportします。

import Metal

MTLDevice は MTLCreateSystemDefaultDevice によって取得できます。

var device: MTLDevice = MTLCreateSystemDefaultDevice()!

3. Metalシェーダ

Metalシェーダは Metal Shading Language という、C++をベースとする独自言語を使用します。

#include <metal_stdlib>
using namespace metal;

struct Object {
    uint group;
    float x;
    float y;
    float angle;
};

float getDistance(float2 pos1, float2 pos2) {
    float2 d = pos1 - pos2;
    return sqrt(d.x*d.x+d.y*d.y);
}

kernel void simulate(const device Object *in [[ buffer(0) ]],
                     const device uint &count [[ buffer(1) ]],
                     const device float &width [[ buffer(2) ]],
                     const device float &height [[ buffer(3) ]],
                     device Object *out [[ buffer(4) ]],
                     uint id [[ thread_position_in_grid ]]) {
    Object object = in[id];
    out[id] = object;

    float2 objectPos = float2(object.x,object.y);

    int total = 0;
    int group = 0;
    for (uint i=0; i<count; i++){
        if (i == id){
            continue;
        };

        Object other = in[i];

        float2 otherPos = float2(other.x,other.y);
        float dist = getDistance(objectPos,otherPos);

        int near = step(dist,20);
        total += near;
        group += (object.group == other.group) * near;
    }
    
    int stay = step(1.0,total)*step(total*0.5, group);
    if (stay) {
        return;
    }
    
    float translationX = sin(object.angle);
    float translationY = cos(object.angle);

    float x = object.x + translationX;
    x -= (x > width) ? width : 0;
    x += (x < 0) ? width : 0;

    float y = object.y + translationY;
    y -= (y > height) ? height : 0;
    y += (y < 0) ? height : 0;

    out[id].x = x;
    out[id].y = y;
}

simulate 関数の引数は以下のとおりです。
第一引数: Object in エージェントの配列
第二引数: uint count エージェントの数
第二引数: float width 画面の幅
第二引数: float height 画面の高さ
第二引数: Object out 結果の出力先
第二引数: uint id エージェントのインデックス
thread_position_in_grid でエージェントのインデックスが取得できるので、各エージェントは in[id] で取得できます。
※ thread_position_in_grid についてはまたあとで説明が出てきます。
あとは先述のルールに従って計算処理を書いていき、 out に計算結果を反映させます。

4. MTLComputePipelineState

Metalシェーダを作成したので、これをコンパイルする必要があります。
MTLComputePipelineState プロトコルは、コンパイルされたプログラムへの参照をエンコードするために使用される軽量オブジェクトのインターフェイスを定義します。
MTLComputePipelineState の生成は負荷が高いのでメンバ変数として定義し初期化時に生成します。

var pipelineState: MTLComputePipelineState!

let defaultLibrary = device.makeDefaultLibrary()
let function = defaultLibrary.makeFunction(name: "simulate")!
pipelineState = try! device.makeComputePipelineState(function: function)

device.makeDefaultLibrary() で取得した MTLLibrary オブジェクトを使ってコンパイルしたシェーダへアクセスすることができます。
シェーダの関数名を指定し MTLFunction を作成したら、 device.makeComputePipelineState(function:) で MTLFunction をコンパイルしたコードに変換します。

5. MTLCommandQueue

最後に MTLCommandQueue を生成します。
MTLCommandQueue プロトコルは、GPUが実行するコマンド・バッファの順序付きリストをキューに入れることができるオブジェクトのインタフェースを定義します。

これもメンバ変数として定義しておきます。

var commandQueue: MTLCommandQueue!

commandQueue = device.makeCommandQueue()

演算処理の実行

並列演算を動かすための準備ができたので次は実行するための処理を書いていきます。
実行用に func execute() -> [Object] というメソッドを用意して以下の処理を書いていきます。

1. MTLCommandBuffer の作成
2. 実行コマンドをエンコード
3. MTLCommandBuffer をコミット
4. 計算結果を返す

1. MTLCommandBuffer の作成

実行したいコマンドのリストとなります。

let commandBuffer = commandQueue?.makeCommandBuffer()

2. 実行コマンドをエンコード

let commandBuffer = commandQueue.makeCommandBuffer()

let encoder = commandBuffer?.makeComputeCommandEncoder()
encoder?.setComputePipelineState(pipelineState)

let buffer = device?.makeBuffer(bytes: objects, length: objects.byteLength, options: [])
var count = UInt32(objects.count)
let countBuffer = device?.makeBuffer(bytes: &count, length: MemoryLayout.size(ofValue: count), options: [])
var width = Float(UIScreen.main.bounds.size.width)
var height = Float(UIScreen.main.bounds.size.height)
let widthBuffer = device?.makeBuffer(bytes: &width, length: MemoryLayout.size(ofValue: width), options: [])
letheightBuffer = device?.makeBuffer(bytes: &height, length: MemoryLayout.size(ofValue: height), options: [])
let outBuffer = device?.makeBuffer(bytes: objects, length: objects.byteLength, options: [])!

encoder?.setBuffer(buffer, offset: 0, index: 0)
encoder?.setBuffer(countBuffer, offset: 0, index: 1)
encoder?.setBuffer(widthBuffer, offset: 0, index: 2)
encoder?.setBuffer(heightBuffer, offset: 0, index: 3)
encoder?.setBuffer(outBuffer, offset: 0, index: 4)

let groupsize = MTLSize(width: 64, height: 1, depth: 1)
let numgroups = MTLSize(width: (objects.count + groupsize.width - 1) / groupsize.width, height: 1, depth: 1)

encoder?.dispatchThreadgroups(numgroups, threadsPerThreadgroup: groupsize)
encoder?.endEncoding()

コマンドをエンコードするためのエンコーダ（MTLComputeCommandEncoder）を作成します。
MTLComputeCommandEncoder オブジェクトを作成したら、このオブジェクトを使用して、並列演算処理コマンドをエンコードします。
先程作成した、実行される関数を含む MTLComputePipelineState オブジェクトをセットします。
入力データと出力先を保持するリソースを指定します。（指定するリソースは simulate 関数の引数で定義した5つです。）
dispatchThreadgroups(_: threadsPerThreadgroup:) メソッドを呼び出して、グリッド用の指定された数のスレッドグループとスレッドグループあたりのスレッド数で計算関数をエンコードします。

Metalシェーダの項で thread_position_in_grid が登場しましたが、dispatchThreadgroups(_: threadsPerThreadgroup:) で指定するサイズによってグリッドの範囲が決まり、thread_position_in_grid でそのグリッド内のインデックスが取得できるということだそうです。

endEncoding()を呼び出してエンコードを終了します。

3. MTLCommandBuffer をコミット

コマンドを実行します。

commandBuffer?.commit()
commandBuffer?.waitUntilCompleted()

4. 計算結果を返す

3で実行すると演算処理が走り、outBufferから結果を取得し返します。

let data = Data(bytesNoCopy: outBuffer.contents(), count: objects.byteLength, deallocator: .none)
objects = data.withUnsafeBytes {
    [Object](UnsafeBufferPointer<Object>(start: $0, count: data.count/MemoryLayout<Object>.size))
}
return objects

あとは UIView 計算した座標を反映させれば完了です。

実行結果

最初はランダムに配置されていますが、シミュレートを重ねると右のように同じ人種同士で固まるという結果が得られました。

→

周囲の50%を同人種で満たしていればOKということは残り約半数は異人種でもいいということになるので結構寛容的に思われますが、それでもこのようにはっきり同人種同士で固まることになるのですね。
試しに許容の水準を30%、10%に下げてみてもある程度分居される結果になりました。

30%	10%

|

シェリングは現象を隣人に対する寛容性あるいは我慢の強弱に起因するものであると仮定していますが、確かに個々が寛容的であっても分居は発生し、そのレベルは許容の水準によるものが影響しているようにみえます。

ネイティブでやった場合の比較

最後に、Metalを使うとどれくらい良いのか、ネイティブとの比較を載せておきます。

実行速度

60回シミュレーションを実行した平均はそれぞれ以下の通りです

エージェント数	Metal	ネイティブ
100	0.002582641443	0.007742651304
500	0.003791083892	0.03941152493
1000	0.004310701291	0.1502834062
1500	0.005150139332	0.3368620038
2000	0.00574699839	0.5972672562

CPU

エージェント数2000のCPUです

Metal	ネイティブ

どちらも圧倒的にMetalのほうがパフォーマンスが良いですね。

まとめ

マルチエージェントシミュレーションをMetalで試してみましたが、Metalシェーダ（Metal Shading Language）を書くというハードルはありますが、Metalフレームワークを使うと割と簡単に動かせました。
マルチエージェントシミュレーションについては今回は簡単なモデルで試しましたが、他にもいろいろモデルがるのでシェーダの練習がてらシミュレーションしてみるのも良いかもしれません。
自分でも身近にある現象と理由を考察してシミュレーションしてみると面白そうですね。（私にはハードル高いけど...）

参考

https://qiita.com/yuky_az/items/3f3de1034e0600be490d
http://mas.kke.co.jp/modules/tinyd4/index.php?id=13