电子测t与仪器学报2 0 0 ( ) 增刊 用于传感器有限元法计算中的稀疏方程组解法 曲国明 庄庆德 陈书旺 ( i j 北科技大 学， 石家 庄 0 5 0 0 5 4 ) 【 摘要l 本文报告了在用有限元法计算传感器特性的研究中， 高阶线性稀疏方程组的几种解法的比 较， 就占 用内存容i t 、运算速度、编程的难易等方面，提出了一种较好的解决方法。 关镇词:传感鹅有限元稀琉方程 侄 宜1 l 引言 有限元方法是一种将连续问 题离散化. 近似地求解数学物理问 题的 方法。 利用该方法分析计 算热型传感器内墓板上的 温度分布， 对传感器性能参数的优化和传感器的设计提供了一种重要的 手段。 一般在计算机上用有限元法计算热型传感器墓板上稳定温度场的步蕊是: 先将一定形状的传 感器墓板划分成单元并编制节点号和单元号， 确定单元号与边界号的关系， 单元的导热系数和边 界条件以 及节点坐标等有关数据: 然后把这些数据输入计算机，由 计算机按有限元变分计算公式 自 动生成沮度刚度矩阵K和右端向 量P ，最后求解关于温度的N阶线性有限元方程组: K】 ， = P ( 1 ) 其中K为N x N的 对称正定刚度矩阵. N为区域内划分的节点总数; T为墓板上N个节点 温 度组成的列向f矩阵:P为己知的右端列向f矩阵。 显然， 要求得高精度的近似解. 必须采用十分辅细的有限单元划分， 从而大大增加了区 域内 的节点总数， 必然会遇到求解几百阶至几万阶， 甚至几十万阶的大型有限元方程组问题。 在复杂 结构的有限 元分析中，这往往占 去总计算I的四分之三以 上，同时总刚度矩阵 K的存储常占 去 大里的微机内存。 因 此如何选择高效的有限元方程组的解法， 及选择占 用内 存较少的存储温度刚 度矩阵的方法。 对用有限元法求解传感器的 温度分布有着十分重要的愈义， 亦是迫切需要解决的 问题。下面就这两个问硬分别加以 讨论. 2温度刚度矩阵的基本特征 为得到温度刚度矩阵的有效存储方法，必须注愈温度有限元方程组的如下墓本特征: ( ) 对称正定 性和半正定性 经典的基于 极小势能原理得到的有限元方程组， 其温度总刚度矩阵总是 对称正定性的。据此 我们只需存储其下三角阵中的元素.对角线上部的元素依对称原则即可得到. ( 2 ) 离度稀硫性 当 要提高计算结果的精度或遇到大型复杂组合结构有限元时， 有限 元方程组的阶数可达到几 万或数十万阶。 尽管在传感器墓板上有数万个节点， 但与某一特定节点相邻的节点数一般不超过 本论文得到河北省自 然科学基金 ( 6 0 0 2 4 9 ) 资助。 曲 国明:.徽授:庄庆恤:教授;陈书旺:讲师. 4 6 7 电 子 测 盆 与 仪 器 学 报 _ _ 一 一一一一一一一一一一一一一一一一一2 丝塑 到万分之几.故尽量少存储刚度矩阵的零元素，必将大大节省占 用计算机内 存空间。 ( 3 )非零元素分布的规则性 合理编排所划分单元的节点号， 可使非零元素在刚度矩阵内的排列具有很强的规则性。 一般 分布在主对角线附 近相当 狭小的一定宽 度内。 常见的分布有等带宽、 变带宽、 子结构等带宽和子 结构变带宽等. 根据温度总刚度矩阵及其非零元素分布特征， 可以构造出不同的存储方法，以 适应不同的求 解方法。 3有限元方程组的存储方法 ( 1 ) 等带宽 压 缩存 储 在有限元方程组 ( 1 ) 式中，温度刚度矩阵K为N阶带型矩阵，即K满足: K ij二 0x 0 ( 1 i+ l ) ( i - 1 -.5 i 5 i + 1 ) ( 2 ) 称为半带宽为I w = 2 w + 1 , ( 2 ) 式可表示为: K二 ，、 K，、 0 K ， K . . KN - I CI KN . N - 1 KN - 1 ,N KN . N ( 3 ) O 对温度刚度矩阵只需存储带区内的非零元素，这既可节省内 存又减少了许多物入工作t。用 一个N行I w列的矩阵A 来表示K 矩阵. 且A的前w + 1行和后w + l 行的元素均以左对齐顺序放 在A 中相应行. 而最右边的空余部分均填入零。即矩阵K 用人 来表示格式为: !胜州州|J a , 丁 a翻 啼 、叭、 K 1 K. . , K , ， ， Q K，、 ，、 。 。 ， ( 4) 一- ，托沁汁 翻:a-.l:奴: a - - a 时L K，、 K， 、 ， ， K “ 一 ，0 all:勘:阮:如 = ( 2 ) 变带宽下 半带按行一维压绷存储 该方法有效地利用了 温度刚度矩阵的 对称性， 只存储沮度刚度矩阵的一半， 同时利用了非零 元素的 稀疏性和其分布的变带宽特征。 具体方法是: 把总刚度矩阵的下三角元素按行， 从第一个 非零元素到该行对角元素止的所有元素存放在一维数组C中。另外，用一维数组I D ( i ) 依次记录 各 主对 角元 素 排在 数组C 中 的 序号. 对于 大 型 方 程组用 数 组N A ( i) 记 录各 行第 一 个非 零 元素 的 列 号。 在计算程 序中 可 依N A ( i) 推算I D ( i 因 为I D ( l ) 必然为1 而且I D ( i) = 1 D ( i- 1 ) + i 一 N A ( i ) 十1 . 按此法矩阵K中的 元素K Y ( i j ) 处 在数组C中的 第I D ( i) 一 i + j 的位1上， 第i 行的第一个非 零元素的列标号是 N A M= i 一 F I D ( i ) 一 I D ( i- 1 ) 1 +1 ( 5 ) 且当 jN A ( i )时，a i j二0 . (L3 下半带紧凑存储 此法充分 利用了 温度总刚度矩阵的 对称性和稀疏性。 其体做法是:按行序号把矩阵中每行从 4 6 8 _. 一 . . . . . 电子测t与仪器学报z o o ( ) 增刊 第一个非 零元素到对角 元素， 按列序号顺序放于 数组C中; 非 零元素的 信息 存放于另一 数组I D ( i) 中， 它按行顺序存放每行对角线元素在数组中的序号，以 及每行非零元素的列标号。 按此法存储 时刚度矩阵中的第i 行第一个非零元素的列标号是 N L=I 以 N + I D ( i 一 1 ) + 1 )( b ) 而当JN L时， a ; j = D . 4 有限元方程组的解法 求解有限元方程组的方法主要有两大类:即直接法和迭代法。 每一类又有多种不同的算法。 但不论采用何种算法都必须既要保证解的 精度又要尽量少占用计算机内存。 提高运算速度， 特别 考虑到有限 元方程组的阶数一般都很高， 因此， 在求解过程中必然伴随着有效的数据组织与管理: 算法的选择上必须考虑温度刚度矩阵的 存储结构， 以实现合理的数据调度和管理。 我们在P C机 上使用V H 语言. 对常见的一些算法， 就其占 用内存容A和运算速度及编程的难易程度进行了比 较。给出如下结果。 ( 1 )高 斯 一赛德尔 迭代法 考虑到A度刚度矩阵K具有主对角线占 优的特点，故其迭代格式为 T W O 0 P _ D j可 W O 一 艺 气 叮 阅 / K 、 (i = 议. . . 冈 ( 7) 其中T ;(m j)为 第M + 1 次 迭 代 近 似 解， T (m 为 第M次 迭 代 近 似 解. 初 值T ( o ) 既 可 取 零 亦 可 取任惫值。 该方法优点是程序简单， 系数矩阵 的内存少 ( 2 ) 。 缺点是由于迭代次数较多， 全选主元高斯消去法 K的一切零元素都不必存储， 也不参加运算。 故占用机器 致使求解时间过长，运算速度过慢。 该方法分为以下几步:( w取值从1 开始直到N + 1 ) ( a ) 全选主元。 即从系数矩阵的第w行、 第w列开始的子矩阵中， 选取绝对值最大的元素， 并将其交换到主元素位置上。 0(iJ ( b ) 归一化 ( c ) 消元 ( d ) 回代 K ， / K . ， K 川 屹j - K. K . j -K . j = w + 一 留 W 十1 N + 1 ) N +1 ) ( 8 ) T . = K . ， K . ， 一 K . , K ; . ; . r - K . ， + .( w=N - 1 ， 1 ) 此方法在温度刚度矩阵K采用在第3 节中 所述的压缩存储方法时. 所占 用计算机内 存比 迭 代法增加的并不多，但却大大减少了 运算时间. 提高了运算速度。 ( 3 ) L D L T 分 解 法 这是一种与变带宽下半带按行存储格式相适应的直接求解方法。 ( a ) 分 解计算:K二 L D 口 计算过程如下: L ij ( 9) ( 1 0 ) 0. 山 山n口 厂!1月.fL 卜 一|曰 ILQ rlesl几ee.L 电卜 其 4 6 9 - - 一 钟 ， ， 碑 ， ， 种 种门 . . 电子测t与仪器学报2 0 0 0 增刊 i =1 久刘 i =m m , 二 ， , i- 1 ( 1 1 ) 如如 械轼喊 L tr = K ,j I 峨夕 其中m ; 为第i 行第一个非零元素的列标号. ( b ) 正代和回代 用上 述 矩阵分 解 的 方 法。 使 方程 组K T = P 变成以 D L T I = P 可 等价于 求 解下 面 两个 三 角 型方 程组:L Y 二 P 和L T T = D - Y .计 算公 式分 别是: Y = P , 一 氛又 79 = Y I t4 ， 一 芝 从 1兀 ( i = 议 . ， N ) ( 1 2 ) ( i = N , N - 1 ， 二 ， 1 ) 用该解法求解大型稀疏方程组，既可高效地节省微机内 存. 亦可大大提高徽机运算速度。 5结语 以上， 分析了用有限元方法计算传感器特性时， 求解稀琉线形方程组的几种算法及特点。 在 保证运算精度的前提下， 尽A减少占 用徽机内 存和提高运算速度， 实现合理的数据调度与管理是 评价各种解法优劣的关键。 结果表明， 求解小型稀疏方程组用迭代法比 较方便， 而对于中型和大 型 稀疏 方 程组 用 压 缩 存储 法存 储 温度总 刚 度 矩 阵. 同 时 选用L D L T 分 解 法求 解效率 较高。 用上 述方 法实际 地 计算了 平板结构的 热型传感器， 计算结果 请参阅参考文献川。 参考文献 D 庄庆 德 ， 陈书 旺 ， 曲 国 明有限 元 法计 算 热型 传 感 器的 特性( 特 发表) 1 2 1 王 助 成 ， 邵 敏， 有 限 单 元 法基 本原 理与 教 值 方 法， 清华 大 学出 版 社. 1 9 8 8 口 孔样谁，有限单元法 在 传热学中的 应用. 科学出 版社出 版， 1 9 8 6 S o l u t i o n o f S p a r s e Ma t r i x E q u a t i o n s A p p l i e d i n C a l c u l a t i o n o f S e n s o r U s i n g F i n i t e E l e m e n t Me t h o d Q u G u o m i n g , Z h