如何使用scipy.stats.maxwell从粒子速度数据拟合获取系统温度

1. 拟合返回参数说明 #

scipy中所有连续概率分布的fit方法默认返回三类参数：形状参数、位置参数loc、缩放参数scale。麦克斯韦分布没有可变的形状参数，因此调用maxwell.fit()只会返回两个值：

• loc：分布的水平偏移参数，物理意义上麦克斯韦速率分布的起始点为0，因此拟合时必须固定该参数为0
• scale：分布的缩放参数，是唯一和物理量直接相关的拟合参数

2. 与系统温度的对应关系 #

三维麦克斯韦速率分布的物理表达式为：

$f(v) = \sqrt{\frac{2}{\pi}} \cdot \left( \frac{m}{kT} \right)^{\frac{3}{2}} \cdot v^2 \cdot e^{-\frac{mv2}{2kT}}$
其中$m$为单个粒子质量，$k$为玻尔兹曼常数，$T$为热力学温度

和scipy内置的maxwell分布概率密度公式（loc=0时）对比，可以得到缩放参数与温度的转换关系：
$T = \frac{m \cdot scale^2}{k}$
各物理量单位要求：

• $m$：单个粒子质量，单位为kg
• $k$：玻尔兹曼常数，固定为1.380649e-23 J/K
• $scale$：拟合得到的缩放参数，单位和输入速率一致（通常为m/s）
• $T$：输出热力学温度，单位为开尔文（K）

3. 正确拟合代码示例 #

import numpy as np
from scipy.stats import maxwell
import matplotlib.pyplot as plt

# 你的粒子速率数据，必须为三维速度的模长，单位m/s
velocity_data = # 替换为你的实际数据数组

# 固定loc=0拟合，避免无意义的偏移参数
loc_fit, scale_fit = maxwell.fit(velocity_data, floc=0)

# 计算系统温度
k = 1.380649e-23
m = 5.3e-26 # 示例为氧分子质量，替换为你的粒子实际质量，单位kg
T = m * (scale_fit ** 2) / k
print(f"拟合得到系统温度：{T:.2f} K")

# 拟合结果可视化验证
plt.hist(velocity_data, bins=50, density=True, alpha=0.6, label="原始速率分布")
v_arr = np.linspace(0, velocity_data.max(), 1000)
plt.plot(v_arr, maxwell.pdf(v_arr, loc=loc_fit, scale=scale_fit), "r-", label="拟合麦克斯韦分布")
plt.xlabel("速率 (m/s)")
plt.ylabel("概率密度")
plt.legend()
plt.show()

注意事项 #

• 输入数据必须是三维速度的模长（速率），不能是单方向的速度分量，单分量服从正态分布而非麦克斯韦分布
• 拟合时必须加floc=0参数固定偏移量，否则拟合得到的偏移参数会导致结果完全不符合物理定义

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