第1章 绪论
1.1 数值分析方法
1.1.1 数值分析方法来源于工程实际
在科学技术及工程实际问题中,有许许多多的实际问题可以归结为微分方程(常微分方程或偏微分方程)及边值问题。通常能用解析法求出精确解的问题一般都是方程性质简单、求解域几何形状规则的问题,如图1-1所示悬臂梁弹性变形条件下各横截面上的位移、应变及应力计算,可以采用弹性力学的方法解出。但对于大多数问题,由于控制方程一般具有非线性性质,或由于求解区域几何形状比较复杂,则不能得出精确解。如图1-2所示具有复杂边界形状的物体承受简单载荷作用,即使在弹性变形条件下也难以求出精确解。
图1-1 承受集中载荷作用的带孔洞悬臂梁(解析法)
图1-2 具有复杂边界形状的物体承受集中载荷作用(数值解法)
科学家和工程师寻求另一种求解问题的方法——数值解法,这种方法一般来说分为两个阶段。首先将连续的求解域划分为有限个网格及节点,选取适当的途径将微分方程及其定解条件转化为网格节点上相应的代数方程组,即建立离散方程组;然后在计算机上求解离散方程组,得到节点上的解,进一步利用插值函数计算网格内部任一点的解,从而得到微分方程定解问题在整个求解域上的近似解。近几十年来随着计算机技术的飞速发展,数值分析方法已经成为解决工程实际问题的主要工具。
1.1.2 数值分析方法分类
(1)有限差分法
有限差分法(Finite Difference Method,简称FDM)是数值解法中最经典的方法。它是将求解域划分为差分网格,用有限个网格节点代替连续的求解域,并在网格节点上用差分代替微分,推导出含有离散点上有限个未知数的差分方程组。求差分方程组(即代数方程组)的解,就是微分方程定解问题的数值近似解,这是一种直接将微分问题变为代数问题的近似数值解法。
有限差分法特别适合求解建立于空间坐标系的流体流动问题,在流体力学领域,至今仍占据支配地位。当用于求解域几何形状复杂时,精度降低或求解困难。
(2)有限体积法
有限体积法(Finite Volume Method,FVM)又称为控制体积法(Control Volume Method,CVM)。其基本思路是:将计算区域划分为网格,并使网格点周围有一个互不重复的控制体积,将微分方程对每一个控制体积积分,从而得出一组离散方程,其中的未知量就是网格点上的场函数值。为了求出控制体积的积分,必须假定场函数在网格点之间的变化规律。
有限体积法是近年发展非常迅速的一种离散化方法,其特点是计算效率高,目前在流体力学领域得到广泛应用,大多数商用计算流体动力学(CFD)软件都采用这种方法。
(3)有限元法
有限元法(Finite Element Method,FEM)的出现是数值分析方法研究领域内重大突破性进展。其基本思想是将一个连续的求解域分成有限个适当形状的离散化单元,在单元上分片构造场函数,然后根据变分原理或加权余量法将问题的控制方程(微分方程)转化为所有单元上的有限元方程,将所有单元总体合成,形成总体有限元方程(代数方程组),代入边界条件求出各个节点上的场函数值,再进一步利用插值函数求出求解域任一点上的未知量。
有限元法具有广泛的适应性,特别适用于几何及物理条件比较复杂的非线性问题,而且便于实现程序的标准化。在固体力学领域,绝大多数都采用有限元法,而在流体力学领域由于计算效率相对较低,因此得到一定的应用。
就离散方法而言,有限体积法可视作有限元法和有限差分法的中间物。有限元法必须假定场函数在网格节点之间的变化规律(插值函数),并将其作为近似解。有限差分法只考虑网格点上场函数的数值而不考虑场函数在网格点之间如何变化。有限体积法计算时只寻求场函数的节点值,但在寻求控制体积的积分时,必须假定场函数在网格点之间的分布。在有限体积法中,插值函数值只用于计算控制体积的积分,得出离散方程后就可以忘掉插值函数;如果需要的话,可以对微分方程中的不同项采用不同的插值函数。
1.1.3 有限元法的发展历史
有限元法的基本思想由特朗(Courant)在1943年提出,他第一次尝试把定义在三角形区域上的分片连续函数和最小位能原理相结合来求解圣维南扭转问题——变分问题的里兹方法。工程上真正利用计算机采用有限元法解决工程实际问题的是特纳(Turner)和克拉夫(Clough)所做的工作。1956年他们在分析飞机结构时,第一次给出了用三角形单元求解平面应力问题的正确解答。三角形单元的单元特性由弹性理论方程确定,采用直接刚度法推导单元刚度阵。1960年,Clough进一步处理了平面弹性问题,并第一次提出了“有限单元法”的名称。1963年Melosh认识到,有限元法的数学基础是变分原理,是一种基于变分原理的分片Ritz法,这就奠定了有限元的数学理论基础。
早期的有限元法建立在虚功原理和最小势能原理基础上,随着认识的加深,各国学者建立了基于不同变分原理的有限元法。如基于最小势能原理及其修正形式的位移元、位移杂交元和非协调元等,基于最小余能原理及其修正形式的应力元、应力杂交元等,基于Hellinger-Reissner二场广义变分原理及其修正形式或Hu-Washizu三场广义变分原理及其修正形式的混合元、混合杂交元等。由于多变量有限元法的参数匹配以及稳定性和收敛性理论的复杂性,在工程应用中目前仍以位移为未知量的位移法为主。20世纪70年代初,有限元基本理论和方法已发展成熟,随后的研究致力于高精度单元、板壳单元、非线性问题的迭代求解方法、适用于新型材料的有限元方法、多尺度有限元法和多场耦合等问题的研究。现在有限元法的应用已从弹性力学平面问题扩展到空间问题、轴对称问题、板壳问题等,从静力平衡问题扩展到稳定问题、动力学问题、传热学及热力学问题等。分析的对象从弹性材料扩展到弹塑性、刚塑性、黏弹性、黏塑性和复合材料等,从固体力学扩展到流体力学、传热学、流固耦合、热力耦合、电磁场、多场耦合问题等。在工程分析中的作用已从分析和校核扩展到优化设计并和计算机辅助设计技术相结合。可以预计,随着现代力学、计算数学和计算机技术等学科的发展,有限元法及数值模拟技术必将在国民经济建设及科学技术发展中发挥更大的作用,其自身也将得到进一步的发展和完善。
在塑性成形数值模拟过程中,主要应用非线性有限元方法。Yamada等在1968年推导了小变形弹塑性本构关系矩阵;Hibbit等在1970年提出了基于Lagrange描述的大变形弹塑性有限元列式;Mcmeeking等提出了基于Euler描述的大变形弹塑性有限元列式;Argris、Wempner及Belytschko等最早提出了基于共旋格式的弹塑性大变形算法,并广泛应用于非线性有限元计算中。至今,大变形弹塑性有限元法仍在不断发展完善。在实际问题中,很多问题可忽略弹性变形而采用刚塑性本构关系。1973年,Lee和Kobayashi提出了以Lagrange乘子法引入不可压缩条件的刚塑性有限元法;Zienkiewicz等在1975年提出了采用罚函数法处理不可压缩条件的刚塑性有限元法;Mori等在1977年提出了可压缩刚塑性有限元法。相对于弹塑性有限元法,刚塑性有限元避开了几何非线性问题,在增量较大的情况下以小变形方法来处理大变形问题,建模简单且计算效率高,收到了工程界的欢迎。在20世纪90年代以前,非线性有限元的解法主要是静力隐式算法,加州大学Berkeley分校的学者Hughes等人对这种方法的进步做出了杰出贡献,在一定程度上解决了这种方法求解复杂非线性问题时出现的收敛性方面的困难。在计算方法上,Costantino在1964年发展了最早的显式有限元程序。显式积分方法不需要计算刚度矩阵及其分解运算,不存在收敛性问题。显式有限元的重大进展是由John Hallquist博士在Lawrence Livermore实验室做出的,他在1976年发布了Dyna程序,该程序在世界上各大学和实验室广泛应用,成为当前显式求解程序的基础。由于不存在收敛性问题,目前显式分析程序广泛应用于复杂问题的非线性有限元分析中。