如何在while循环中实现delta_omega累积求和及计算绘图

需求说明 #

原代码通过while循环逐个计算delta_omega数组中每个频率对应的原子捕获通量f，现需修改为：每次计算使用delta_omega的累积求和值（第1次用第1个元素，第2次用前2个元素之和，依此类推），最终绘制累积和与对应f值的关系图。

修改后的完整代码 #

import numpy as np
import matplotlib.pyplot as plt  

# -------------------------- 请先定义以下常量（根据你的实验参数填写） --------------------------
m = 100  # delta_omega数组的长度
k = 1.0  # 波矢
v_cap = 100.0  # 初始速度
mu_eff = 1.0  # 有效磁矩
B = 0.1  # 磁场强度
hbar = 1.0545718e-34  # 约化普朗克常数
s0 = 1.0  # 饱和参数
L = 0.5  # 相互作用长度
m_Rb = 85.4678e-3 / 6.02214076e23  # 铷原子质量
z0 = 1.0  # 减速区长度
t_b = 0.01  # 某段时间参数
a_max = 10.0  # 最大加速度
T = 300.0  # 温度
k_b = 1.380649e-23  # 玻尔兹曼常数
vp = 10.0  # 最概然速度
# -------------------------------------------------------------------------------------------

delta_omega = np.linspace(-900*10**6, -100*10**6, m)  # Hz - 频率范围
cum_delta_omega = np.cumsum(delta_omega)  # 计算累积和数组
f_values = []  # 存储每次计算的f值
i = 0

while i < len(cum_delta_omega):
    # 使用累积和数组的第i个元素进行计算
    delta = cum_delta_omega[i] - (k * v_cap) + (mu_eff * B) / hbar
    p_ee = (s0 * L / 2) / (1 + s0 + (2 * delta / L)**2)  # 激发态布居数
    R = L * p_ee  # 散射率
    F = hbar * k * R  # 原子受到的散射力
    a = F / m_Rb  # 假设加速度恒定
    
    vf_slower = np.sqrt(v_cap**2 - (2 * a * z0))  # 减速区末端速度
    t_d = (v_cap - vf_slower) / a  # 减速区耗时
    
    # -------- 减速区之后 --------
    da = 0.1  # 减速区末端到MOT中心的距离(m)
    vf_MOT = np.sqrt(vf_slower**2 - (2 * a * da))  # MOT中心处粒子速度(m/s)
    t_a = da / vf_MOT  # 减速区后到MOT中心的耗时(s)
    
    r0 = 0.01  # MOT捕获半径
    vr_max = r0 / (t_b + t_d + t_a)  # 最大横向速度
    
    vz_max = np.sqrt(v_cap**2 + 2 * a_max * z0)  # 最大轴向速度(m/s)
    
    # -------- 捕获原子通量 --------
    P = 10**(4.312 - (4040 / T))  # 液相蒸气压（固相用4.857替换4.312）
    A = 5e-4  # 炉口面积
    n = P / (k_b * T)  # 原子数密度
    
    f_oven = ((n * A) / 4) * np.sqrt(2 / np.pi) * np.sqrt((2 * k_b * T) / m_Rb)
    f = f_oven * (1 - np.exp(-vr_max**2 / vp**2)) * (1 - np.exp(-vz_max**2 / vp**2))
    
    f_values.append(f)
    i += 1

# 绘制累积和与f值的关系图
plt.figure(figsize=(10,6))
plt.plot(cum_delta_omega, f_values, 'b-', linewidth=2)
plt.xlabel('delta_omega 累积和 (Hz)')
plt.ylabel('原子捕获通量 f')
plt.title('累积频率偏移与原子捕获通量的关系')
plt.grid(True)
plt.show()

关键修改点 #

• 新增cum_delta_omega = np.cumsum(delta_omega)：计算delta_omega的累积和数组，满足每次计算的输入要求。
• 新增f_values = []：用于存储每次循环计算得到的f值，避免原代码中f被覆盖的问题。
• 循环内使用cum_delta_omega[i]替代原delta_omega[i]：确保每次计算基于累积和值。
• 绘图部分改为plt.plot(cum_delta_omega, f_values)：绘制累积和与对应f值的关系曲线。
• 补充了常量定义的注释提示：原代码中未定义的实验参数需用户根据实际情况填写。

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