# プログラム名 / Program name: perturbation_energy_1D_well.py
# 内容 / Purpose:
# 一次摂動論に基づいて、無限ポテンシャル井戸中の粒子に摂動ポテンシャル U(x) = alpha * x を加えたときの
# エネルギー補正（1次）を数値積分で計算する。
# This script numerically calculates the first-order energy correction due to a perturbation potential
# U(x) = alpha * x in a 1D infinite potential well using first-order perturbation theory.

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import quad

# --- 定数定義 / Define constants ---
hbar = 1.0    # プランク定数 (簡略化のために1に規格化) / Reduced Planck constant (normalized)
m = 1.0       # 質量 (規格化) / Mass (normalized)
L = 1.0       # ポテンシャル井戸の幅 / Width of the potential well
alpha = 1.0   # 線形摂動ポテンシャルの係数 / Coefficient for linear perturbation
n = 1         # 対象とする状態 (n=1: 基底状態) / Quantum number (n=1: ground state)

# --- 非摂動波動関数 / Unperturbed wavefunction ---
def psi_n(x, n, L):
    """非摂動の波動関数 / Unperturbed wavefunction"""
    return np.sqrt(2 / L) * np.sin(n * np.pi * x / L)

# --- 摂動ポテンシャル / Perturbation potential ---
def U(x):
    """摂動ポテンシャル: U(x) = alpha * x / Linear perturbation"""
    return alpha * x

# --- 一次エネルギー補正の積分関数 / Integrand for first-order correction ---
def integrand(x, n, L):
    """ψ_n(x)^2 * U(x) を返す / Return psi^2 * U"""
    return psi_n(x, n, L)**2 * U(x)

# --- 数値積分による一次補正の計算 / Compute first-order correction via numerical integration ---
E1_n, _ = quad(integrand, 0, L, args=(n, L))  # 積分区間: [0, L]

# --- 結果表示 / Show result ---
print(f"一次摂動によるエネルギー補正 E1({n}) = {E1_n:.5f} [規格化単位 / normalized units]")

# --- プロット: 波動関数と摂動ポテンシャル / Plot: wavefunction and perturbation potential ---
x_vals = np.linspace(0, L, 1000)  # xの範囲 / Range of x
psi_vals = psi_n(x_vals, n, L)   # 波動関数の値 / ψ(x)
U_vals = U(x_vals)               # ポテンシャルの値 / U(x)

plt.plot(x_vals, psi_vals**2, label=r"$|\psi_n(x)|^2$")  # 確率密度 / Probability density
plt.plot(x_vals, U_vals, label=r"$U(x) = \alpha x$")     # 摂動ポテンシャル / Perturbation
plt.title("Wavefunction and Perturbation Potential")
plt.xlabel("x")
plt.ylabel("Amplitude")
plt.legend()
plt.grid(True)
plt.show()

# プログラム名 / Program name: heterostructure_wavefunction_plot.py
# 内容: 有効質量近似による電子波動関数 (包絡関数 × ブロッホ関数) を可視化
# Purpose: Visualize the electron wavefunction in a heterostructure using effective mass approximation

import numpy as np
import matplotlib.pyplot as plt

# --- パラメータ定義 / Define parameters ---
L = 1.0                     # ポテンシャル井戸の幅 / Width of quantum well
x = np.linspace(-1.5*L, 1.5*L, 1000)

# --- 包絡関数: 井戸内の波動関数 (n=1) / Envelope function: particle in a box ---
def chi(x, L):
    return np.where(
        (x >= -L/2) & (x <= L/2),
        np.sqrt(2/L) * np.sin(np.pi * (x + L/2) / L),
        0.0
    )

# --- ブロッホ関数 (周期性を持つ正弦波) / Bloch function (rapid oscillation) ---
def bloch(x, k0=30*np.pi):
    return np.cos(k0 * x)

# --- 総波動関数: Psi = chi × psi / Total wavefunction ---
envelope = chi(x, L)
bloch_func = bloch(x)
psi = envelope * bloch_func

# --- プロット / Plot ---
plt.figure(figsize=(10, 6))
plt.plot(x, envelope, label="Envelope $\chi(x)$", color='blue', linestyle='--')
plt.plot(x, bloch_func, label="Bloch $\psi(x)$", color='red', alpha=0.4)
plt.plot(x, psi, label="Total $\Psi(x)$", color='black')
plt.title("Electron Wavefunction in Heterostructure (Effective Mass Approximation)")
plt.xlabel("Position x")
plt.ylabel("Amplitude")
plt.axvline(-L/2, color='gray', linestyle=':')
plt.axvline(L/2, color='gray', linestyle=':')
plt.grid(True)
plt.legend()
plt.tight_layout()
plt.show()