[Rust/Python] 周期境界条件

実数 $\mathbb{R}$ に周期 $b > 0$ の周期境界条件を課す ($\mathbb{R} / b \mathbb{R}$) ことを考えます. つまり, ふたつの数 $x$ と $x+b$ を同一視します. このとき任意の数 $x \in \mathbb{R}$ に対応する座標の一意的な表現として

• $0 \leq x < b$
• $-\frac{b}{2} \leq x < \frac{b}{2}$
  という二通りの表現方法が可能です. 本記事では任意の実数 $x$ が与えられたときにこれらの制約を満足する値へと変換するコードを確認します.

Rust 実装

コード全体と動作確認はこちら: PlayGround

0 ≦ x < b

Rust では剰余 x%b は x と同符号になります. 従って x が正ならば x%b でよいのですが, そうでなければ b を足す必要があります.

#[inline]
fn pbc1(x: f64, b: f64) -> f64 {
    let dx = x%b;
    if dx >= 0. { dx } else { dx + b }
}

最初 f64::is_sign_negative を使っていたのですが, -0.0 の扱いで困った結果, 不等号を使うべきだという結論に到達しました.

-b/2 ≦ x < b/2

この場合 (x/b + 0.5).floor() を計算すると x を欲しい区間に移動するために何回 b で動かせばよいのかがわかります.

#[inline]
fn pbc2(x: f64, b: f64) -> f64 {
    x - (x/b + 0.5).floor()*b
}

Python 実装

以下のコードは入力 x が float でも np.ndarray でも使えます. 動作確認はこちら: Paiza

0 ≦ x < b

Python では剰余は $0 \leq x%b < b$ を満たすようになっている (いま $b > 0$ としています) ので一発です.

def pbc1(x, b):
    return x%b

-b/2 ≦ x < b/2

これは上の Rust 実装をそのまま Python に書き換えるだけです.

def pbc2(x, b):
    return x - np.floor(x/b + 0.5)*b