以下は、Pythonスクリプト benzene.py を解析し、Sphinx (MyST) で問題なくビルド可能なMarkdownドキュメントです。

# ``benzene.py``: ベンゼン分子のHückel近似による電子構造計算と状態密度可視化

## 概要

``benzene.py`` は、ベンゼン分子（6炭素環）の電子構造をHückel近似法を用いて解析し、その結果を可視化するPythonスクリプトです。このスクリプトは、Hückel行列を定義し、その固有値問題を解くことで分子軌道のエネルギー準位と波動関数を計算します。計算されたエネルギー準位はソートされ、その縮重度が特定されます。その後、ガウス型ブロードニングを適用して状態密度（DOS）を計算します。最終的に、``matplotlib`` ライブラリを使用して、計算された電子準位と状態密度の両方を一つの図にプロットし、ベンゼン分子の電子状態を視覚的に表現します。

## 特徴

*   Hückel近似によるベンゼン分子の電子構造計算。
*   分子軌道のエネルギー準位（固有値）の計算。
*   エネルギー準位の縮重度の特定。
*   ガウス型ブロードニングを用いた状態密度（DOS）の計算。
*   電子準位とDOSを``matplotlib``で統合的に可視化。
*   Pythonの標準ライブラリと科学計算ライブラリ（``numpy``, ``scipy``, ``matplotlib``）のみを使用。

## 使い方

### 前提条件

このスクリプトを実行するには、Python 3と以下のライブラリが必要です。

*   ``numpy``
*   ``matplotlib``
*   ``scipy``

### インストール

必要なライブラリは ``pip`` を使用してインストールできます。

```bash
pip install numpy matplotlib scipy

実行方法

スクリプトはコマンドラインから直接実行できます。

python benzene.py

実行例

上記コマンドを実行すると、ターミナルにHückel行列と計算されたエネルギー準位が出力され、同時にベンゼン分子の電子準位と状態密度を示すグラフウィンドウが表示されます。

コンソール出力例

Hückel行列（ベンゼン）:
array([[0, 1, 0, 0, 0, 1],
       [1, 0, 1, 0, 0, 0],
       [0, 1, 0, 1, 0, 0],
       [0, 0, 1, 0, 1, 0],
       [0, 0, 0, 1, 0, 1],
       [1, 0, 0, 0, 1, 0]])

エネルギー準位（ソート後）:
[-2. -1. -1.  1.  1.  2.]

グラフ出力例

実行すると、以下のようなグラフが表示されるmatplotlibウィンドウが開きます。

左側のサブプロット: Hückel準位がエネルギー軸に沿って水平線で表示されます。
右側のサブプロット: 状態密度（DOS）がエネルギーに対するプロットとして表示され、ガウスブロードニングが適用されています。

出力

標準出力: Hückel行列の定義、計算されたエネルギー準位（ソート済み）。
グラフウィンドウ: matplotlibによって生成された、ベンゼン分子のHückel準位と状態密度を示す図。この図はアプリケーションが閉じるまで表示されます。

コード解析

グローバルパラメータ

スクリプトの冒頭で、Hückel近似計算と状態密度計算に関連するグローバルパラメータが定義されています。

H_BENZENE: numpy.ndarray ベンゼン分子（6炭素環）のHückel行列です。alpha = 0、beta = 1 と仮定した隣接行列形式で表現されます。これは6x6の対称行列です。

H_BENZENE = np.array([[0, 1, 0, 0, 0, 1],
                      [1, 0, 1, 0, 0, 0],
                      [0, 1, 0, 1, 0, 0],
                      [0, 0, 1, 0, 1, 0],
                      [0, 0, 0, 1, 0, 1],
                      [1, 0, 0, 0, 1, 0]])

NUM_POINTS: int 状態密度（DOS）を計算する際のエネルギー範囲の分割点数です。DOS曲線の滑らかさに影響します。デフォルト値は 1000 です。
SIGMA: float DOS計算におけるガウス型ブロードニングの幅（標準偏差）です。この値が大きいほど、DOS曲線は滑らかになります。デフォルト値は 0.1 です。

関数

main() このスクリプトの主要な処理を実行する関数です。Hückel行列の対角化、エネルギー準位の計算、縮重度の特定、状態密度（DOS）の計算、そしてその結果のプロットを行います。この関数は、スクリプトが直接実行された場合にのみ呼び出されます (if __name__ == "__main__": ブロック内)。

主要処理フロー

main() 関数内の主要な処理フローは以下の通りです。

Hückel行列の表示: グローバル変数 H_BENZENE の内容が pprint を用いて標準出力に表示されます。
固有値・固有ベクトルの計算: numpy.linalg.eig() を使用して、H_BENZENE の固有値と固有ベクトルを計算します。固有値は分子軌道のエネルギー準位に対応します。
エネルギー準位のソートと表示: 計算された固有値の実数部のみを取り出し、numpy.sort() を用いて昇順にソートします。ソートされたエネルギー準位は標準出力に表示されます。
縮重度の計算: 浮動小数点数の微小な差を考慮し、numpy.round() でエネルギー準位を小数点以下8桁で丸めます。その後、collections.Counter() を使用して、各エネルギー準位の出現回数（縮重度）を計算します。
状態密度（DOS）の計算:
- エネルギー範囲は、計算された最小および最大のエネルギー準位からそれぞれ 1 離れた範囲を NUM_POINTS で分割して定義されます。
- 定義されたエネルギー範囲に対し、各エネルギー準位とその縮重度に基づいてガウス関数を適用し、DOSを加算していきます。これにより、各離散的な準位がガウス関数でブロードニングされます。
- 計算されたDOSは、scipy.ndimage.gaussian_filter1d() を用いてさらに平滑化されます。平滑化の sigma は、定義された SIGMA を基にエネルギー範囲に合わせてスケールされます。
結果のプロット: matplotlib.pyplot を使用して、計算結果を可視化します。
- matplotlib.gridspec.GridSpec を用いて、1行2列のサブプロットレイアウトを作成します。左側のサブプロットの幅は右側のサブプロットの3分の1に設定されます。
- 電子準位図 (左側のサブプロット): ax0 として作成され、各エネルギー準位が水平線でプロットされます。Y軸はエネルギー（beta単位）を示します。
- 状態密度（DOS）図 (右側のサブプロット): ax1 として作成され、DOS曲線が線グラフでプロットされ、その下側が灰色で塗りつぶされます。X軸はDOS（任意単位）を示し、Y軸は左側のサブプロットと共有されます。
- plt.tight_layout() でプロット領域が最適化され、plt.show() で作成された図が表示されます。

注意事項

このスクリプトはベンゼン分子に特化しており、H_BENZENE 行列が固定されています。他の分子の計算には、Hückel行列の変更が必要です。
Hückel近似では、alpha = 0、beta = 1 の単位系が用いられています。プロットのY軸の単位もこれを反映しています。
縮重度の計算は浮動小数点数の丸め（小数点以下8桁）に依存します。計算精度によっては異なる結果が生じる可能性があります。
グラフはプログラムの実行中に表示され、手動で閉じる必要があります。