Untitled

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.animation import FuncAnimation
from mpl_toolkits.mplot3d import Axes3D
import matplotlib as mpl

# 设置中文字体
plt.rcParams['font.sans-serif'] = ['SimHei']
plt.rcParams['axes.unicode_minus'] = False

class HeatConduction3D:
    def __init__(self, nx=20, ny=20, nz=20, dx=0.01, dy=0.01, dz=0.01, 
                 k=0.5, density=1000, specific_heat=4186, 
                 T_hot=100, T_init=20):
        self.nx, self.ny, self.nz = nx, ny, nz
        self.dx, self.dy, self.dz = dx, dy, dz
        
        # 材料物理参数
        self.k = k                  # 热导率 (W/(m·K))
        self.density = density      # 密度 (kg/m³)
        self.specific_heat = specific_heat  # 比热容 (J/(kg·K))
        
        # 计算热扩散系数
        self.alpha = k / (density * specific_heat)
        
        # 根据热导率调整时间步长
        dx_min = min(dx, dy, dz)
        if k > 100:
            self.dt = 0.15 * dx_min**2 / (6 * self.alpha)  # 高导热材料
        elif k > 10:
            self.dt = 0.25 * dx_min**2 / (6 * self.alpha)  # 中等导热材料
        else:
            self.dt = 0.35 * dx_min**2 / (6 * self.alpha)  # 低导热材料
        
        self.T_hot = T_hot
        self.T_init = T_init
        self.current_step = 0
        
        # 初始化温度场
        self.T = np.ones((nx, ny, nz)) * T_init
        self.T[0, :, :] = T_hot  # 左面加热
        
        # 创建3D图形
        self.fig = plt.figure(figsize=(8, 6))
        self.ax = self.fig.add_subplot(111, projection='3d')
        
        # 预先创建体素网格
        self.voxel_grid = np.ones((nx-1, ny-1, nz-1), dtype=bool)
        
        # 设置初始视角
        self.ax.view_init(elev=20, azim=45)
        
        # 创建颜色映射
        self.cmap = plt.cm.plasma
        self.norm = mpl.colors.Normalize(vmin=self.T_init, vmax=self.T_hot)
        
        # 添加实际时间计数器
        self.real_time = 0.0  # 实际时间（秒）

    def step(self):
        """使用向量化运算计算温度场演化"""
        T_new = np.copy(self.T)
        T_new[1:-1, 1:-1, 1:-1] += self.alpha * self.dt * (
            (self.T[2:, 1:-1, 1:-1] + self.T[:-2, 1:-1, 1:-1] +
             self.T[1:-1, 2:, 1:-1] + self.T[1:-1, :-2, 1:-1] +
             self.T[1:-1, 1:-1, 2:] + self.T[1:-1, 1:-1, :-2] -
             6 * self.T[1:-1, 1:-1, 1:-1]) / (self.dx**2)
        )
        
        # 应用边界条件
        T_new[0, :, :] = self.T_hot
        T_new[-1, :, :] = T_new[-2, :, :]
        T_new[:, 0, :] = T_new[:, 1, :]
        T_new[:, -1, :] = T_new[:, -2, :]
        T_new[:, :, 0] = T_new[:, :, 1]
        T_new[:, :, -1] = T_new[:, :, -2]
        
        self.T = T_new
        self.real_time += self.dt  # 更新实际时间
        self.current_step += 1

    def update(self, frame):
        """更新动画帧"""
        # 增加每帧的计算步数
        for _ in range(5):  # 从2增加到5
            self.step()
        
        # 清除之前的图像
        self.ax.clear()
        
        # 创建颜色数组
        colors = self.cmap(self.norm(self.T[:-1, :-1, :-1]))
        
        # 绘制体素（降低透明度以提高性能）
        self.ax.voxels(self.voxel_grid,
                      facecolors=colors,
                      alpha=0.6)
        
        # 设置坐标轴范围和标签
        self.ax.set_xlim(0, self.nx-1)
        self.ax.set_ylim(0, self.ny-1)
        self.ax.set_zlim(0, self.nz-1)
        self.ax.set_xlabel('X 方向')
        self.ax.set_ylabel('Y 方向')
        self.ax.set_zlabel('Z 方向')
        self.ax.set_title(f'三维温度分布 (时间步：{self.current_step})\n热导率：{self.k} W/(m·K)')
        
        # 添加颜色条
        if not hasattr(self, 'cbar'):
            sm = plt.cm.ScalarMappable(cmap=self.cmap, norm=self.norm)
            self.cbar = self.fig.colorbar(sm, ax=self.ax, label='温度 (°C)')
        
        return []

def main():
    # 创建不同材料的模拟实例
    # 示例材料的热导率（W/(m·K)）：
    # 铜: 385
    # 铝: 237
    # 铁: 80.2
    # 不锈钢: 16.2
    # 水: 0.6
    # 空气: 0.026
    
    # 选择材料的热导率
    k = float(input("请输入材料的热导率 (W/(m·K)): "))
    
    # 统一网格大小
    nx, ny, nz = 20, 20, 20  # 统一使用20x20x20的网格
    
    # 根据热导率仅调整计算参数
    if k > 100:
        frames = 50
        interval = 20
        steps_per_frame = 12
    elif k > 10:
        frames = 80
        interval = 25
        steps_per_frame = 8
    else:
        frames = 120
        interval = 30
        steps_per_frame = 4
    
    heat_model = HeatConduction3D(
        nx=nx,
        ny=ny,
        nz=nz,
        dx=0.005,
        dy=0.005,
        dz=0.005,
        k=k,
        density=1000,
        specific_heat=4186
    )
    
    # 修改custom_update函数中的视角设置
    def custom_update(frame):
        # 根据热导率使用不同的计算步数
        for _ in range(steps_per_frame):
            heat_model.step()
        
        # 清除之前的图像
        heat_model.ax.clear()
        
        # 创建颜色数组
        colors = heat_model.cmap(heat_model.norm(heat_model.T[:-1, :-1, :-1]))
        
        # 绘制体素
        heat_model.ax.voxels(heat_model.voxel_grid,
                          facecolors=colors,
                          alpha=0.6)
        
        # 设置坐标轴范围和标签
        heat_model.ax.set_xlim(0, heat_model.nx-1)
        heat_model.ax.set_ylim(0, heat_model.ny-1)
        heat_model.ax.set_zlim(0, heat_model.nz-1)
        heat_model.ax.set_xlabel('X 方向')
        heat_model.ax.set_ylabel('Y 方向')
        heat_model.ax.set_zlabel('Z 方向')
        
        # 更新标题，添加实际时间信息
        heat_model.ax.set_title(
            f'三维温度分布\n'
            f'热导率：{heat_model.k:.1f} W/(m·K)\n'
            f'模拟时间：{heat_model.real_time:.3f} 秒'
        )
        
        # 添加颜色条
        if not hasattr(heat_model, 'cbar'):
            sm = plt.cm.ScalarMappable(cmap=heat_model.cmap, norm=heat_model.norm)
            heat_model.cbar = heat_model.fig.colorbar(sm, ax=heat_model.ax, label='温度 (°C)')
        
        # 保持固定视角
        heat_model.ax.view_init(elev=20, azim=45)
        
        # 添加计时器文本（在图形右上角）
        plt.figtext(
            0.98, 0.98,
            f'计算步数：{heat_model.current_step}\n'
            f'每帧步数：{steps_per_frame}\n'
            f'时间步长：{heat_model.dt:.2e} 秒',
            ha='right', va='top',
            bbox=dict(facecolor='white', alpha=0.8, edgecolor='none')
        )
        
        return []
    
    # 创建动画
    anim = FuncAnimation(
        heat_model.fig, 
        custom_update,
        frames=frames,
        interval=interval,
        blit=False
    )
    
    # 启用鼠标交互
    plt.rcParams['axes.prop_cycle'] = plt.cycler('color', ['r'])
    heat_model.ax.mouse_init()
    
    plt.show()

if __name__ == "__main__":
    main()