`vib_irreps.py`: 分子の振動既約表現計算スクリプト

概要

このスクリプト vib_irreps.py は、指定されたXYZファイルと点群シンボルに基づき、分子の振動モードの既約表現を計算します。

詳細説明

本スクリプトは、分子の対称性を表す点群（Schoenflies 表記、例: C2v, D3d, Td, Oh, C3h, C4h, Ch, I, Ih）と分子の構造情報を含むXYZファイルを引数として受け取ります。

コードのdocstringによれば、pymatgen が利用可能な環境では、ヘルマン・モーガン表記を内部で使用することがありますが、それ以外の場合や特定の点群（I, Ih, C3h など）については、既知のクラスサイズと一般位置の仮定に基づく「抽象モード」にフォールバックします。ただし、本スクリプトのコード内では pymatgen は直接インポートされていません。

計算は以下のステップで進行します。

まず、全自由度に対応する指標 (Γ_3N) を、固定された原子のカウント（サイトベース）または一般位置の仮定から構築します。
次に、分子の並進運動 (Γ_T) と回転運動 (Γ_R) の指標を差し引き、純粋な振動モードに対応する指標 (Γ_v) を導出します。
最終的に、この Γ_v を指標の直交性定理を用いて、各既約表現への分解係数（多重度）を決定します。

コマンドライン引数

スクリプトは以下のコマンドライン引数をサポートします。

--pg:
- 必須引数です。
- 目的: 計算対象の点群シンボル（Schoenflies 表記）を指定します。
- 例: Td, Oh, D6h, I
- サポートされている点群は、tkpointgroup ライブラリによって定義されています。
--xyz:
- 必須引数です。
- 目的: 分子の構造情報を含むXYZファイルのパスを指定します。
- 例: CH4.xyz, SF6.xyz, C6H6.xyz, C60.xyz
--tol:
- オプション引数です。
- 型: 浮動小数点数 (float)。
- デフォルト値: 1e-5
- 目的: 固定された原子を識別するための許容誤差（単位: オングストローム、Å）を指定します。

依存ライブラリ

このスクリプトは以下の非標準ライブラリに依存しています。

numpy: 数値計算、特に多次元配列操作に使用されます。
tkpointgroup: 点群に関する操作、特性表、指標計算などを行うカスタムライブラリです。コード内のコメントには「共通ライブラリ」と記載されています。

使用例

python vib_irreps.py --pg Td  --xyz CH4.xyz
python vib_irreps.py --pg Oh  --xyz SF6.xyz
python vib_irreps.py --pg D6h --xyz C6H6.xyz
python vib_irreps.py --pg I   --xyz C60.xyz

入力

XYZファイルの形式

load_coords_from_xyz() 関数によって読み込まれるXYZファイルは、以下の一般的な形式を想定しています。

ファイルの最初の行には原子の総数が含まれると予想されますが、必須ではありません。
2行目はコメント行と見なされます。
3行目以降は、各原子の要素名と3次元座標（x, y, z）が記述された行が続きます。

例:

5
methane molecule
C   0.000000   0.000000   0.000000
H   0.629130   0.629130   0.629130
H  -0.629130  -0.629130   0.629130
H  -0.629130   0.629130  -0.629130
H   0.629130  -0.629130  -0.629130

ファイル内の空行やコメント行、または座標情報を含まない行はスキップされます。
もし最初の行が整数でない場合（例えばコメント行で始まる場合）は、座標情報の読み込み開始位置が自動的に調整されます。
各行の要素名はスキップされ、最後の3つの数値がそれぞれx, y, z座標として解析されます。カンマ区切りの数値もスペース区切りに変換されて処理されます。

その他の入力

コマンドライン引数 --pg で指定される点群シンボル。
コマンドライン引数 --tol で指定される許容誤差。

出力

スクリプトは計算結果を標準出力に整形して表示します。

計算結果のヘッダー:
- 指定された点群シンボル、対応する Hermann-Mauguin 表記、および計算モード（abstract または normal）が出力されます。
- 例: Point group: Td (HM = F-43m) mode=normal
クラス情報:
- 点群の各対称クラスのラベルが一覧表示されます。
- 例: Classes : E 8C3 3C2 6S4 6Sd

指標 (キャラクター) の表示:

以下の各指標について、各対称クラスごとのキャラクター値が示されます。
- Γ_3N: 全自由度に対応する指標。
- Γ_T: 並進運動に対応する指標。
- Γ_R: 回転運動に対応する指標。
- Γ_v: 純粋な振動モードに対応する指標 (Γ_3N - Γ_T - Γ_R)。

例:

Γ_3N chars:    E:  15.0    8C3:   0.0    3C2:  -1.0    6S4:  -1.0    6Sd:   3.0
Γ_T  chars:    E:   3.0    8C3:   0.0    3C2:  -1.0    6S4:  -1.0    6Sd:   1.0
Γ_R  chars:    E:   3.0    8C3:   0.0    3C2:  -1.0    6S4:   1.0    6Sd:  -1.0
Γ_v  chars:    E:   9.0    8C3:   0.0    3C2:   1.0    6S4:  -1.0    6Sd:   3.0

振動既約表現の結果:
- Γ_v を各既約表現に分解した結果（多重度）が、pg.pretty_irreps() 関数によって整形された形式で表示されます。
- 例: Vibrational irreps: 1A1 + 1E + 2F2

エラー出力

以下の場合に、エラーメッセージを標準出力に出力し、スクリプトは終了します。

指定された点群がサポートされていない場合。
- 例: Error: point group 'XYZ' not supported. Supported: C1, Ci, Cs, C2, C2h, C2v, C3, C3h, C3v, ...
XYZファイルから原子座標が一つも解析できなかった場合。
- 例: Error: no coordinates parsed from XYZ.

関数リファレンス

`load_coords_from_xyz(path: str) -> numpy.ndarray`

概要: 指定されたXYZファイルから原子の座標を読み込みます。
詳細説明: ファイル内の空行やコメント行、または座標情報を含まない行をスキップし、数値データのみを解析して numpy.ndarray として返します。最初の数行はヘッダー（原子数、コメント行）であると想定し、適切な行から座標の読み込みを開始します。もし最初の行が整数でない場合（例えばコメント行で始まる場合）は、行の読み込み開始位置を自動的に調整します。
引数:
- path (str): 読み込むXYZファイルのパス。
戻り値:
- numpy.ndarray: 読み込まれた原子座標の配列 (N, 3)。ファイルが空の場合、空の配列を返します。

`main() -> None`

概要: スクリプトの主処理を実行し、分子の振動既約表現を計算して表示します。
詳細説明: コマンドライン引数から点群とXYZファイルパスを受け取り、分子の振動モードの既約表現を計算して結果を表示します。処理は、点群が「抽象モード」（I, Ih, C3h など）であるかどうかに応じて分岐します。
- 通常モードの場合: tkpointgroup から実際の対称操作行列 (ops) を取得し、それに基づいて Γ_3N, Γ_T, Γ_R の指標を計算・集約し、 Γ_v を導出して分解します。
- 抽象モードの場合: 事前に定義されたクラスサイズ (pg.ABSTRACT_CLASS_SIZES) を使用して、 Γ_3N （特に 'E' クラスで 3N 、他は 0）、 Γ_T, Γ_R の指標を計算し、 Γ_v を導出して分解します。計算された各指標のキャラクターと最終的な振動既約表現の結果を標準出力に整形して表示します。
引数:
- なし (コマンドライン引数から情報を取得します)。
戻り値:
- なし (計算結果を標準出力に表示します)。

その他の情報

スクリプトの実行が完了した後、ターミナルがすぐに閉じないように、main() 関数内で以下の行が呼び出されています。

input("\nPress ENTER to terminate>>\n")

この行により、ユーザーが ENTER キーを押すまでプログラムは待機状態となります。

vib_irreps.py: 分子の振動既約表現計算スクリプト

概要