crystal_convert_cell.py 技術ドキュメント

プログラムの動作

crystal_convert_cell.py は、指定された結晶構造の単位格子を、異なる基本ベクトル系に変換し、変換前後の単位格子と原子配置を3Dで可視化するPythonスクリプトです。

主な機能:

  • 結晶情報の取得: FCC (面心立方), BCC (体心立方), Hex (六方晶), Rhomb (菱面体晶) といった一般的な結晶構造の初期格子パラメータと原子サイト情報を、内部定義またはコマンドライン引数に基づいて設定します。

  • 単位格子の変換: 指定された変換行列 (conversion_mode) を用いて、元の単位格子の格子ベクトルと、その単位格子内の原子の分数座標を新しい単位格子の基準で計算します。

  • 計算結果の表示: 変換前後の格子パラメータ (\(a, b, c, \alpha, \beta, \gamma\))、格子ベクトル、計量テンソル、単位格子の体積、および原子サイトの分数座標を標準出力に詳細に表示します。

  • 3D可視化: Matplotlibライブラリを利用して、元の単位格子と変換後の単位格子を同じ3D空間上に描画します。これにより、両単位格子の幾何学的関係と原子配置の変化を視覚的に理解できます。元の単位格子は黒線と指定された色の原子で、変換後の単位格子は青線と灰色原子で描画されます。

解決する課題:

結晶学において、ある単位格子を別の単位格子（例: 面心立方格子をプリミティブセル）に変換する際、その幾何学的関係や原子配置がどのように変化するかを直感的に理解することは重要です。このプログラムは、抽象的な変換操作を具体的な3Dモデルとして提示することで、この理解を深めるのに役立ちます。特に、変換前後の単位格子の違いや原子の相対的な位置関係を視覚的に確認できるため、教育目的や研究の初期段階での概念理解に有効です。

原理

このプログラムは、結晶学における単位格子の定義と変換の基本的な数学的原理に基づいています。

  1. 単位格子の定義: 結晶構造は、3つの基本並進ベクトル \(\mathbf{a}_1, \mathbf{a}_2, \mathbf{a}_3\) によって定義されます。これらのベクトルは、格子定数 \(a, b, c\)（ベクトルの長さ）と、ベクトル間の角度 \(\alpha, \beta, \gamma\) から計算されます。

  2. 格子定数と格子ベクトルの相互変換:

    • 格子定数から格子ベクトルへの変換: tkcrystalbase モジュール内の cal_lattice_vectors 関数によって行われます。一般的に、\(\mathbf{a}_1\) を \(x\) 軸上に、\(\mathbf{a}_2\) を \(xy\) 平面上に配置し、以下のように計算されます。 $\( \mathbf{a}_1 = (a, 0, 0) \)\( \)\( \mathbf{a}_2 = (b \cos \gamma, b \sin \gamma, 0) \)\( \)\( \mathbf{a}_3 = (c \cos \beta, c \frac{\cos \alpha - \cos \beta \cos \gamma}{\sin \gamma}, c \frac{\sqrt{1 - \cos^2 \alpha - \cos^2 \beta - \cos^2 \gamma + 2 \cos \alpha \cos \beta \cos \gamma}}{\sin \gamma}) \)$

    • 格子ベクトルから格子定数への変換: tkcrystalbase モジュール内の cal_lattice_parameters_from_lattice_vectors 関数によって行われます。ベクトルの長さと内積から角度を計算します。 $\( a = |\mathbf{a}_1|, \quad b = |\mathbf{a}_2|, \quad c = |\mathbf{a}_3| \)\( \)\( \cos \alpha = \frac{\mathbf{a}_2 \cdot \mathbf{a}_3}{|\mathbf{a}_2| |\mathbf{a}_3|}, \quad \cos \beta = \frac{\mathbf{a}_3 \cdot \mathbf{a}_1}{|\mathbf{a}_3| |\mathbf{a}_1|}, \quad \cos \gamma = \frac{\mathbf{a}_1 \cdot \mathbf{a}_2}{|\mathbf{a}_1| |\mathbf{a}_2|} \)$

  3. 単位格子の変換: 新しい基本ベクトル \(\mathbf{a}'_1, \mathbf{a}'_2, \mathbf{a}'_3\) は、元の基本ベクトル \(\mathbf{a}_1, \mathbf{a}_2, \mathbf{a}_3\) の線形結合で表されます。この関係は変換行列 \(\mathbf{T} = (t_{ij})\) を用いて次のように記述されます。 $\( \begin{pmatrix} \mathbf{a}'_1 \\ \mathbf{a}'_2 \\ \mathbf{a}'_3 \end{pmatrix} = \mathbf{T} \begin{pmatrix} \mathbf{a}_1 \\ \mathbf{a}_2 \\ \mathbf{a}_3 \end{pmatrix} \)\( または、各ベクトルについて \)\( \mathbf{a}'_i = \sum_{j=1}^3 t_{ij} \mathbf{a}_j \)\( ここで、プログラム内の `tij` はこの変換行列 \)\mathbf{T}$ に対応します。

  4. 原子サイトの分数座標変換: 元の単位格子における原子の分数座標を \(\mathbf{x} = (x, y, z)\) とし、その実空間位置を \(\mathbf{r} = x\mathbf{a}_1 + y\mathbf{a}_2 + z\mathbf{a}_3\) とします。新しい単位格子における分数座標を \(\mathbf{x}' = (x', y', z')\) とすると、\(\mathbf{r} = x'\mathbf{a}'_1 + y'\mathbf{a}'_2 + z'\mathbf{a}'_3\) です。 このとき、新しい分数座標 \(\mathbf{x}'\) と元の分数座標 \(\mathbf{x}\) の間には以下の関係があります。 $\( \begin{pmatrix} x' \\ y' \\ z' \end{pmatrix} = (\mathbf{T}^{-1})^T \begin{pmatrix} x \\ y \\ z \end{pmatrix} \)\( ここで \)(\mathbf{T}^{-1})^T\( は変換行列 \)\mathbf{T}$ の逆行列の転置です。プログラム内の convert_fractional_coordinates_by_tij 関数は、この原理に基づいて原子の分数座標を変換します。tij を引数として受け取り、内部で逆行列の転置を計算して使用していると推測されます。

  5. その他の機能:

    • reduce01: 原子座標を \([0, 1)\) の範囲に正規化するために使用されます（例: x - floor(x)）。

    • add_site: 原子サイトのリストに新しいサイトを追加する際、既存のサイトとの距離が rmin より小さい場合は同一と見なし、重複を避けます。これは新しい単位格子内の原子を効率的に集計するために使われます。

    • cal_volume: 格子ベクトルから単位格子の体積を計算します。

    • cal_metrics: 格子定数から計量テンソル \(g_{ij} = \mathbf{a}_i \cdot \mathbf{a}_j\) を計算します。

必要な非標準ライブラリとインストール方法

crystal_convert_cell.py を実行するためには、以下の非標準ライブラリが必要です。

  1. numpy: 数値計算、特にベクトルや行列の操作に利用されます。

  2. matplotlib: 結晶構造の3D可視化（グラフ描画）に利用されます。

また、本プログラムが依存するカスタムモジュール tkcrystalbase.py が必要です。

インストール方法:

numpy と matplotlib は、Pythonのパッケージマネージャ pip を使ってインストールできます。

pip install numpy matplotlib

tkcrystalbase.py は、この crystal_convert_cell.py スクリプトと同じディレクトリに配置するか、Pythonのモジュール探索パス (PYTHONPATH) が通っているディレクトリに配置する必要があります。通常は、crystal_convert_cell.py と同じディレクトリに置いておくのが最も簡単です。

必要な入力ファイル

このプログラムは、直接的な入力ファイルを必要としません。すべての設定はスクリプト内の変数定義またはコマンドライン引数として与えられます。

コマンドライン引数による入力:

  • crystal_name: 変換対象の結晶の種類を指定します。

  • conversion_mode: 使用する変換行列のモードを指定します。

  • kRatom: 描画時に原子半径に乗じる係数（Float型）。

内部定義による初期値:

スクリプトの冒頭で、以下のパラメータが初期値として設定されています。コマンドライン引数が与えられない場合、これらの値が使用されます。

  • crystal_name: 'FCC'

  • conversion_mode: 'FCCPrim'

  • lattice_parameters: [ 5.62, 5.62, 5.62, 90.0, 90.0, 90.0] (FCCの例)

  • crystal_sites: [['Na', 'Na1', ..., np.array([0.0, 0.0, 0.0])], ...] (FCCのNa原子サイトの例)

  • kr: 1000.0

  • rmin: 0.1 (同一原子サイトと判断する距離)

  • nrange: [[-0.1, 1.1], [-0.1, 1.1], [-0.1, 1.1]] (描画範囲)

  • figsize: (12, 12) (Matplotlibの図のサイズ)

  • converted_atom_color: 'gray' (変換後原子の色)

生成される出力ファイル

crystal_convert_cell.py は、プログラム自身がファイルを出力することはありません。

  • 標準出力 (stdout): プログラムの実行中に、以下のような情報がターミナルに表示されます。

    • 指定された結晶名と変換モード

    • 元の単位格子の格子パラメータ、格子ベクトル、計量テンソル、体積、および構成原子のサイト情報

    • 変換後の単位格子の変換行列、新しい格子ベクトル、格子パラメータ、計量テンソル、体積、および構成原子のサイト情報

    • 描画されるすべての原子サイトのリスト

  • Matplotlib 3Dプロットウィンドウ: プログラムの最後に、Matplotlibによって生成された3Dグラフが表示されます。このウィンドウには、元の単位格子（黒線、元の原子色）と、変換された新しい単位格子（青線、灰色原子）が同時に描画されます。ユーザーはこのウィンドウから手動で画像を保存できます（例: PNG, JPEGなどの形式）。

コマンドラインでの使用例 (Usage)

crystal_convert_cell.py は、以下の形式でコマンドラインから実行します。

python crystal_convert_cell.py crystal_name conversion_mode kRatom

  • crystal_name: 変換する結晶の種類を指定します。以下のいずれかを指定できます。

    • 'FCC' (面心立方格子)

    • 'BCC' (体心立方格子)

    • 'Hex' (六方晶)

    • 'Rhomb' (菱面体晶)

  • conversion_mode: 適用する変換行列のモードを指定します。以下のいずれかを指定できます。

    • 'FCCPrim' (FCCからプリミティブセルへ)

    • 'BCCPrim' (BCCからプリミティブセルへ)

    • 'ACenterPrim' (A面心からプリミティブセルへ)

    • 'BCenterPrim' (B面心からプリミティブセルへ)

    • 'CCenterPrim' (C面心からプリミティブセルへ)

    • 'RhombHex' (菱面体晶から六方晶へ)

    • 'HexRhomb' (六方晶から菱面体晶へ)

    • 'HexOrtho' (六方晶から直方晶へ)

  • kRatom: 描画される原子の半径に乗じる係数です。浮動小数点数を指定します。値を大きくすると原子が大きく表示されます。

コマンドラインでの具体的な使用例

例1: FCC格子をプリミティブセルに変換する

この例では、面心立方格子 (FCC) をそのプリミティブセルに変換し、両方の単位格子と原子配置を可視化します。

コマンド:

python crystal_convert_cell.py FCC FCCPrim 1000.0

実行結果の説明:

  • 標準出力: FCC格子の初期情報（格子パラメータ、格子ベクトル、サイト）と、FCCPrim 変換行列 ([[0.5, 0.5, 0], [0, 0.5, 0.5], [0.5, 0, 0.5]]) が適用された後の新しいプリミティブセルの格子情報（格子パラメータ、格子ベクトル、サイト）が出力されます。変換後の単位格子の体積は、元のFCC単位格子の体積の \(1/4\) になります。これは、FCC単位格子が4つの化学式単位（またはプリミティブセル）を含むためです。

  • 3Dプロット: 元のFCC単位格子が、中心と各面心に原子を持つ立方体として表示されます（黒い線、赤い原子）。この中に、変換によって得られたプリミティブな菱面体状の単位格子が青い線と灰色の原子で表示されます。これにより、FCC格子がどのように複数のプリミティブセルで構成されているかを視覚的に確認できます。

例2: BCC格子をプリミティブセルに変換する

この例では、体心立方格子 (BCC) をそのプリミティブセルに変換します。

コマンド:

python crystal_convert_cell.py BCC BCCPrim 1000.0

実行結果の説明:

  • 標準出力: BCC格子の初期情報と、BCCPrim 変換行列 ([[-0.5, 0.5, 0.5], [0.5, -0.5, 0.5], [0.5, 0.5, -0.5]]) が適用された後の新しいプリミティブセルの格子情報が出力されます。変換後の単位格子の体積は、元のBCC単位格子の体積の \(1/2\) になります。これは、BCC単位格子が2つの化学式単位を含むためです。

  • 3Dプロット: 元のBCC単位格子が、角と中心に原子を持つ立方体として表示されます（黒い線、赤い原子）。変換されたプリミティブな菱面体状の単位格子が青い線と灰色の原子で表示され、BCCのプリミティブな繰り返し単位が可視化されます。

例3: 六方晶を直方晶に変換する

この例では、六方晶単位格子を直方晶単位格子に変換します。

コマンド:

python crystal_convert_cell.py Hex HexOrtho 1000.0

実行結果の説明:

  • 標準出力: 六方晶の格子情報と、HexOrtho 変換行列 ([[1.0, 0.0, 0.0], [1.0, 2.0, 0.0], [0.0, 0.0, 1.0]]) が適用された後の直方晶の格子情報が出力されます。この変換は、通常、六方晶系の格子をより大きな直方晶の超格子に変換する場合に用いられます。

  • 3Dプロット: 元の六方晶単位格子が、六角柱の形状で表示されます（黒い線、赤い原子）。変換された直方晶単位格子が青い線と灰色の原子で表示され、両者の幾何学的関係が示されます。