[发明专利]基于改进三重互易边界元法的瞬态温度计算方法

权利要求书

1.一种基于改进三重互易边界元法的瞬态温度计算方法，包括以下步骤：

S1：基于计算时长T_T和最小时间步长Δt，采用倍数时刻方案确定计算时刻序列；

所述步骤S1具体步骤为：

S11：假设初始时刻为0，针对时刻t_i，采用基于Laplace变换的频域法进行求解，计算t_i时刻的频域参数s_ij为

β＝log2，整数N_L表示逆变换精度参数，一个时刻t_i对应2N_L个频域参数s_ij，每一个频域参数s_ij对应一个需要求解的频域方程，针对时刻t_i与t_l＝2t_i，根据上式可知，时刻t_i与t_l的频域参数s_ij与s_lj有一半重合，考虑到成2倍关系的时刻中t_l时刻可节约一半的计算量，提出倍数时刻方案，即取计算时刻

t_i＝2ⁱ·Δt i＝0,1,2,......；t_i≤T_T

式中Δt表示最小时间步长；

S2：基于Laplace变换得到瞬态热传导问题的频域方程，并基于分析时刻t和Laplace数值逆变换方法确定频域方程的频域参数；

所述步骤S2具体步骤为：

瞬态热传导问题的控制方程为：

式中a为热扩散系数，a＝k/(cρ)，k表示材料的热传导系数，c和ρ分别表示材料的比热和密度；Ω表示模型所占的几何区域，x为Ω中的任意点，称为场点，u(x,t)表示在t时刻点x处的温度，w(x,t)表示热源密度函数，对方程两边均进行Laplace变换，得到

式中分别表示Laplace变换后的温度函数、热源密度函数，u0(x)为初始温度分布函数，s为频域参数；

S3：基于修正亥姆霍兹基本解得到频域方程的边界-域积分方程；

S4：基于改进三重互易公式将边界-域积分方程包含的域积分转化为等效边界积分，得到纯边界积分方程；

所述步骤S4具体步骤如下：

S41：针对瞬态热传导问题中热源及初始温度分布引起边界积分方程中的域积分，基于三重互易法，记域积分函数为：采用满足如下方程的b1(x)对域函数b(x)进行近似，

b₁(x_nm)＝b(x_nm)

其中δ(x,x_m)为狄拉克δ函数，x_nm为所有边界节点和域内插值点，nm＝1,2,3,…,N+M，其中N为边界节点的个数，M为域内插值点的个数，x_m表示域内插值点，m＝1,2,3…M，b₁(x)、b₂(x)、b₃(x_m)为待求的未知函数；

S42：对S41中的近似函数b₁(x)满足的方程，基于Laplace基本解及边界元法进行求解b₁(x)、b₂(x)、b₃(x_m)，即得到如下边界积分方程：

其中为Laplace方程基本解，为高阶基本解，满足对于θ＝1和2有y_m为域内源点，从域内插值点x_m中依次选取；

S43：对S42中的边界积分方程，对第一、第二个方程中的源点y依次选N个边界节点，并对第三方程中的源点y依次选M个域内插值点，将共得到2N+M个方程，并基于边界元法理论，得到如下2N+M维线性代数方程组

G₁D₁+H₂B₂-G₂D₂-G_2DB_3D＝H₁B₁

H₁B₂-G₁D₂-G_1DB_3D＝0

G₁、G₂、H₁、H₂为方程系数矩阵，为N维方阵，其第行ψ¹列元素分别为：

其中为第ψ¹个边界节点的形函数；系数矩阵为M×N维矩阵，其第行ψ²列元素为：

G_1D、G_2D为N×M维系数矩阵，其第行ψ³列元素分别为：

为M维系数方矩，其第行ψ⁴列元素分别为：

其中B₁、B₂、D₁、D₂为N维列向量，分别为函数b₁(x)、b₂(x)、d₁(x)、d₂(x)在N个边界节点上的取值，其中B₁、B_1D为已知；M维列向量B_1D、B_3D为函数b₁(x)和b₃(x)在M个域内插值点上的取值；

令B₂＝0，得到求解总方程组如下

S44：对S43中2N+M维的线性代数方程组，基于自由度凝集方法进行降维求解，先记右端向量为：

将方程组拆散重组，进行自由度凝集，得低维矩阵表达形式：

进一步推导另外两个未知向量的表达式，并简记为：

B_3D＝TB_1D+Τ¹B₁

D₂＝Τ²B_3D

D₁＝Τ³B₁+Τ⁴D₂+Τ⁵B_3D

其中

通过上述降维的表达式，求解得到b₁(x)高阶函数d₁(x)、d₂(x)、b₃(x)在边界节点及域内插值点上的取值向量D₁、D₂、B_3D；

S5：基于边界元法理论表面网格与域内插值点，求解频域边界积分方程；

S6：对频域边界积分方程的解进行Laplace逆变换，求解时域温度分布，得到计算时刻序列中所有时刻的温度值；

S7：判断所得所有温度值中相邻时刻的最大温度变化值是否小于Δu_max，Δu_max为相邻时刻间允许最大温度变化值，若否，则返回步骤S1，基于倍数时刻方案进行新的计算时刻插入，重复S2-S6步骤计算新插入时刻的温度值，直至满足所有相邻时刻的最大温度变化值小于Δu_max。