q1_2新生成

org/chatgpt2/__init__.py

@@ -69,35 +69,45 @@ plt.grid(True)
 plt.show()
 
 # -----------------------------
-# 3. 净辐射功率计算（小问2）
+# 问题2
+# -----------------------------
+T_film = 420  # 薄膜温度 K
+T_sky = 280   # 天空温度 K
+alpha = 0.1   # 吸收可见光比例
+I_sun_total = 1000  # W/m²，太阳总辐射
+I_sun = np.ones_like(wl_fine) * I_sun_total / len(wl_fine)  # 均匀分布
+tau_atm = 0.9  # 大气透射率假设
+
+# -----------------------------
+# 4. Planck 黑体辐射函数
 # -----------------------------
-# 黑体辐射谱 (Planck)
 def planck_spectrum(wl, T):
     wl_m = wl * 1e-6  # μm → m
     return (2*h*c**2 / wl_m**5) / (np.exp(h*c/(wl_m*k*T)) - 1)
 
-T_film = 300  # K
+# -----------------------------
-T_sky = 280   # K
+# 5. 计算净辐射功率
-# 假设太阳吸收率 alpha = 0.1
+# -----------------------------
-alpha = 0.1
+P_net_list = []
-# 假设太阳总辐射 1000 W/m²
-I_sun_total = 1000
-I_sun = np.ones_like(wl_fine) * I_sun_total / len(wl_fine)
-# 假设大气透射率 tau = 0.9
-tau_atm = 0.9
-
-plt.figure(figsize=(8,5))
 for d in thicknesses:
     epsilon = emission_dict[d]
     P_emit = np.trapz(epsilon * planck_spectrum(wl_fine, T_film), wl_fine)
     P_sun = np.trapz(alpha * I_sun, wl_fine)
     P_atm = np.trapz(epsilon * planck_spectrum(wl_fine, T_sky) * tau_atm, wl_fine)
     P_net = P_emit - P_sun - P_atm
-    print(f"d={d} μm, P_net = {P_net:.2f} W/m²")
+    P_net_list.append(P_net)
-    plt.bar(d, P_net, width=0.3)
 
-plt.xlabel("Thickness (μm)")
+# -----------------------------
+# 6. 绘图：净辐射功率 vs 厚度
+# -----------------------------
+plt.figure(figsize=(8,5))
+plt.plot(thicknesses, P_net_list, marker='o', linestyle='-', color='blue', label='Net Radiative Cooling')
+for x, y in zip(thicknesses, P_net_list):
+    plt.text(x, y, f"{y:.1f}", ha='center', va='bottom')
+
+plt.xlabel("Film Thickness (μm)")
 plt.ylabel("Net Radiative Cooling Power (W/m²)")
-plt.title("PDMS Thin Film Net Radiative Cooling")
+plt.title("PDMS Thin Film Radiative Cooling Performance")
 plt.grid(True)
+plt.legend()
 plt.show()

org/chatgpt2/q1_2.py

@@ -0,0 +1,182 @@
+import numpy as np
+import matplotlib.pyplot as plt
+from scipy.interpolate import CubicSpline
+
+
+# -----------------------------
+# 1. 从data.txt读取分块格式数据（先wl+n，再wl+k）
+# -----------------------------
+def read_split_data(file_path):
+    """
+    解析分块数据格式：
+    第一部分：wl	n（多行数据）
+    第二部分：wl	k（多行数据）
+    返回：sorted_wl, n, k（波长已排序，保证n和k一一对应）
+    """
+    with open(file_path, 'r', encoding='utf-8') as f:
+        lines = [line.strip() for line in f if line.strip() and not line.startswith('#')]
+
+    # 分割n数据和k数据（以"wl	k"为分界）
+    split_idx = None
+    for i, line in enumerate(lines):
+        if line == 'wl	k':  # 找到k数据的表头
+            split_idx = i
+            break
+
+    # 提取n数据（表头后到split_idx前）
+    n_lines = lines[1:split_idx]  # 跳过"wl	n"表头
+    wl_n = []
+    n_list = []
+    for line in n_lines:
+        wl, n_val = line.split()
+        wl_n.append(float(wl))
+        n_list.append(float(n_val))
+
+    # 提取k数据（split_idx后）
+    k_lines = lines[split_idx + 1:]  # 跳过"wl	k"表头
+    wl_k = []
+    k_list = []
+    for line in k_lines:
+        wl, k_val = line.split()
+        wl_k.append(float(wl))
+        k_list.append(float(k_val))
+
+    # 转换为numpy数组
+    wl_n = np.array(wl_n)
+    n_list = np.array(n_list)
+    wl_k = np.array(wl_k)
+    k_list = np.array(k_list)
+
+    # 确保n和k的波长完全一致（否则插值会出错）
+    assert np.allclose(wl_n, wl_k), "n和k的波长列表不一致！"
+
+    # 排序（按波长递增，避免插值异常）
+    sorted_idx = np.argsort(wl_n)
+    sorted_wl = wl_n[sorted_idx]
+    sorted_n = n_list[sorted_idx]
+    sorted_k = k_list[sorted_idx]
+
+    return sorted_wl, sorted_n, sorted_k
+
+
+# 读取数据
+wl_all, n_all, k_all = read_split_data('/Users/spasolreisa/IdeaProjects/asiaMath/data.txt')
+
+# 三次样条插值（覆盖全波段，保证计算精度）
+cs_n = CubicSpline(wl_all, n_all)  # 折射率n的插值函数
+cs_k = CubicSpline(wl_all, k_all)  # 消光系数k的插值函数
+
+# 定义待研究的PDMS薄膜厚度（μm），可按需调整
+thicknesses = [0.5, 1.0, 1.5, 2.0]
+
+
+# -----------------------------
+# 2. 核心物理模型：薄膜发射率计算
+# -----------------------------
+def fresnel_reflectance(n1, k1, n2, k2):
+    """
+    菲涅尔反射率（垂直入射近似，考虑消光系数k）
+    输入：介质1的n/k，介质2的n/k
+    输出：垂直入射时的反射率R（0-1）
+    """
+    m1 = n1 + 1j * k1  # 介质1的复折射率
+    m2 = n2 + 1j * k2  # 介质2的复折射率
+    return np.abs((m1 - m2) / (m1 + m2)) ** 2
+
+
+def thin_film_emissivity(n_film, k_film, d, wl, n_air=1.0, k_air=0.0):
+    """
+    薄膜发射率计算（考虑多光束干涉和吸收）
+    输入：
+        n_film/k_film: 薄膜的折射率/消光系数
+        d: 薄膜厚度（μm）
+        wl: 波长（μm）
+        n_air/k_air: 空气的折射率/消光系数（默认n=1, k=0）
+    输出：
+        epsilon: 发射率（0-1）
+    """
+    # 1. 计算上下表面的菲涅尔反射率
+    R12 = fresnel_reflectance(n_air, k_air, n_film, k_film)  # 空气→薄膜
+    R23 = fresnel_reflectance(n_film, k_film, n_air, k_air)  # 薄膜→空气
+
+    # 2. 计算相位差和吸收衰减
+    delta = 2 * np.pi * n_film * d / wl  # 干涉相位差（无吸收时）
+    alpha = 4 * np.pi * k_film * d / wl  # 吸收导致的振幅衰减系数
+
+    # 3. 多光束干涉反射率（考虑吸收）
+    numerator = R12 + R23 * np.exp(-alpha) + 2 * np.sqrt(R12 * R23 * np.exp(-alpha)) * np.cos(2 * delta)
+    denominator = 1 + R12 * R23 * np.exp(-alpha) + 2 * np.sqrt(R12 * R23 * np.exp(-alpha)) * np.cos(2 * delta)
+    R_total = numerator / denominator
+
+    # 4. 多光束干涉透射率（考虑吸收）
+    T_total = (1 - R12) * (1 - R23) * np.exp(-alpha) / denominator
+
+    # 5. 基尔霍夫定律：ε = 1 - R - T（局部热平衡）
+    epsilon = 1 - R_total - T_total
+    return epsilon
+
+
+# -----------------------------
+# 3. 计算并绘制不同厚度的发射率光谱
+# -----------------------------
+# 定义计算波长范围（覆盖数据的有效波长区间，步长0.01μm保证平滑）
+wl_min = wl_all.min()
+wl_max = wl_all.max()
+wl_fine = np.linspace(wl_min, wl_max, int((wl_max - wl_min) / 0.01) + 1)
+
+# 创建绘图
+plt.figure(figsize=(12, 7))
+plt.rcParams['font.sans-serif'] = ['Arial']  # 统一字体
+
+# 存储不同厚度的发射率结果（用于后续分析）
+emission_dict = {}
+
+for d in thicknesses:
+    # 插值得到当前波长下的n和k
+    n_film = cs_n(wl_fine)
+    k_film = cs_k(wl_fine)
+
+    # 计算发射率
+    epsilon = thin_film_emissivity(n_film, k_film, d, wl_fine)
+    emission_dict[d] = epsilon
+
+    # 绘制光谱曲线
+    plt.plot(wl_fine, epsilon, linewidth=2, label=f'Thickness = {d} μm')
+
+# -----------------------------
+# 4. 图表美化与标注（突出辐射制冷关键波段）
+# -----------------------------
+# 标注大气透明窗口（8-13μm，辐射制冷核心波段）
+if wl_min <= 13 and wl_max >= 8:  # 仅当数据覆盖该波段时才标注
+    plt.axvspan(8, 13, alpha=0.15, color='red', label='Atmospheric Window (8-13 μm)')
+
+# 坐标轴与标题
+plt.xlabel('Wavelength (μm)', fontsize=14)
+plt.ylabel('Emissivity ε(λ)', fontsize=14)
+plt.title('PDMS Thin Film Spectral Emissivity for Different Thicknesses', fontsize=16, fontweight='bold')
+
+# 网格、图例与范围设置
+plt.grid(True, alpha=0.3, linestyle='--')
+plt.legend(fontsize=12, loc='best')
+plt.ylim(0, 1.05)  # 发射率范围0-1.05（留出余量）
+plt.xlim(wl_min, wl_max)
+
+# 保存图片（高分辨率）
+plt.tight_layout()
+plt.savefig('PDMS_emissivity_spectrum.png', dpi=300, bbox_inches='tight')
+plt.show()
+
+# -----------------------------
+# 5. 输出关键波段（8-13μm）的平均发射率（若数据覆盖该波段）
+# -----------------------------
+if wl_min <= 13 and wl_max >= 8:
+    print("=== 8-13 μm 大气窗口平均发射率 ===")
+    wl_window = np.linspace(8, 13, 500)  # 大气窗口内的波长点
+    for d in thicknesses:
+        n_window = cs_n(wl_window)
+        k_window = cs_k(wl_window)
+        epsilon_window = thin_film_emissivity(n_window, k_window, d, wl_window)
+        avg_epsilon = np.mean(epsilon_window)
+        print(f"厚度 {d} μm: 平均发射率 = {avg_epsilon:.4f}")
+else:
+    print("数据未覆盖8-13μm大气窗口，跳过该波段平均发射率计算。")