Rustエンジニア必見：ndarrayで数値計算の効率を上げる方法

By /

Rustエンジニア必見：ndarrayで数値計算の効率を上げる方法

Rustは、その安全性とパフォーマンスによって、システムプログラミングやWeb開発など幅広い分野で注目されています。しかし、PythonやMATLAB、Rといった言語に比べて、大規模な数値計算や科学計算の分野ではまだ普及途上と言えるでしょう。その理由の一つに、標準ライブラリだけでは多次元配列の操作や高度な数値計算機能が不足している点が挙げられます。

しかし、Rustのエコシステムは急速に発展しており、数値計算を効率的に行うための強力なライブラリが登場しています。その代表格が、多次元配列ライブラリであるndarrayです。ndarrayは、PythonのNumPyにインスパイアされており、Rustの強みである安全性とパフォーマンスを活かしながら、柔軟かつ効率的な多次元配列操作を提供します。

この記事では、Rustで数値計算を行う際にndarrayがいかに強力であるか、そしてどのようにすればndarrayを最大限に活用して計算効率を向上させられるかを、詳細な説明と豊富なコード例とともに解説します。

この記事を読むことで、以下のことが学べます。

  • ndarrayの基本的な使い方と概念
  • ndarrayがなぜ高速なのか、その内部構造と設計思想
  • 効率的な数値計算を実現するための具体的なテクニック（ビュー、ブロードキャスト、並列化など）
  • 線形代数や統計計算など、高度な数値計算タスクへのndarrayの応用
  • パフォーマンス測定の方法と、よくある落とし穴

それでは、Rustでの数値計算の世界へ、ndarrayとともに踏み込んでいきましょう。

1. はじめに：Rustにおける数値計算の現状とndarrayの重要性

Rustは、メモリ安全性を保証しつつCやC++に匹敵するパフォーマンスを発揮できるという独自の強みを持っています。この特性は、数値計算、特に大規模なデータセットや計算量の多いアルゴリズムを扱う場合に非常に魅力的です。しかし、標準ライブラリのVec<T>は1次元配列としては優秀ですが、多次元配列の表現や操作、行列演算といった機能は提供していません。PythonのNumPyのような、多次元配列を効率的に扱うための基盤がRustには不可欠でした。

そこで登場したのがndarrayクレートです。ndarrayは、固定および動的な次元を持つ多次元配列を扱うためのライブラリであり、NumPyユーザーにとって非常に馴染みやすいAPIを提供します。これにより、RustでもPythonやMATLABのように直感的に多次元配列を操作し、効率的な数値計算を記述することが可能になりました。

単に機能を提供するだけでなく、ndarrayはRustの所有権システムや借用ルールと調和するように設計されています。これにより、データ競合や無駄なデータコピーを防ぎつつ、安全かつ高性能なコードを書くことができます。また、並列計算ライブラリrayonとの連携や、最適化された外部の線形代数ライブラリ（BLASやLAPACK）との統合により、さらなるパフォーマンス向上を実現しています。

この先では、まずndarrayの基本的な使い方から始め、その設計がどのように効率性に貢献しているのか、そして具体的なテクニックを学ぶことで、Rustにおける数値計算の効率を飛躍的に向上させる方法を探っていきます。

2. ndarrayの基本：多次元配列の生成、操作、形状

ndarrayの最も基本的な構成要素は、Array<T, D>構造体です。これは、要素型Tと次元数D（Ix1、Ix2、Ix3などの型エイリアス、あるいはIxDynによる動的な次元数）を持つ多次元配列を表します。

2.1. ndarrayのインストール

Cargo.tomlに依存関係を追加します。

toml
[dependencies]
ndarray = "0.15" # 最新バージョンを指定してください

2.2. 配列の生成

ndarrayで配列を生成する方法はいくつかあります。

  • マクロ array!: 簡単な配列を手っ取り早く作成する場合に便利です。ネストされたブラケットを使って次元を表現します。

    “`rust
    use ndarray::array;

    // 1次元配列 (Vector)
    let a = array![1, 2, 3, 4, 5];
    println!(“a: {}”, a); // [1, 2, 3, 4, 5]

    // 2次元配列 (Matrix)
    let b = array![[1.0, 2.0, 3.0], [4.0, 5.0, 6.0]];
    println!(“b:\n{}”, b);
    / 出力例:
    b:
    [[1.0, 2.0, 3.0],
    [4.0, 5.0, 6.0]]
    /

    // 3次元配列
    let c = array![[[1, 2], [3, 4]], [[5, 6], [7, 8]]];
    println!(“c:\n{}”, c);
    “`

  • Array::zeros(), Array::ones(), Array::from_elem(): 指定した形状で、全ての要素が0、1、または特定の値の配列を作成します。形状はタプルで指定します。

    “`rust
    use ndarray::{Array, Ix2};

    // 3×4のゼロ行列
    let zeros = Array::::zeros((3, 4));
    println!(“zeros:\n{}”, zeros);
    / 出力例:
    zeros:
    [[0.0, 0.0, 0.0, 0.0],
    [0.0, 0.0, 0.0, 0.0],
    [0.0, 0.0, 0.0, 0.0]]
    /

    // 2×2のイチ行列
    let ones = Array::::ones((2, 2)); // 次元は推論させることも可能
    println!(“ones:\n{}”, ones);

    // 4×1の全てが7の行列
    let sevens = Array::::from_elem((4, 1), 7);
    println!(“sevens:\n{}”, sevens);
    “`

  • Array::from_shape_vec(): 1次元のVec<T>を指定した形状の多次元配列に変換します。要素数と形状の要素数の積が一致している必要があります。

    “`rust
    use ndarray::{Array, Ix3};

    let data = vec![1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12];
    let shape = (2, 3, 2); // 2層, 3行, 2列
    let arr = Array::from_shape_vec(shape, data).expect(“形状と要素数が一致しません”);
    println!(“arr:\n{}”, arr);
    /* 出力例:
    arr:
    [[[1, 2],
    [3, 4],
    [5, 6]],

    [[7, 8],
    [9, 10],
    [11, 12]]]
    */
    “`

  • Array::from_shape_fn(): 形状と、各要素の座標を受け取り値を計算するクロージャを指定して配列を作成します。

    “`rust
    use ndarray::Array;

    // 3×3の単位行列
    let identity = Array::from_shape_fn((3, 3), |(i, j)| if i == j { 1.0 } else { 0.0 });
    println!(“identity:\n{}”, identity);
    “`

2.3. 要素へのアクセス

要素へのアクセスは、インデックスを使って行います。多次元配列の場合、タプルを使って各次元のインデックスを指定します。

“`rust
use ndarray::array;

let mut arr = array![[1, 2, 3], [4, 5, 6]];

// 要素の取得 (参照)
let element = arr[(0, 1)]; // 1行目 (0-indexed), 2列目 (0-indexed)
println!(“Element at (0, 1): {}”, element); // 2

// 要素の変更
arr[(1, 2)] = 10; // 2行目, 3列目の要素を10に変更
println!(“Modified arr:\n{}”, arr);
/ 出力例:
Modified arr:
[[1, 2, 3],
[4, 5, 10]]
/
“`

単一のインデックスでアクセスできるのは、1次元配列の場合か、ndim()が1の場合（例えば形状(N,)や(1, N)、(N, 1)など）です。

2.4. 配列の形状と次元数

配列の形状（shape）と次元数（ndim）は重要な情報です。

“`rust
use ndarray::array;

let arr = array![[[1, 2], [3, 4]], [[5, 6], [7, 8]]];

println!(“Shape: {:?}”, arr.shape()); // Shape: [2, 2, 2]
println!(“Number of dimensions: {}”, arr.ndim()); // Number of dimensions: 3
println!(“Total number of elements: {}”, arr.len()); // Total number of elements: 8
“`

shape()メソッドは、各次元のサイズを含むスライス&[usize]を返します。

2.5. 形状操作 (Reshaping)

配列の要素数を変えずに形状を変更できます。

“`rust
use ndarray::array;
use ndarray::Array;

let arr = array![1, 2, 3, 4, 5, 6];

// 2×3の行列にリシェイプ
let reshaped = arr.into_shape((2, 3)).expect(“Reshape failed”);
println!(“Reshaped:\n{}”, reshaped);
/ 出力例:
Reshaped:
[[1, 2, 3],
[4, 5, 6]]
/

// 3×2の行列にリシェイプ (再度元の配列から行うか、コピー/所有権移動が必要)
// into_shape()は元の配列の所有権を奪うため、再利用するには元の配列に戻すか、クローンが必要です。
let arr2 = array![1, 2, 3, 4, 5, 6];
let reshaped2 = arr2.into_shape((3, 2)).expect(“Reshape failed”);
println!(“Reshaped2:\n{}”, reshaped2);

// 次元数を減らす（flattening）
let arr3 = array![[1, 2], [3, 4]];
let flattened = arr3.into_shape(4).expect(“Reshape failed”); // 1次元に
println!(“Flattened: {}”, flattened); // [1, 2, 3, 4]

// 次元数を増やす（特定の場所に次元を追加）
let arr4 = array![1, 2, 3];
let expanded = arr4.insert_axis(ndarray::Axis(1)); // 1次元目の後 (インデックス1) にサイズ1の次元を追加
println!(“Expanded:\n{}”, expanded);
/ 出力例:
Expanded:
[[1],
[2],
[3]]
/
println!(“Shape of expanded: {:?}”, expanded.shape()); // Shape: [3, 1]

let expanded2 = arr4.insert_axis(ndarray::Axis(0)); // 0次元目 (先頭) にサイズ1の次元を追加
println!(“Expanded2:\n{}”, expanded2);
/ 出力例:
Expanded2:
[[[1, 2, 3]]]
/
println!(“Shape of expanded2: {:?}”, expanded2.shape()); // Shape: [1, 3]
“`

into_shapeメソッドは、可能な限りデータをコピーせずに形状を変更しようとします。ただし、新しい形状が元のデータのメモリレイアウトと互換性がない場合（例えば、C-orderの配列をFortran-orderとして解釈する場合など）、データの再配置（コピー）が必要になることがあります。デフォルトでは、into_shapeはデータコピーを避けるために形状変更の制約を課しますが、into_shape_any()を使えば、コピーが発生する可能性を許容して任意の形状に変更できます。

3. ndarrayが効率的な理由

ndarrayがRustで効率的な数値計算を可能にする理由はいくつかあります。

3.1. メモリレイアウトとストライド

ndarrayの配列は、要素がメモリ上にどのように配置されているかをストライドという概念で管理します。ストライドとは、ある次元のインデックスを1つ進める際に、メモリ上の次の要素にアクセスするためにスキップする必要のあるバイト数（あるいは要素数）です。

例えば、形状(3, 4)の2次元配列を考えます。

  • C-order (行優先): 要素は行ごとに連続してメモリに配置されます。1つの行内の次の要素に進むストライドは要素サイズ分、次の行に進むストライドは「行の要素数 × 要素サイズ」となります。ストライドは[4 * element_size, 1 * element_size]のようになります（要素数単位なら[4, 1]）。

  • Fortran-order (列優先): 要素は列ごとに連続してメモリに配置されます。1つの列内の次の要素に進むストライドは要素サイズ分、次の列に進むストライドは「列の要素数 × 要素サイズ」となります。ストライドは[1 * element_size, 3 * element_size]のようになります（要素数単位なら[1, 3]）。

ndarrayは、配列がC-orderかFortran-orderか、あるいはそれ以外の任意のストライドを持つかを内部で管理できます。特に、into_shapeや.transpose()のような操作は、可能な限りストライド情報だけを変更して、実際のデータコピーを避けます。これにより、大きな配列の形状変更を高速に行うことができます。

“`rust
use ndarray::array;
use ndarray::Array;

let arr = array![[1, 2, 3], [4, 5, 6]]; // デフォルトはC-order

println!(“Array:\n{}”, arr);
println!(“Shape: {:?}”, arr.shape()); // [2, 3]
println!(“Strides: {:?}”, arr.strides()); // [3, 1] (要素数単位)

let transposed = arr.t(); // 転置はビューであり、データはコピーされない

println!(“Transposed:\n{}”, transposed);
println!(“Shape: {:?}”, transposed.shape()); // [3, 2]
println!(“Strides: {:?}”, transposed.strides()); // [1, 3] (要素数単位 – Fortran-orderに対応)
“`

転置された配列transposedは、元の配列arrと同じデータバッファを共有していますが、ストライド情報が異なるため、アクセスパターンが列優先になります。これがビュー（View）の重要な概念です。

3.2. ビューとスライス：ゼロコピー操作

ndarrayのパフォーマンスの鍵の一つは、データコピーを最小限に抑えるビューとスライス機能です。

  • ビュー (.view(), .view_mut()): 配列全体やその一部への参照を提供します。元のデータバッファを共有するため、ビューを作成したり操作したりしてもデータのコピーは発生しません。view()は不変参照（&）に、view_mut()は可変参照（&mut）に相当します。

  • スライス (.slice(), .slice_mut()): ビューの特殊なケースで、特定の軸に沿って配列の一部を切り出します。インデックス範囲やステップを指定できます。

“`rust
use ndarray::{array, Slice};

let mut arr = array![[1, 2, 3], [4, 5, 6], [7, 8, 9]];

// 最初の2行を切り出す (ビュー)
let rows = arr.slice(s![0..2, ..]); // Slice![行の範囲, 列の範囲]
println!(“First two rows:\n{}”, rows);
/ 出力例:
First two rows:
[[1, 2, 3],
[4, 5, 6]]
/

// 2列目を切り出す (ビュー – この場合、形状は (3,) になる)
let col = arr.slice(s![.., 1]);
println!(“Second column: {}”, col); // [2, 5, 8]
println!(“Shape of column: {:?}”, col.shape()); // [3]

// 可変ビューを取得し、要素を変更
let mut mutable_view = arr.slice_mut(s![1, ..]); // 2行目全体
mutable_view[0] = 100; // 元の配列の要素が変更される
println!(“Modified original array:\n{}”, arr);
/ 出力例:
Modified original array:
[[1, 2, 3],
[100, 5, 6],
[7, 8, 9]]
/
“`

スライスに使うs![]マクロは、NumPyのスライス記法（start:stop:step）に似ています。

  • i : 特定のインデックス i
  • i..j : インデックス i から j-1 まで
  • i.. : インデックス i から最後まで
  • ..j : 最初からインデックス j-1 まで
  • .. : 最初から最後まで (全範囲)
  • i..; step : インデックス i から最後まで、ステップ step で
  • ..; step : 最初から最後まで、ステップ step で

ビューとスライスは非常に強力ですが、注意点もあります。ビューは元のデータの借用であるため、ビューが存在する間は元の配列を可変的に操作したり、所有権を移動したりすることはできません。ビューが必要なくなった場合は、スコープから外れるか、明示的に破棄される必要があります。

3.3. 要素型とジェネリクス

ndarrayはジェネリクス(Array<T, D>)を使って実装されているため、Tには数値型（f32, f64, i32, u8など）だけでなく、CloneトレイトとDebugトレイトなどを実装した任意の型を使用できます。これにより、型ごとに異なる配列実装を持つ必要がなく、コードの再利用性が高まります。

パフォーマンスの観点からは、ジェネリクスはコンパイル時に具体的な型に特化される（モノモーフィゼーション）ため、実行時のオーバーヘッドがありません。数値計算では、通常f32やf64が使用されます。

3.4. 内部実装と外部ライブラリ連携 (BLAS/LAPACK)

ndarrayの多くの数値計算操作（特に線形代数演算）は、内部で最適化された外部ライブラリに委譲されます。最も一般的なのは、BLAS (Basic Linear Algebra Subprograms) と LAPACK (Linear Algebra PACKage) です。これらは、長年にわたって高度に最適化されてきた数値計算ライブラリであり、行列乗算や線形システムの解法などを非常に高速に実行できます。

ndarrayは、これらのライブラリのRustバインディングであるndarray-linalgクレートを介して連携します。デフォルトでは、ndarray-linalgはRustで実装されたフォールバックを使用しますが、openblasやnetlibなどのBLAS/LAPACK実装をリンクするようにコンパイルオプションを設定することで、劇的なパフォーマンス向上を得られます。

Cargo.tomlでBLAS/LAPACK実装を指定する例：

toml
[dependencies]
ndarray = "0.15"
ndarray-linalg = { version = "0.14", features = ["openblas"] } # "openblas"の部分を他の実装名に変えられます

これにより、ndarray-linalg経由で呼び出される行列乗算（.dot()）などが、OpenBLASのような最適化されたライブラリによって実行されるようになります。

3.5. 並列化 (rayonとの連携)

Rustの強力な並列処理ライブラリであるrayonとndarrayは密接に連携しています。ndarrayは、rayonを使って多次元配列の要素に対する操作を簡単に並列化できるアダプタを提供しています。

  • .par_iter(): 配列の要素に対するイテレータを並列化します。
  • .par_azip!(): 複数の配列をまとめて（要素ごとに）並列処理します。最も一般的で効率的な並列化方法です。
  • .par_map(), .par_fold()など: rayonが提供する他の並列処理パターンを利用できます。

“`rust
use ndarray::array;
use rayon::prelude::*;

let mut arr = array![[1.0, 2.0, 3.0], [4.0, 5.0, 6.0]];

// 各要素を2倍にする (並列化)
arr.par_map_inplace(|x| x = 2.0);
println!(“Parallel map:\n{}”, arr);

// 複数の配列の要素を並列に足し合わせる
let a = array![1.0, 2.0, 3.0];
let b = array![4.0, 5.0, 6.0];
let mut c = array![0.0, 0.0, 0.0];

// azip! は複数の配列を同時にイテレートするマクロ
// par_azip! はそれを並列化する
ndarray::par_azip!((a in &a, b in &b, c in &mut c) c = a + *b);
println!(“Parallel azip add: {}”, c); // [5.0, 7.0, 9.0]
“`

par_azip!マクロは、複雑な要素ごとの並列操作を非常に簡潔に記述できます。これは、特に大規模な配列に対する演算で大きなパフォーマンス gainsをもたらします。

4. 効率的な数値計算のための具体的なテクニック

ndarrayの基本的な構造と効率の理由を理解したところで、次は具体的なテクニックを見ていきましょう。

4.1. ベクトル化演算とブロードキャスト

ベクトル化演算とは、ループを使って個々の要素にアクセスする代わりに、配列全体やその一部に対して直接演算を行うことです。ndarrayは、NumPyと同様に豊富なベクトル化演算子とメソッドを提供しています。

“`rust
use ndarray::array;

let a = array![[1.0, 2.0], [3.0, 4.0]];
let b = array![[5.0, 6.0], [7.0, 8.0]];

// 要素ごとの加算
let c = &a + &b; // 参照を使う
println!(“Element-wise addition:\n{}”, c);

// 要素ごとの乗算
let d = &a * &b;
println!(“Element-wise multiplication:\n{}”, d);

// スカラーとの演算
let e = &a * 2.0;
println!(“Scalar multiplication:\n{}”, e);

// 行列乗算 (dot product)
use ndarray::linalg::dot;
let f = dot(&a, &b); // または a.dot(&b)
println!(“Matrix multiplication:\n{}”, f);
“`

ブロードキャストは、形状が異なる配列間で要素ごとの演算を行うための機能です。ndarrayはNumPyと同様のブロードキャスト規則に従います。簡単に言えば、次元数が少ない方の配列の形状を、演算が行われる次元に合わせて拡張します。

“`rust
use ndarray::array;

let a = array![[1, 2, 3], [4, 5, 6]]; // 形状 (2, 3)
let b = array![[10], [20]]; // 形状 (2, 1)

// b を形状 (2, 3) にブロードキャストして加算
// [[10, 10, 10], [20, 20, 20]] のように拡張されるイメージ
let c = &a + &b;
println!(“Broadcast addition:\n{}”, c);
/ 出力例:
Broadcast addition:
[[11, 12, 13],
[24, 25, 26]]
/

let d = array![100, 200, 300]; // 形状 (3,)

// d を形状 (2, 3) にブロードキャストして加算
// [[100, 200, 300], [100, 200, 300]] のように拡張されるイメージ
let e = &a + &d;
println!(“Broadcast addition 2:\n{}”, e);
/ 出力例:
Broadcast addition 2:
[[101, 202, 303],
[404, 505, 606]]
/
“`

ブロードキャスト可能な形状かどうかは、次元を末尾から比較し、各次元のサイズが一致するか、どちらかが1であるか、あるいはどちらかの次元が存在しない（次元数が少ない）場合に可能です。ブロードキャストはデータのコピーを伴わない効率的な操作です。

4.2. インデックス操作とビュー/スライス

前述の通り、.slice()や.slice_mut()を使ったスライスは非常に重要です。その他にも、特定の要素や部分配列に効率的にアクセスするためのメソッドがあります。

  • .select(): 特定のインデックスリストに基づいて、特定の軸に沿って要素を選択します。
  • .diag(): 配列の対角要素をビューとして取得します。
  • .row(), .column(), .rows(), .columns(): 行や列をビューとして取得する便利なメソッドです。
  • .take(): 特定のインデックスリストに基づいて要素をコピーして取得します。

“`rust
use ndarray::{array, Axis, s};

let arr = array![[1, 2, 3], [4, 5, 6], [7, 8, 9]];

// 0行目と2行目を選択 (コピーされる)
let selected_rows = arr.select(Axis(0), &[0, 2]);
println!(“Selected rows:\n{}”, selected_rows);
/ 出力例:
Selected rows:
[[1, 2, 3],
[7, 8, 9]]
/

// 対角要素を取得 (ビュー)
let diagonal = arr.diag();
println!(“Diagonal: {}”, diagonal); // [1, 5, 9]

// 1行目を取得 (ビュー)
let row1 = arr.row(1);
println!(“Row 1: {}”, row1); // [4, 5, 6]

// 0列目と2列目を取得 (ビュー)
let cols = arr.columns(s![0..;2]);
println!(“Columns 0 and 2:\n{}”, cols);
/ 出力例:
Columns 0 and 2:
[[1, 3],
[4, 6],
[7, 9]]
/

// 特定の要素をインデックスリストで取得 (コピーされる)
let elements = arr.into_iter().take(3).collect::<Vec<_>>(); // イテレータと組み合わせる例
println!(“First 3 elements: {:?}”, elements); // [1, 2, 3]
“`

ビュー (slice, diag, row, column, rows, columns) はデータコピーが発生しないため、大規模な配列の一部を頻繁に参照する場合に非常に効率的です。一方、selectやtakeは通常データをコピーするため、必要な場合にのみ使用します。

4.3. データ型の選択 (f32 vs f64)

数値計算で浮動小数点数を使う場合、f32（単精度）とf64（倍精度）のどちらを選ぶかはパフォーマンスと精度のトレードオフになります。

  • f64 (倍精度): より広い範囲の数値と高い精度を提供します。科学技術計算や、高い精度が要求されるアプリケーションで一般的に使用されます。メモリ使用量はf32の2倍になります。
  • f32 (単精度): f64より少ないメモリ使用量と、多くの場合高速な計算が可能です。特にSIMD命令（単一の命令で複数のデータを処理するCPU機能）が利用できる場合、f32はf64よりも大幅に高速になることがあります。ディープラーニングなど、多少の精度低下が許容されるが大量の計算が必要な分野でよく使われます。

パフォーマンスを重視する場合、特に大きな配列を扱う際にはf32を検討する価値があります。ただし、数値的に不安定なアルゴリズムや、金融計算など厳密な精度が求められる分野ではf64が必須となるでしょう。ndarrayはどちらの型でも効率的に動作します。

4.4. アロケーションの削減とインプレース操作

Rustではメモリのアロケーション（確保）とデアロケーション（解放）にはコストがかかります。ndarrayを使った計算でパフォーマンスボトルネックになるのは、不要な配列の生成、つまりアロケーションです。

効率的なコードを書くためには、以下の点を考慮します。

  • ビューやスライスを活用する: データのコピーを避けることで、新しい配列のアロケーションを防ぎます。
  • インプレース操作 (_inplaceメソッド): 可能であれば、新しい配列を生成する代わりに、既存の配列をその場で（インプレースで）変更するメソッド（例: map_inplace, par_map_inplace, fill）を使用します。
  • 結果を格納する配列を事前に用意する: 計算結果を格納するための配列を一度だけアロケートし、繰り返し計算に再利用します。
  • 所有権の移動と借用を適切に管理する: into_owned()や.to_owned()はデータコピーを発生させます。本当に所有権が必要か、借用（ビュー）で十分かを検討します。

“`rust
use ndarray::array;

let mut a = array![[1.0, 2.0], [3.0, 4.0]];

// 新しい配列を生成する例 (コピー/アロケーションが発生)
let b = &a * 2.0;

// インプレースで変更する例 (アロケーションなし)
a.map_inplace(|x| x = 2.0); // a自体が変更される
println!(“Inplace multiplication:\n{}”, a);

// 結果を格納する配列を事前に用意する例
let mut result = array![[0.0, 0.0], [0.0, 0.0]];
let c = array![[5.0, 6.0], [7.0, 8.0]];

// result に a + c の結果を格納 (計算によってはインプレースまたは結果配列を指定できるものもある)
// この例では &a + &c が一時配列を作るため、result.assign() がコピーを行う
// ndarray::azip! や par_azip! を使えば、より効率的に複数配列間のインプレース演算ができる
ndarray::azip!((r in &mut result, a in &a, c in &c) r = a + *c);
println!(“Result in pre-allocated array:\n{}”, result);
“`

特に大規模なループ内で配列を生成するコードは、アロケーションコストが積み重なりパフォーマンスを低下させる原因となるため注意が必要です。

4.5. BLAS/LAPACKの活用 (ndarray-linalg)

前述のように、ndarray-linalgクレートと組み合わせてBLAS/LAPACK実装を利用することは、線形代数演算（行列乗算、固有値分解、線形システム解法など）のパフォーマンスを向上させる最も効果的な方法の一つです。これらの外部ライブラリは、プロセッサのキャッシュ効率やSIMD命令などを最大限に活用するように高度に最適化されています。

ndarray-linalgの主な機能：

  • 行列乗算 (.dot())
  • 逆行列 (.inv())
  • 固有値・固有ベクトル計算 (.eig())
  • 線形システム Ax=b の解法 (.solve())
  • 特異値分解 (SVD)
  • LU分解、QR分解、Cholesky分解など

これらの操作は、デフォルトのRust実装よりもBLAS/LAPACKを使う方が圧倒的に高速です。線形代数計算がワークロードの大部分を占める場合は、ndarray-linalgを導入し、BLAS実装を適切にリンクすることを強く推奨します。

“`rust
use ndarray::{array, Array};
use ndarray_linalg::Inverse; // Inverseトレイトをインポートして.inv()を使う

let a = array![[1.0, 2.0], [3.0, 4.0]];

// 逆行列を計算
match a.inv() {
Ok(a_inv) => println!(“Inverse of a:\n{}”, a_inv),
Err(e) => println!(“Error computing inverse: {:?}”, e),
}

use ndarray_linalg::Solve; // Solveトレイトをインポートして.solve()を使う

let b = array![5.0, 11.0];
// Ax = b を解く
match a.solve(&b) {
Ok(x) => println!(“Solution x:\n{}”, x), // x = [1.0, 2.0] となるはず
Err(e) => println!(“Error solving linear system: {:?}”, e),
}
“`

ndarray-linalgのドキュメントを参照し、必要な線形代数操作がサポートされているか確認してください。

5. 高度なトピックと応用

5.1. ブロードキャストの詳細と効率的な使い方

ブロードキャストは非常に便利ですが、その仕組みを理解しておくと、予期しない形状のエラーを防ぎ、効率的なコードを書くのに役立ちます。

ブロードキャスト可能な条件（NumPy互換）：
2つの配列の形状を末尾から比較します。各次元について、以下のいずれかの条件が満たされる必要があります。
1. 次元のサイズが等しい。
2. どちらかの次元のサイズが1である。
3. どちらかの配列にその次元が存在しない（次元数が少ない）。

最も一般的な例は、行列とベクトルの演算です。例えば、形状(M, N)の行列と形状(N,)のベクトルを加算する場合、ベクトルは形状(1, N)として扱われ、M回複製されてから要素ごとの加算が行われるように見えます（実際には複製は行われず、ストライドとインデックス計算で対応します）。また、形状(M, N)の行列と形状(M, 1)のベクトルを加算する場合も同様です。

“`rust
use ndarray::array;

let matrix = array![[1, 2, 3], [4, 5, 6]]; // (2, 3)
let row_vec = array![10, 20, 30]; // (3,) – 行ベクトルとして扱われる
let col_vec = array![[100], [200]]; // (2, 1) – 列ベクトルとして扱われる

let result1 = &matrix + &row_vec; // (2,3) + (3,) -> (2,3) + (1,3) ブロードキャスト -> (2,3)
println!(“Matrix + Row Vec:\n{}”, result1);

let result2 = &matrix + &col_vec; // (2,3) + (2,1) -> (2,3) ブロードキャスト -> (2,3)
println!(“Matrix + Col Vec:\n{}”, result2);
“`

ブロードキャストはメモリコピーなしで行われるため、パフォーマンス上のメリットが大きいです。しかし、意図しないブロードキャストはバグの原因になることもあるため、特に複雑な形状の配列を扱う際は、期待するブロードキャストが起こるか形状をよく確認することが重要です。形状不一致の場合はコンパイル時または実行時にエラーになります。

5.2. Generalized Universal Functions (GUFuncs)

GUFuncsは、NumPyで導入された概念で、異なる形状を持つ配列に対して、要素ごとではなく「サブ配列ごと」に適用される操作を一般化するものです。例えば、行列のバッチに対して行列乗算を行う場合などが該当します。ndarrayにも同様の概念と機能が導入されています。

これは、azip!マクロのより高度な形式と考えることができます。特に、複数の入力配列と出力配列があり、それぞれの次元の一部だけを使って計算を行い、残りの次元は「バッチ」として扱いたい場合に有用です。

GUFuncsを直接使うことは少ないかもしれませんが、ndarrayの内部で並列化やブロードキャストと組み合わせて効率的な操作を実現するために使われています。カスタムの複雑な操作を実装する際に、この概念を理解しておくと、より効率的な設計につながる可能性があります。

5.3. 疎行列 (ndarray-sparse)

数値計算では、多くの要素がゼロである「疎行列」を扱うことがよくあります。密行列と同じ方法で格納すると、メモリが無駄になり、ゼロとの演算に計算時間が費やされます。疎行列を効率的に扱うためには、非ゼロ要素とその位置だけを格納する専用のデータ構造が必要です。

ndarray-sparseクレートは、ndarrayと連携して疎行列を扱うための機能を提供します。CSR (Compressed Sparse Row) やCSC (Compressed Sparse Column) のような一般的な疎行列形式をサポートしています。

toml
[dependencies]
ndarray = "0.15"
ndarray-sparse = "0.2" # 最新バージョンを指定

“`rust
use ndarray::{Array, array};
use ndarray_sparse::{SparseArray, spmatrix};
use ndarray_sparse::coo::CooMatrix; // 座標リスト形式

// 疎行列の生成 (COO形式)
let mut coo_mat: CooMatrix = CooMatrix::new((3, 3));
coo_mat.add_element(0, 1, 1.0);
coo_mat.add_element(1, 0, 2.0);
coo_mat.add_element(2, 2, 3.0);

// CSR形式に変換 (計算に適した形式)
let csr_mat = coo_mat.to_csr();
println!(“Sparse Matrix (CSR):\n{:?}”, csr_mat);

// 密行列との乗算
let dense_vec = array![1.0, 1.0, 1.0];
// 疎行列とベクトル/行列の積はndarray-sparseの機能を使う
let result = csr_mat.dot(&dense_vec);
println!(“Sparse matrix-vector product: {:?}”, result); // [1.0, 2.0, 3.0]
“`

疎行列を扱う場合は、ndarray-sparseのような専用ライブラリを使用することで、メモリ使用量と計算時間の両方を大幅に削減できます。

5.4. 科学計算ライブラリとの連携

ndarrayは、Rustで他の科学計算ライブラリと連携するための基盤としても機能します。

  • ndarray-stats: 配列に対する統計計算（平均、分散、パーセンタイルなど）を提供します。
  • ndarray-ode: 常微分方程式 (ODE) ソルバー。ndarray配列を状態ベクトルとして扱います。
  • rust-pca: 主成分分析 (PCA) ライブラリ。ndarray配列を入力として受け取ります。

これらのライブラリは、ndarray配列を入力または出力として受け取ることで、ndarrayエコシステム内で統一的にデータを取り扱えるように設計されています。

5.5. GPU計算

現状、ndarray自体は直接GPU上で計算を実行する機能を持っていません。ndarrayの配列はCPUメモリ上に存在します。しかし、RustでGPU計算を行うためのライブラリ（例: wgpu, vulkano, cuda-sysなど）は存在しており、ndarrayはこれらのライブラリとデータをやり取りする際の橋渡し役として機能できます。

例えば、CPUで準備したデータをndarray配列に格納し、それをGPUライブラリが理解できる形式（例えば、GPUメモリ上のバッファ）に変換してGPUに転送する、といったワークフローが考えられます。計算結果をGPUからCPUに戻す際も、同様にGPUバッファからndarray配列にデータをコピーします。

将来的には、ndarrayとGPUライブラリ間の連携がよりスムーズになる可能性もありますが、現状では手動でのデータ転送が必要です。GPUを使った計算は、特に大規模な並列計算（行列乗算、畳み込みなど）で大きなパフォーマンス向上をもたらしますが、GPUメモリへのデータ転送コストも考慮する必要があります。

6. パフォーマンス測定とベンチマーク

コードのパフォーマンスを評価し、最適化のボトルネックを特定するためには、正確なベンチマークが不可欠です。「なんとなく速い」ではなく、具体的な数値で比較することが重要です。Rustでは、criterionクレートがベンチマークの標準的なツールとして広く使われています。

6.1. criterionクレート

criterionは、統計的手法を用いてベンチマーク結果のばらつきを考慮し、信頼性の高いパフォーマンスレポートを生成します。

Cargo.tomlに依存関係を追加します。

toml
[dev-dependencies]
criterion = { version = "0.4", features = ["html_reports"] } # 最新バージョンを指定

ベンチマークコードは通常 benches/ ディレクトリに配置します。benches/my_benchmark.rs のようなファイルを作成します。

“`rust
// benches/my_benchmark.rs
use criterion::{black_box, criterion_group, criterion_main, Criterion};
use ndarray::{Array, Array2};
use ndarray::array;
use rand::prelude::*;

// 行列乗算のベンチマーク関数
fn matrix_multiply_benchmark(c: &mut Criterion) {
let mut rng = rand::thread_rng();
let size = 100; // 行列サイズ

let a: Array2<f64> = Array::random((size, size), &mut rng);
let b: Array2<f64> = Array::random((size, size), &mut rng);

c.bench_function("matrix_multiply", |bencher| {
    bencher.iter(|| {
        // ベンチマーク対象のコード
        // black_box() はコンパイラによる最適化（計算結果の除去など）を防ぐ
        black_box(a.dot(&b));
    });
});

}

// Vecを使った単純な行列乗算 (比較用)
fn naive_matrix_multiply(c: &mut Criterion) {
let mut rng = rand::thread_rng();
let size = 100;

// Vec<Vec<f64>> で行列を表現
let a: Vec<Vec<f64>> = (0..size)
    .map(|_| (0..size).map(|_| rng.gen::<f64>()).collect())
    .collect();
let b: Vec<Vec<f64>> = (0..size)
    .map(|_| (0..size).map(|_| rng.gen::<f64>()).collect())
    .collect();

c.bench_function("naive_matrix_multiply", |bencher| {
    bencher.iter(|| {
        let mut c = vec![vec![0.0; size]; size];
        for i in 0..size {
            for j in 0..size {
                for k in 0..size {
                    c[i][j] += a[i][k] * b[k][j];
                }
            }
        }
        black_box(c);
    });
});

}

// ベンチマークグループの定義
criterion_group!(benches, matrix_multiply_benchmark, naive_matrix_multiply);
criterion_main!(benches);
“`

ベンチマークを実行するには、cargo bench コマンドを使います。結果は target/criterion/report/index.html にHTMLレポートとして出力されます。

この例のように、ndarrayの.dot()メソッド（もしBLASが有効になっていればさらに高速）と、Vec<Vec<f64>>を使った3重ループの素朴な実装を比較することで、ndarrayがどれだけ効率的かを定量的に評価できます。

6.2. ボトルネックの特定

ベンチマーク結果を見て、最も時間がかかっている部分（ボトルネック）を特定します。数値計算では、以下のような点がボトルネックになりやすいです。

  • アロケーション: 不要な配列の生成。
  • データコピー: ビューで済むところでコピーをしてしまう。
  • キャッシュミス: メモリ上のデータへのアクセスパターンが非効率的。
  • 並列化不足: 並列化可能な部分がシングルスレッドで実行されている。
  • 非効率なアルゴリズム: 使用しているアルゴリズム自体が計算量の大きいもの。
  • 外部ライブラリの未利用: BLAS/LAPACKで高速化できる線形代数計算をRustネイティブコードで行っている。

これらのボトルネックが見つかったら、前に述べた効率化テクニック（ビュー、インプレース操作、並列化、BLAS活用など）を適用してコードを改善し、再度ベンチマークを実行して効果を確認するという iterative なプロセスで最適化を進めます。

7. 実践例

いくつかの簡単な実践例を通して、ndarrayを使った効率的な数値計算の具体的なコードを見てみましょう。

7.1. 行列乗算 (BLAS活用)

これはベンチマーク例でも触れましたが、BLASを有効にしたndarray-linalgを使った行列乗算は非常に高速です。

“`rust
use ndarray::{array, Array2};
use ndarray_linalg::Dot; // Dotトレイトをインポートして.dot()を使う

let a: Array2 = array![[1.0, 2.0], [3.0, 4.0]];
let b: Array2 = array![[5.0, 6.0], [7.0, 8.0]];

// BLASが有効なら、BLASのdgemm関数などが呼ばれる
let c = a.dot(&b);
println!(“Matrix product:\n{}”, c);
/ 出力:
[[19.0, 22.0],
[43.0, 50.0]]
/
``Cargo.tomlでndarray-linalg = { version = “…”, features = [“openblas”] }` のように指定することで、OpenBLASなどの最適化されたライブラリを利用できます。

7.2. 統計計算 (平均、標準偏差)

ndarray-statsクレートを使うと、統計計算を簡単に行えます。

toml
[dependencies]
ndarray = "0.15"
ndarray-stats = "0.5" # 最新バージョンを指定

“`rust
use ndarray::{array, Array1, Array2, Axis};
use ndarray_stats::DeviationExt; // .mean(), .std() などのトレイト

let data = array![[1.0, 2.0, 3.0], [4.0, 5.0, 6.0], [7.0, 8.0, 9.0]]; // 3×3行列

// 配列全体の平均
let mean_total = data.mean().unwrap();
println!(“Mean of all elements: {}”, mean_total); // 5.0

// 各列の平均 (軸0に沿って計算)
let mean_per_column = data.mean_axis(Axis(0)).unwrap();
println!(“Mean per column: {}”, mean_per_column); // [4.0, 5.0, 6.0]

// 各行の平均 (軸1に沿って計算)
let mean_per_row = data.mean_axis(Axis(1)).unwrap();
println!(“Mean per row: {}”, mean_per_row); // [2.0, 5.0, 8.0]

// 配列全体の標準偏差 (標本標準偏差 – ddof=1)
let std_total = data.std(1.0); // ddof=1.0 で標本標準偏差
println!(“Standard deviation of all elements: {}”, std_total); // 約 2.58

// 各行の標準偏差
let std_per_row = data.std_axis(Axis(1), 1.0);
println!(“Standard deviation per row: {}”, std_per_row); // [1.0, 1.0, 1.0]
``ndarray-stats`は、これらの統計計算を効率的に（可能な場合は並列化も考慮して）実装しています。

7.3. 画像処理 (フィルタリング)

画像は通常、高さ×幅×チャンネル（例: RGB）の3次元配列として表現できます。ndarrayはこのような多次元データを扱うのに適しています。簡単な例として、画像の各ピクセル値を2倍にする操作を考えます。

“`rust
use ndarray::Array3;
use rayon::prelude::*; // 並列化にrayonを使う

// ダミーの画像データ (高さ=2, 幅=3, チャンネル=3)
// u8型は画像データでよく使われる
let mut image: Array3 = Array3::from_shape_fn((2, 3, 3), |(h, w, c)| {
((h * 3 + w) * 3 + c + 1) as u8
});
println!(“Original image:\n{:?}”, image);

// 各ピクセル値を2倍にする (飽和演算付き) – 並列化
image.par_map_inplace(|x| x.saturating_mul(2));
println!(“Processed image:\n{:?}”, image);
“`

より複雑なフィルタリング（例: 畳み込み）も、ndarrayのビューと配列演算を使って効率的に実装できます。例えば、カーネル行列を用意し、.slice()で画像の各部分を切り出し、行列乗算や要素ごとの積和を計算するといった手法です。ndarrayはスライディングウィンドウ操作なども効率的に行えるメソッドを提供しています。

7.4. 簡単なニューラルネットワーク層 (活性化関数)

ニューラルネットワークの計算は多次元配列演算の典型例です。ここでは、ReLU活性化関数（max(0, x)）を配列の各要素に適用する例を示します。

“`rust
use ndarray::{array, Array};
use rayon::prelude::*;

let mut input = array![[-1.0, 0.5, -0.1], [2.0, -0.5, 1.0]]; // 形状 (2, 3)

// ReLU活性化関数を適用 (インプレース、並列化)
input.par_map_inplace(|x| x = x.max(0.0));
println!(“ReLU output:\n{}”, input);
/ 出力:
[[0.0, 0.5, 0.0],
[2.0, 0.0, 1.0]]
*/
“`

行列乗算（a.dot(&w)）による線形変換と、要素ごとの活性化関数適用（par_map_inplaceなど）は、ニューラルネットワークの基本的な層を構成する主要な計算です。ndarrayはこれらの計算を効率的に実行するための優れたツールです。

8. よくある落とし穴とデバッグ

ndarrayは強力ですが、Rustの借用システムや多次元配列特有の概念（形状、軸、ストライドなど）が組み合わさることで、最初は戸惑うこともあるかもしれません。

8.1. 形状不一致のエラー

最もよく遭遇するエラーの一つは、形状が一致しない配列間で演算を行おうとしたときです。ndarrayは多くのケースでコンパイル時または実行時に形状の不一致を検出してエラーを報告します。

“`rust
use ndarray::array;

let a = array![[1, 2], [3, 4]]; // (2, 2)
let b = array![[5, 6, 7], [8, 9, 10]]; // (2, 3)

// let c = &a + &b; // 形状が異なるためコンパイルエラーまたは実行時パニック
// error[E0277]: cannot add &ndarray::Array<i32, ndarray::Dim<[usize; 2]>> to &ndarray::Array<i32, ndarray::Dim<[usize; 2]>>
// 異なる形状のArrayに対するAddトレイトが実装されていないため

let d = array![10, 20]; // (2,)

// let e = &a + &d; // この加算はブロードキャスト可能なのでエラーにならない
// [[11, 12], [33, 34]] となる (dが [[10, 10], [20, 20]] にブロードキャストされる)
“`

エラーメッセージだけでは分かりにくい場合もありますが、「shape mismatch」や「dimensions not compatible」といったキーワードに注目し、arr.shape()やarr.ndim()を使って配列の形状を確認することがデバッグの第一歩です。

8.2. 借用チェッカーとの戦い (.into_owned(), .to_owned())

ビューやスライスは元の配列の借用であるため、Rustの借用ルールが適用されます。可変ビュー (&mut) が存在する場合、同じ配列への他の参照や可変参照、あるいは所有権を奪う操作はできません。

“`rust
use ndarray::array;
use ndarray::{Array, Array2, s};

let mut arr = array![[1, 2], [3, 4]];

let view1 = arr.slice_mut(s![0, ..]); // arrへの可変参照

// let view2 = arr.slice_mut(s![1, ..]); // エラー: arrは既に可変参照されている
// let owned_arr = arr.to_owned(); // エラー: arrは既に可変参照されている
“`

計算結果を新しい配列として取得したいが、元の配列のビューしか持っていない場合など、所有権が必要になることがあります。そのような場合は、.to_owned()や.into_owned()メソッドを使ってデータをコピーし、新しい独立した配列を作成します。

“`rust
use ndarray::{array, s};

let arr = array![[1, 2], [3, 4]];
let view = arr.slice(s![0, ..]); // arrへの不変参照

// view自体は参照なので、要素を変更したり所有権を移動したりできない
// let mut view_mut = view.slice_mut(s![0]); // エラー: 不変ビューからは可変ビューを取得できない

// viewのデータをコピーして独立した1次元配列にする
let owned_vec = view.to_owned(); // Array1 型になる
let mut owned_vec_mut = owned_vec.into_owned(); // owned_vec の所有権を移動して可変に
owned_vec_mut[0] = 100;
println!(“Owned vec: {}”, owned_vec_mut); // [100, 2] (元のarrは変更されない)

println!(“Original arr:\n{}”, arr);
“`

to_owned()やinto_owned()はデータコピーを伴うため、パフォーマンスに影響を与える可能性があります。本当に所有権が必要か、あるいは計算結果を格納するために新しい配列を事前にアロケートし、そこに計算結果を書き込む形にできないか検討することが重要です。

8.3. ビューの有効期間 (Lifetime)

ビューは元の配列の借用であるため、元の配列よりも長く生存することはできません。これはRustの通常のライフタイムルールに従います。

“`rust
use ndarray::{Array, Array2, s};

fn get_first_row<‘a>(arr: &’a Array2) -> ndarray::ArrayView1<‘a, i32> {
arr.row(0) // arrのライフタイム<‘a>を持つビューを返す
}

// main関数など
let arr = Array2::zeros((3, 3));
let view = get_first_row(&arr); // viewはarrのライフタイムに依存する
// arrは view がスコープを抜けるまで有効でなければならない
“`

通常はコンパイラがライフタイムを推論してくれますが、複雑な関数やデータ構造でビューを扱う場合は、明示的なライフタイム注釈が必要になることがあります。

8.4. デバッグ方法

  • 形状と値の確認: println!マクロを使って、配列の形状 (.shape()) や値 ({} または {:?}) を出力します。特に形状不一致エラーの場合は、各配列の形状を確認するのが必須です。
  • ステップ実行: デバッガー (GDB, LLDBなど) を使ってコードをステップ実行し、配列の状態（形状、ストライド、データポインタなど）を確認します。
  • サブ配列の確認: 特定の部分で問題が発生していると思われる場合は、.slice()や.row()などを使ってその部分を切り出し、個別に確認します。

ndarrayはDebugトレイトを実装しているため、println!("{:?}, array)のように出力すると、形状、ストライド、および要素の値が整形されて表示され、デバッグに役立ちます。

9. まとめ

この記事では、Rustで数値計算を効率的に行うための強力なライブラリであるndarrayについて、その基本的な使い方から、なぜ効率的なのか、具体的なテクニック、そして応用例や注意点に至るまで、詳細に解説しました。

ndarrayは、多次元配列を扱うための直感的でNumPyライクなAPIをRustに提供します。その内部では、効率的なメモリレイアウト、データコピーを避けるビューとスライス、BLAS/LAPACKとの連携、そしてrayonによる簡単な並列化といった仕組みが組み合わされており、Rustの強みである安全性とパフォーマンスを最大限に活かせるよう設計されています。

効率的な数値計算を実現するためには、以下の点を意識することが重要です。

  • ベクトル化演算とブロードキャストを積極的に使う: ループよりも高速で記述も簡潔になります。
  • ビューとスライスを使い、データコピーを避ける: 特に大きな配列の一部にアクセスする場合に重要です。
  • インプレース操作を活用し、アロケーションを減らす: 不要なメモリ確保・解放のコストを削減します。
  • 線形代数計算にはndarray-linalgとBLAS/LAPACKを利用する: 高度に最適化された外部ライブラリの力を借ります。
  • rayonと.par_azip!などを使って並列化する: マルチコアCPUの性能を引き出します。
  • 適切なデータ型 (f32 vs f64) を選択する: 精度とパフォーマンスのバランスをとります。
  • 疎行列にはndarray-sparseを利用する: 特化したデータ構造で効率を高めます。
  • ベンチマークでパフォーマンスを定量的に評価し、ボトルネックを特定する: 勘ではなく、データに基づいた最適化を行います。

Rustは、安全性、並行性、そしてパフォーマンスを高いレベルで両立できるユニークな言語です。ndarrayのような成熟したライブラリが登場したことで、Rustはシステムプログラミングの領域を超え、数値計算や科学計算の分野でも実用的な選択肢となりつつあります。

もしあなたがRustを使って数値計算プロジェクトに取り組んでいる、あるいはこれから始めようとしているなら、ぜひndarrayを試してみてください。その効率性と柔軟性は、きっとあなたの開発体験とコードのパフォーマンスを向上させるでしょう。

10. 参考文献とリソース

これらのリソースを活用しながら、ndarrayを使った効率的な数値計算のスキルをさらに磨いていってください。

コメントする

メールアドレスが公開されることはありません。 ※ が付いている欄は必須項目です

上部へスクロール