求助：Python绘制含积分的非分离变量u-T函数图像的方法

问题背景 #

需要绘制化学势u与温度T（1到1000整数）的关系图，原方程为：
$$5 \times 10^{28} = \int_0^\infty \frac{\pi}{2} \cdot \frac{8m}{h^2} \cdot E^{0.5} \cdot \frac{1}{\exp\left(\frac{E-u}{kT}\right)+1} dE$$
无法直接整理为u=f(T)的形式，需通过数值方法对每个T求解满足方程的u值。

原代码的关键错误 #

• 循环逻辑缺失：T是列表，不能直接代入积分函数，需逐个遍历每个T值求解对应u
• 运算顺序错误：积分函数中指数部分(E-u)/k*T的运算优先级出错，正确应为(E-u)/(k*T)
• 目标函数定义错误：Func_Equal_Zero不能直接计算，需针对固定T定义关于u的函数，用于数值求解
• 优化逻辑错误：每个T对应唯一u，无需同时优化T和u，应针对单个T给出u的初始猜测并求解

修正后的代码 #

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

# 定义物理常数
m = 9.11e-31    # 电子质量（kg）
h = 6.63e-34    # 普朗克常数（J·s）
k = 1.38e-23    # 玻尔兹曼常数（J/K）
target = 5e28   # 方程右侧目标值

# 积分函数：E为积分变量，u为化学势，T为当前温度
def integrand(E, u, T):
    prefactor = (np.pi / 2) * (8 * m / (h ** 2))
    return prefactor * (E ** 0.5) / (np.exp((E - u) / (k * T)) + 1)

# 目标函数：计算积分值与目标值的差的平方，用于最小化求解u
def objective(u, T):
    integral, _ = quad(integrand, 0, np.inf, args=(u, T))
    return (integral - target) ** 2

# 初始化温度列表和化学势结果列表
T_list = np.arange(1, 1001)
u_list = []

# 遍历每个温度，求解对应化学势u
for T in T_list:
    # 初始猜测：低温下u接近费米能，这里用1e-19作为初始值（可按需调整）
    u_guess = 1e-19
    # 调用最小化函数求解u
    result = spo.minimize(objective, u_guess, args=(T,), method='Nelder-Mead')
    u_list.append(result.x[0])

# 绘制关系图
plt.figure(figsize=(8, 6))
plt.plot(T_list, u_list, label='Chemical Potential u(T)')
plt.xlabel('Temperature (K)')
plt.ylabel('Chemical Potential (J)')
plt.title('Chemical Potential vs Temperature')
plt.grid(True)
plt.legend()
plt.show()

代码说明 #

1. 积分函数integrand：明确传入E、u、T三个参数，修正了指数部分的运算顺序，确保物理公式正确
2. 目标函数objective：对固定T，计算积分值与目标值的差的平方，通过最小化该值找到满足原方程的u
3. 循环求解逻辑：逐个遍历每个温度，针对每个T调用spo.minimize求解对应u，初始猜测值可根据物理意义调整
4. 绘图逻辑：先收集所有T对应的u值，再统一绘制关系图，确保图像生成正确

内容的提问来源于stack exchange，提问作者Jonathan Locker