anharmonic プログラム仕様

D:/git/sphinx/tkProg/source/tiny_simulations/anharmonic.py

anharmonic_chain_Tshift.py

1次元単原子チェーン（最隣接）の非調和フォノン振動数シフトの温度依存性を計算するモジュール。

詳細説明: モデル（相対変位 Δ_n = u_{n+1} - u_n）:

V = Σ_n [ (K2/2) Δ_n^2 + (K3/3!) Δ_n^3 + (K4/4!) Δ_n^4 ]

摂動論を用いて以下のように近似計算を行います:

ω_q(T) ≈ ω_q + Δω_q^(4)(T) + Δω_q^(3)(T)

ここで、

Δω_q ≈ Re Σ(q, ω_q) / (2 ω_q) となります。

• 4次 (K4) の「tadpole（オタマジャクシ）」型シフト: ⟨|u_k|^2⟩ の単純な熱平均によって生じます。

• 3次 (K3) の「bubble（バブル）」型シフト: kに関する主値積分で、微小パラメータetaによって正則化されます。

単位系: - デフォルトでSI単位系を使用します（aはメートル、Mはkg、K2はN/m、K3はN/m^2、K4はN/m^3）。 - 振動数は内部的にrad/sで計算され、必要に応じてTHzに変換して出力可能です。

ユーザーは、k2_func/k3_func/k4_funcを編集することで、K2/K3/K4を (delta, T) の関数としてカスタマイズできます。

tiny_simulations.anharmonic.C_of_q(q: ndarray, a: float) → ndarray

フーリエ変換の位相因子 C(q) = e^{iqa} - 1 を計算する。

詳細説明: 計算の多くでは実数となる |C(q)|^2 = 4 sin^2(qa/2) が使用されます。

パラメータ:

• q -- np.ndarray: 波数の配列

• a -- float: 格子定数

戻り値:

np.ndarray: C(q) の配列

class tiny_simulations.anharmonic.Params(a: float, M: float, Nk: int, k_cut: float, use_quantum: bool, eta_rel: float)

ベースクラス: object

計算に使用するパラメータを保持するデータクラス。

パラメータ:

• a -- float: 格子定数

• M -- float: 質量

• Nk -- int: k空間のグリッド分割数

• k_cut -- float: 赤外発散を防ぐためのkのカットオフ値

• use_quantum -- bool: 量子分布を使用するかどうか

• eta_rel -- float: バブル型シフトの主値積分に使用する相対的な正則化パラメータ

M: float

Nk: int

a: float

eta_rel: float

k_cut: float

use_quantum: bool

tiny_simulations.anharmonic.absC2(q: ndarray, a: float) → ndarray

位相因子の絶対値の2乗 |C(q)|^2 = 4 sin^2(qa/2) を計算する。

パラメータ:

• q -- np.ndarray: 波数の配列

• a -- float: 格子定数

戻り値:

np.ndarray: |C(q)|^2 の配列

tiny_simulations.anharmonic.bose_n(omega: ndarray, T: float, use_quantum: bool = True) → ndarray

ボース分布関数 n(ω) を計算する。

詳細説明: 量子的な場合は n(ω) = 1/(exp(β ħ ω) - 1) を計算します。 use_quantum=False の場合は、古典極限 n ~ kT/(ħ ω) を計算します。

パラメータ:

• omega -- np.ndarray: 振動数の配列

• T -- float: 温度 [K]

• use_quantum -- bool: 量子分布を使用するかどうか（デフォルトはTrue）

戻り値:

np.ndarray: 計算されたボース分布関数の配列

tiny_simulations.anharmonic.coth_halfbeta_hw(omega: ndarray, T: float, use_quantum: bool = True) → ndarray

双曲線余接関数 coth(βħω/2) の値を計算する。

詳細説明: 温度 T -> 0 の極限では 1 に近づきます。 古典極限では coth ~ 2kT/(ħω) として計算されます。

パラメータ:

• omega -- np.ndarray: 振動数の配列

• T -- float: 温度 [K]

• use_quantum -- bool: 量子分布を使用するかどうか（デフォルトはTrue）

戻り値:

np.ndarray: 計算結果の配列

tiny_simulations.anharmonic.delta_omega_cubic_bubble(q: float, kgrid: ndarray, a: float, M: float, K2: float, K3: float, T: float, use_quantum: bool, eta: float, k_cut: float) → float

3次 (K3) のバブル型シフト Δω_q^(3) を計算する。

詳細説明: 自己エネルギーの実部を用いて、Δω_q^(3) = Re Σ^(3)(q, ω_q) / (2 ω_q) を算出します。

パラメータ:

• q -- float: 対象の波数

• kgrid -- np.ndarray: 積分用のkグリッド配列

• a -- float: 格子定数

• M -- float: 質量

• K2 -- float: 調和力定数

• K3 -- float: 3次の非調和力定数

• T -- float: 温度 [K]

• use_quantum -- bool: 量子分布を使用するかどうか

• eta -- float: 正則化のための微小パラメータ

• k_cut -- float: 赤外発散を避けるためのkのカットオフ値

戻り値:

float: バブル型シフトによる振動数の補正量

tiny_simulations.anharmonic.delta_omega_quartic_tadpole(q: float, kgrid: ndarray, a: float, M: float, K2: float, K4: float, T: float, use_quantum: bool, k_cut: float) → float

4次 (K4) のオタマジャクシ型シフト Δω_q^(4) を計算する。

詳細説明: 実用的な1次元モデルの実装として以下を評価します:

Σ^(4)(q) ∝ K4 |C(q)|^2 * (1/N) Σ_k |C(k)|^2 <|u_k|^2>

ここで、<|u_k|^2> = ħ/(2 M ω_k) coth(βħω_k/2) です。 単純な係数選択を用いて Δω_q^(4) = ReΣ/(2ω_q) を返します。

パラメータ:

• q -- float: 対象の波数

• kgrid -- np.ndarray: 積分用のkグリッド配列

• a -- float: 格子定数

• M -- float: 質量

• K2 -- float: 調和力定数

• K4 -- float: 4次の非調和力定数

• T -- float: 温度 [K]

• use_quantum -- bool: 量子分布を使用するかどうか

• k_cut -- float: 赤外発散を避けるためのkのカットオフ値

戻り値:

float: オタマジャクシ型シフトによる振動数の補正量

tiny_simulations.anharmonic.k2_func(delta: float = 0.0, T: float | None = None) → float

調和力の定数K2を返す。

詳細説明: デフォルトでは定数を返しますが、ユーザーが任意の関数に置き換えることができます。 単位は [N/m]、または無次元の慣例に従う場合は一貫した単位となります。

パラメータ:

• delta -- float: 相対変位

• T -- float | None: 温度 [K]

戻り値:

float: 調和力定数 K2

tiny_simulations.anharmonic.k3_func(delta: float = 0.0, T: float | None = None) → float

3次の非調和力定数K3を返す。

詳細説明: デフォルトでは0を返しますが、バブル型シフトを有効にするためにゼロ以外の値に設定することができます。

パラメータ:

• delta -- float: 相対変位

• T -- float | None: 温度 [K]

戻り値:

float: 3次の非調和定数 K3 [N/m^2]

tiny_simulations.anharmonic.k4_func(delta: float = 0.0, T: float | None = None) → float

4次の非調和力定数K4を返す。

詳細説明: デフォルトの値を返しますが、オタマジャクシ型シフトの挙動を変更するために任意の関数に置き換え可能です。

パラメータ:

• delta -- float: 相対変位

• T -- float | None: 温度 [K]

戻り値:

float: 4次の非調和定数 K4 [N/m^3]

tiny_simulations.anharmonic.main() → None

スクリプトのエントリポイント。

詳細説明: パラメータの設定、kグリッドの作成、各温度に対する振動数シフトの計算を実行し、 結果をCSVファイルとして保存します。

戻り値:

None: 戻り値なし

tiny_simulations.anharmonic.make_kgrid(Nk: int, a: float) → ndarray

第1ブリルアンゾーンにおける一様なkグリッドを作成する。

詳細説明: 範囲は (-pi/a, pi/a] となります。

パラメータ:

• Nk -- int: k空間のグリッド分割数

• a -- float: 格子定数

戻り値:

np.ndarray: 作成されたkの配列

tiny_simulations.anharmonic.omega_T_for_q(q: float, T: float, kgrid: ndarray, p: Params, k2f: Callable[[float, float | None], float], k3f: Callable[[float, float | None], float], k4f: Callable[[float, float | None], float]) → Tuple[float, float, float, float]

指定された波数qと温度Tにおける、各シフト量と温度補正された振動数を計算する。

パラメータ:

• q -- float: 対象の波数

• T -- float: 温度 [K]

• kgrid -- np.ndarray: 積分用のkグリッド配列

• p -- Params: 計算パラメータのデータクラス

• k2f -- Callable[[float, float | None], float]: 調和力定数を返す関数

• k3f -- Callable[[float, float | None], float]: 3次の非調和力定数を返す関数

• k4f -- Callable[[float, float | None], float]: 4次の非調和力定数を返す関数

戻り値:

Tuple[float, float, float, float]: (調和振動数, 3次シフト量, 4次シフト量, 温度補正後の振動数) のタプル

tiny_simulations.anharmonic.omega_harmonic(q: ndarray, a: float, M: float, K2: float) → ndarray

最隣接相互作用モデルによる調和振動数を計算する。

詳細説明: 分散関係 ω_q^2 = (4 K2 / M) sin^2(qa/2) に基づいて計算します。

パラメータ:

• q -- np.ndarray: 波数の配列

• a -- float: 格子定数

• M -- float: 質量

• K2 -- float: 調和力定数

戻り値:

np.ndarray: 調和振動数の配列

tiny_simulations.anharmonic.parse_args() → Namespace

コマンドライン引数をパースする。

戻り値:

argparse.Namespace: パースされた引数オブジェクト

tiny_simulations.anharmonic.re_sigma_cubic_bubble(q: float, kgrid: ndarray, a: float, M: float, K2: float, K3: float, T: float, use_quantum: bool, eta: float, k_cut: float) → float

3次 (K3) のバブル型シフトによる自己エネルギーの実部 Re Σ^(3)(q, ω_q) を計算する。

詳細説明: バブル型の標準的な数式構造を用い、etaで正則化された主値積分 Re[1/(x + iη)] = x/(x^2 + η^2) を用いて評価します。

Σ^(3)(q, ω) = Σ_k |V3|^2 [ (n_k + n_{q-k} + 1)/(ω - ω_k - ω_{q-k} + i0)

• (n_k - n_{q-k})/(ω - ω_k + ω_{q-k} + i0)

• (n_{q-k} - n_k)/(ω + ω_k - ω_{q-k} + i0)

• (n_k + n_{q-k} + 1)/(ω + ω_k + ω_{q-k} + i0) ]

パラメータ:

• q -- float: 対象の波数

• kgrid -- np.ndarray: 積分用のkグリッド配列

• a -- float: 格子定数

• M -- float: 質量

• K2 -- float: 調和力定数

• K3 -- float: 3次の非調和力定数

• T -- float: 温度 [K]

• use_quantum -- bool: 量子分布を使用するかどうか

• eta -- float: 正則化のための微小パラメータ

• k_cut -- float: 赤外発散を避けるためのkのカットオフ値

戻り値:

float: 自己エネルギーの実部

tiny_simulations.anharmonic.v3_squared(q: float, k: ndarray, a: float, M: float, K3: float, omega_q: float, omega_k: ndarray, omega_qmk: ndarray) → ndarray

3次相互作用の要素 |V3(q,k,q-k)|^2 に比例する項を計算する。

詳細説明: V3 ∝ K3 * C(q)C(k)C(q-k) / sqrt( (2M)^3 ω_q ω_k ω_{q-k} ) となります。 複素位相を避けるため、|C|^2 = 4 sin^2(...) の形を使用します。

パラメータ:

• q -- float: 対象の波数

• k -- np.ndarray: 積分用の波数配列

• a -- float: 格子定数

• M -- float: 質量

• K3 -- float: 3次の非調和力定数

• omega_q -- float: 波数qの振動数

• omega_k -- np.ndarray: 波数kの振動数配列

• omega_qmk -- np.ndarray: 波数(q-k)の振動数配列

戻り値:

np.ndarray: |V3|^2 に対応する配列