如何在该场景下正确使用curve_fit拟合求解带隙Eg？

问题现象 #

• 拟合后params与初始猜测值完全一致
• 协方差矩阵cov为无穷大

核心问题分析 #

1. 未定义全局常量：代码中kB（玻尔兹曼常数）和c（前置系数）未定义，导致拟合算法无法正常计算目标函数。
2. 函数缩进错误：return I0的缩进层级错误，函数无法正确返回计算结果。
3. 数值缩放失衡：Eg的量级为1e-19（焦耳），而I0为1e-5，方程T³*exp(-Eg/(kB*T))的数值波动范围极大，拟合算法无法找到有效梯度，只能停留在初始值。

修正方案与代码实现 #

方案1：修复非线性拟合代码 #

补全常量定义、修正缩进，并将c作为拟合参数（若已知c可直接固定），同时估算合理的初始值帮助收敛：

import numpy as np
from scipy.optimize import curve_fit

# 定义玻尔兹曼常数（单位：J/K）
kB = 1.380649e-23

# 原始数据
T = np.array([93, 100, 110, 130, 140, 159, 179, 195, 210, 224, 241, 255, 265, 283, 293])
I0 = np.array([3.89071815e-05, 3.86379100e-05, 3.19573788e-05, 2.75543444e-05,
               3.00226067e-05, 3.59134864e-05, 3.92161213e-05, 4.37321918e-05,
               4.95912369e-05, 5.24583115e-05, 5.47069424e-05, 5.73620696e-05,
               5.80254459e-05, 7.23664434e-05])

# 修正后的拟合函数
def fit_func(T, Eg, c):
    T3 = T ** 3
    arg = -Eg / (kB * T)
    return c * T3 * np.exp(arg)

# 估算c的初始值，帮助拟合收敛
c_init = I0[0] / (T[0]**3 * np.exp(-4.239e-19/(kB*T[0])))
# 执行拟合
params, cov = curve_fit(fit_func, T, I0, p0=[4.239e-19, c_init])

print(f"拟合得到的Eg: {params[0]} J")
print(f"拟合得到的c: {params[1]}")

方案2：对数变换转线性拟合（更稳定） #

对原方程取对数，将非线性问题转化为线性问题，避免数值缩放问题：
原方程：I0 = c*T³*exp(-Eg/(kB*T))
取对数后：ln(I0/T³) = ln(c) - Eg/(kB*T)
令x=1/T，y=ln(I0/T³)，则转化为线性方程y = ax + b，其中a = -Eg/kB，b=ln(c)

import numpy as np
import matplotlib.pyplot as plt

kB = 1.380649e-23
T = np.array([93, 100, 110, 130, 140, 159, 179, 195, 210, 224, 241, 255, 265, 283, 293])
I0 = np.array([3.89071815e-05, 3.86379100e-05, 3.19573788e-05, 2.75543444e-05,
               3.00226067e-05, 3.59134864e-05, 3.92161213e-05, 4.37321918e-05,
               4.95912369e-05, 5.24583115e-05, 5.47069424e-05, 5.73620696e-05,
               5.80254459e-05, 7.23664434e-05])

# 数据变换
x = 1 / T
y = np.log(I0 / (T ** 3))

# 线性拟合
slope, intercept = np.polyfit(x, y, 1)
Eg_fit = -slope * kB
c_fit = np.exp(intercept)

print(f"线性拟合得到的Eg: {Eg_fit} J")
print(f"线性拟合得到的c: {c_fit}")

# 绘图验证拟合效果
plt.scatter(x, y, label='原始数据点')
plt.plot(x, slope*x + intercept, 'r-', label='拟合直线')
plt.xlabel('1/T (K⁻¹)')
plt.ylabel('ln(I0/T³)')
plt.legend()
plt.show()

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