一种高效的递推消元算法.docxVIP

一种高效的递推消元算法 1 精细积分求解 许多科学和工程问题可以概括为：在微组中，可以概括为两个边值问题。这样的解决方案方法具有重要的理论和应用意义。在两点边值问题的解法中,差分类方法的应用最为普遍,具有简便易用的特点,但差分引入的离散误差使得精度与计算量之间始终存在着矛盾。打靶法是一种将边值问题转化为初值问题进行求解的方法,其优点是计算量小,但同样存在离散误差,且当求解区域较长时,容易导致问题的病态。 精细积分法由于具有高精度、高效率等优点,在各类数学物理方程的求解中已经获得了日益广泛的应用。文献研究了齐次微分方程组的两点边值问题,建立了一种可以避免病态的精细积分方法。文献利用“增维”思想将非齐次方程化为齐次方程,从而使精细积分方法对非齐次方程也能适用。精细积分法在求解边值问题时,由于不引入离散误差,因此具有更高的精度,并且对于病态问题具有较好的适应性。 本文在精细积分研究的基础上,给出了两点边值问题精细积分解法的另一种思路。利用相邻两点间的精细积分关系式建立代数方程组,将其代入边界条件后,可将系数矩阵化为块三对角形式。利用这一特点本文给出了一种求解的高效递推方法,其优点是能灵活地处理各类边界条件,也能方便地处理非齐次问题。此外,本文方法对于奇异摄动边值问题的求解也十分有效。 2 算法的基本思想 由于高阶常微分方程(组)一般可化为一阶常微分方程(组),故仅考虑一阶常微分方程组。 设非齐次边值问题为 其中:v={q,p}T为n维向量,n=nq+np;H为n阶系数矩阵;上标“?”表示对x的微商。 边界条件为 将区间[0,xf]均匀划分为m等份,步长τ=xf/m。记:xi=iτ;vi=v(xi),有 其中:为特解;T=exp(Hτ)为传递矩阵。文献给出了T的一种高精度算法,即将步长τ分为2N等份,令?t=τ/2N,则 由于?t非常小,故对exp(H?t)进行Taylor展开时可仅取前几项,如 其中 式(4)相当于执行如下语句操作 当循环结束时 按上述过程可给出通解在计算机上的“精确解”。特解v*i+1的计算可采用特解精细积分法,如果非齐次项属于不能做特解精细积分的情况,则可采用数值积分法。也可通过“增维”的方法将非齐次方程转化为齐次方程。 将式(3)写成矩阵形式,为 将边界条件代入,并将T采用分块记法 其中:T11为nq阶矩阵;T11为nq×np阶矩阵;T21为np×nq阶矩阵;T22为np阶矩阵。式(9)可化为 其中 BL和BR分别表示对角线上的第一个和最后一个矩阵,其形式与边界条件有关。由此得到边值问题式(1)~式(2)的“离散形式”的代数方程组。与差分类方法不同的是,精细积分“离散”过程没有引入离散误差,因而具有较高的精度。式(11)是一个块三对角的方程,通常此类方程组可采用块追赶法求解。 对式(11)可以采用一种高效解法。取m=2M+1,利用第2个方程消去第1个方程的第2项,用第m-1个方程消去第m个方程的第m-1项,用第2i,2(i+1)个方程消去第2i+1个方程的第2i,2(i+1)项(其中i=1,2,,(m-3)/2)。得 其中 式(12)的方程数和未知数比式(11)减少了2M-1个。类似的,重复上述的消元过程,得到消元的递推公式为 其中 经过M次递推消元,得到 解得 注意到u1(M)即为1u,求出u1(M)后便可确定0p。然后由式(3)或者按照消元的逆过程回代便可求出问题的全部解答。 在递推过程中A(k)的左半部分始终为零,而C(k)只有左下1/4子块非零,这种形式对于降低计算量十分有利。 比较块追赶法与本文递推方法的计算量,由于矩阵乘法和矩阵求逆的计算量远大于矩阵与向量乘法的计算量,仅考虑矩阵乘法与矩阵求逆的计算量,采用块追赶法,所需矩阵求逆的次数为2M,矩阵乘法次数为2M+1-1。而按式(12)~式(15),矩阵求逆次数仅为M+2,矩阵乘法次数为6M+2。由此可见当M较大时,本文递推方法较块追赶法能大幅度提高计算效率,特别当方程式(1)为齐次方程时,在前M次递推消元过程中方程组式(13)的自由项无需进行运算,这使运算更加简便。对于取m=2M,递推公式也可类似推出。 3 任意边界条件问题 对于边界条件更为一般的情况,假设左右边界均为混合边界,即v(0)和v(xf)中既有部分q已知,又有部分p已知。 设x=0时的边界条件为 x=xf处的边界条件为 其中N1+N2+N3+N4=n。将边界条件式(16)~式(17)代入式(9),可以得到与式(11)形式上完全相同的代数方程。所不同的是矩阵BL、BR和向量u1、um、y1、ym的形式有所改变,但它们的表达式十分容易确定,不再赘述,故上述消元方法可以适用于任意边界条件问题。 本文方法也可以处理衔接问题,以图1所示的弹性力学问题为例。对于衔接问题,可先在各材料区域内进行“精