ツイスパラメーター

分布生成

空間電荷効果

ブリルアンフロー

加速空洞

Panofsky-Wenzel
theorem

電磁石

真空

研究加速器理論ツイスパラメーター: 分布生成

ツイスパラメーター分布生成

ツイスパラメーターと RMS エミッタンスから，マクロ粒子の二次元正規分布を作る具体的な手順を説明する．

はじめに

ツイスパラメーターから粒子分布を作る作業は，ビームダイナミクス計算の入口である．入力として欲しいのは，一平面のツイスパラメーター \((\alpha,\,\beta)\) と RMS エミッタンス \(\varepsilon_{rms}\) であり，出力として欲しいのは，その統計量を持つ多数のマクロ粒子 \((x_i,\,x_i^\prime)\) である．

楕円の意味，二次元正規分布の基本式，統計量とツイスパラメーターの関係は，すでに楕円，正規分布，統計のページで説明した．このページではそれらを繰り返さず，分布をどう作るかに話を絞る．

結論だけ先に言うと，標準二次元正規分布を作り，目標の共分散行列のコレスキー分解で得た下三角行列を掛ければよい．この方法は実装が簡単で，得られる分布の平均，分散，共分散をそのまま確認できる．平均がゼロでないビームを作りたい場合は，最後に平行移動を加えればよい．

必要な線形代数の知識

基本

考え方の中心はとても単純である．平均 0，共分散行列が単位行列の二次元確率変数を

\begin{align} \boldsymbol{z}_i= \begin{bmatrix} u_i \\ v_i \end{bmatrix}, \qquad \mathrm{Cov}(\boldsymbol{z}_i)= \begin{bmatrix} 1 & 0 \\ 0 & 1 \end{bmatrix} \label{eq:standard_normal_vector} \end{align}

とする．そして，目標の共分散行列を \(\Sigma_{target}\) とし，ある行列 \(L\) が

\begin{align} \Sigma_{target}=LL^T \label{eq:cholesky_basic} \end{align}

を満たしているとする．このとき，

\begin{align} \boldsymbol{w}_i=L\boldsymbol{z}_i \label{eq:linear_transform} \end{align}

で作った新しい確率変数の共分散行列は

\begin{align} \mathrm{Cov}(\boldsymbol{w}_i) = L\,\mathrm{Cov}(\boldsymbol{z}_i)\,L^T = LIL^T = LL^T = \Sigma_{target} \label{eq:cov_linear_transform} \end{align}

になる．つまり，単位円に対応する標準正規分布を，行列 \(L\) で引き伸ばして傾ければ，欲しい分布が得られる．ここで使う \(L\) を，正定値対称行列に対して系統的に求める方法がコレスキー分解である．

ツイスパラメーターの共分散行列のコレスキー分解

ツイスパラメーターと RMS エミッタンスが表す目標の共分散行列は，

\begin{align} \Sigma_{tw} = \varepsilon_{rms} \begin{bmatrix} \beta & -\alpha \\ -\alpha & \gamma \end{bmatrix}, \qquad \beta\gamma-\alpha^2=1 \label{eq:twiss_covariance} \end{align}

である．通常は \(\beta>0\) と \(\varepsilon_{rms}>0\) であり，さらに \(\beta\gamma-\alpha^2=1\) なので，この行列は正定値である．したがってコレスキー分解ができる．

下三角行列を \(L_{tw}\) とすると，

\begin{align} L_{tw} = \sqrt{\varepsilon_{rms}} \begin{bmatrix} \sqrt{\beta} & 0 \\ -\cfrac{\alpha}{\sqrt{\beta}} & \sqrt{\gamma-\cfrac{\alpha^2}{\beta}} \end{bmatrix} \label{eq:twiss_cholesky_general} \end{align}

であり，ツイスパラメーターの拘束条件を使うと

\begin{align} L_{tw} = \sqrt{\varepsilon_{rms}} \begin{bmatrix} \sqrt{\beta} & 0 \\ -\cfrac{\alpha}{\sqrt{\beta}} & \cfrac{1}{\sqrt{\beta}} \end{bmatrix} \label{eq:twiss_cholesky_simple} \end{align}

と簡単になる．実装上はこちらの形が便利である．\(\gamma\) を独立に入力するより，\(\gamma=(1+\alpha^2)/\beta\) として計算した方が，パラメーターの不整合を避けやすい．

二次元正規分布の作成方法

ここでは，一平面 \((x,\,x^\prime)\) の分布を作る．\(y\) 平面や \((\phi,\,\Delta E)\) に対しても，全く同じ手順を使えばよい．

手続き (1) 標準二次元正規分布

最初に，互いに独立な標準正規乱数 \(u_i\) と \(v_i\) を \(N\) 個ずつ生成する．これを

\begin{align} \boldsymbol{z}_i= \begin{bmatrix} u_i \\ v_i \end{bmatrix}, \qquad u_i \sim \mathcal{N}(0,\,1), \qquad v_i \sim \mathcal{N}(0,\,1) \label{eq:standard_pair} \end{align}

と書く．この分布は原点中心の真ん丸なガウス分布であり，共分散行列は単位行列である．言い換えると，位置方向にも角度方向にも同じ広がりを持ち，相関がない状態から出発することになる．

ここで重要なのは，作りたいのが楕円の輪郭ではなく，楕円状のガウス分布だという点である．したがって，半径と角度から楕円周上の点を作る方法ではなく，まず正規乱数を使って雲のような点群を作る必要がある．

手続き (2) コレスキー分解

次に，入力された \(\alpha\)，\(\beta\)，\(\varepsilon_{rms}\) から式 \eqref{eq:twiss_cholesky_simple} の下三角行列 \(L_{tw}\) を作る．この行列が，標準正規分布の真ん丸な雲を，目的の広がりと傾きを持つ雲へ変形する役目を持つ．

行列の形をもう一度書くと，

\begin{align} L_{tw} = \sqrt{\varepsilon_{rms}} \begin{bmatrix} \sqrt{\beta} & 0 \\ -\cfrac{\alpha}{\sqrt{\beta}} & \cfrac{1}{\sqrt{\beta}} \end{bmatrix} \label{eq:twiss_cholesky_for_use} \end{align}

である．1 行目は位置 \(x\) の広がりを決め，2 行目は \(x\) と \(x^\prime\) の相関，および角度方向の広がりを決める．\(\alpha\) が 0 なら傾きはなくなり，\(\alpha\) が負なら右上がり，正なら右下がりの相関を作る．

手続き (3) 座標変換

最後に，各粒子に対して

\begin{align} \begin{bmatrix} x_i \\ x_i^\prime \end{bmatrix} = L_{tw} \begin{bmatrix} u_i \\ v_i \end{bmatrix} \label{eq:final_transform} \end{align}

を計算する．成分で書けば，

\begin{align} x_i &= \sqrt{\beta\varepsilon_{rms}}\,u_i \label{eq:generated_x} \\ x_i^\prime &= \sqrt{\cfrac{\varepsilon_{rms}}{\beta}}\,(v_i-\alpha u_i) \label{eq:generated_xp} \end{align}

である．これで，標準正規分布のサンプルは，目標のツイスパラメーターを持つ分布へ変換される．

この式から，生成された分布の二次モーメントはすぐに確かめられる．\(u_i\) と \(v_i\) は独立で，平均 0，分散 1 なので，

\begin{align} \langle x^2 \rangle &= \beta\varepsilon_{rms} \\ \langle xx^\prime \rangle &= -\alpha\varepsilon_{rms} \\ \langle x^{\prime 2} \rangle &= \cfrac{1+\alpha^2}{\beta}\varepsilon_{rms} = \gamma\varepsilon_{rms} \label{eq:generated_moments} \end{align}

となる．したがって，共分散行列は確かに式 \eqref{eq:twiss_covariance} になる．言い換えると，散布図を描いたときの RMS 楕円は，入力したツイスパラメーターそのものになる．

実装だけが必要なら，次のような簡単な Python プログラムで十分である．sample_twiss.py などの名前で保存すれば，そのまま実行できる．

#!/usr/bin/env python3
import numpy as np

def sample_twiss(alpha, beta, eps_rms, n, x0=0.0, xp0=0.0, seed=1234):
    if beta &lt;= 0.0:
        raise ValueError("beta must be positive.")
    if eps_rms &lt;= 0.0:
        raise ValueError("eps_rms must be positive.")

    rng = np.random.default_rng(seed)
    u = rng.standard_normal(n)
    v = rng.standard_normal(n)

    x = np.sqrt(beta * eps_rms) * u + x0
    xp = np.sqrt(eps_rms / beta) * (v - alpha * u) + xp0
    return np.column_stack((x, xp))


def calc_twiss(particles):
    x = particles[:, 0] - np.mean(particles[:, 0])
    xp = particles[:, 1] - np.mean(particles[:, 1])

    xx = np.mean(x * x)
    xxp = np.mean(x * xp)
    xpxp = np.mean(xp * xp)
    eps_rms = np.sqrt(xx * xpxp - xxp ** 2)

    alpha = -xxp / eps_rms
    beta = xx / eps_rms
    gamma = xpxp / eps_rms
    return alpha, beta, gamma, eps_rms


def main():
    alpha = 1.65
    beta = 0.051
    eps_rms = 2.894352e-5
    n = 10000

    particles = sample_twiss(alpha, beta, eps_rms, n)
    np.savetxt("particles.txt", particles, fmt="%.6e", header="x xp")

    alpha2, beta2, gamma2, eps2 = calc_twiss(particles)
    print(f"alpha   = {alpha2:.6f}")
    print(f"beta    = {beta2:.6f}")
    print(f"gamma   = {gamma2:.6f}")
    print(f"eps_rms = {eps2:.6e}")
    print("particles.txt written")


if __name__ == "__main__":
    main()

この例は numpy だけで動作する．particles.txt には 1 列目が \(x\)，2 列目が \(x^\prime\) として保存される．平均位置 \(x_0\) や平均角度 \(x_0^\prime\) を持たせたい場合は，sample_twiss() の引数に値を与えればよい．

実装上の注意

生成されるのは楕円の境界線ではない．得られるのは楕円状のガウス分布であり，RMS 楕円はその統計量を表す補助線である．
粒子数が有限だと統計量は少し揺らぐ．\(N\) が小さいと，計算し直した \(\alpha\)，\(\beta\)，\(\varepsilon_{rms}\) は入力値と完全一致しない．粒子数を増やせばずれは小さくなる．
入力の整合性を確認する．少なくとも \(\beta>0\) と \(\varepsilon_{rms}>0\) は必須である．\(\alpha\)，\(\beta\)，\(\gamma\) を別々に与えるなら，\(\beta\gamma-\alpha^2=1\) も確認したい．
独立な 3 平面なら同じ手続きを 3 回行えばよい．\(x\)，\(y\)，\((\phi,\,\Delta E)\) を独立とみなす場合は，各平面ごとに同じ計算を行えばよい．平面間の相関まで入れたい場合は，6 次元の共分散行列を直接扱う必要がある．

まとめ

ツイスパラメーターから二次元正規分布を作る手順は，次の 3 段階である．

独立な標準正規乱数から，原点中心で相関のない二次元分布を作る．
ツイスパラメーターが表す共分散行列をコレスキー分解し，下三角行列 \(L_{tw}\) を求める．
\(\boldsymbol{w}_i=L_{tw}\boldsymbol{z}_i\) により座標変換し，必要なら最後に平均値を加える．

この方法の利点は，理屈が明快で，実装も短く，生成後の統計量の確認も容易な点にある．より詳しい導出や，実際のプログラムと入出力ファイルまで含めた説明は，同じディレクトリに置いた PDFを参照されたい．

ページ情報

参考資料

山本昌志, 「ツイスパラメーターとマクロ粒子分布」, twiss_parameter.pdf, 2026年2月15日．
G. Strang, Linear Algebra and Its Applications．

更新履歴

2026年03月07日

twiss_parameter.pdf と twiss_parameter.tex をもとに，分布生成の手順を中心にページを再構成