有限元法-PPT教学课件

第3章 有限元法基础 3.6 有限元方程组的求解 用有限元 方法 快递客服问题件处理详细方法山木方法pdf计算方法pdf华与华方法下载八字理论方法下载 求解电磁场问 题 快递公司问题件快递公司问题件货款处理关于圆的周长面积重点题型关于解方程组的题及答案关于南海问题 ，在导出有限元方程组后，即归结为求解一个方程组的问题，方程组阶数等于未知节点参数的节点数。 γ＝常数时：线性方程组 γ≠常数时：非线性方程组 3.6.1 有限元线性方程组解法 高斯-若尔当解法 注意有限元方程组系数矩阵的特点： (1) 对称性 (2) 稀疏性与正定性 第3章 有限元法基础 下面主要介绍牛顿-拉夫逊迭代法。 3.6.2 有限元非线性方程组解法 线性化方法：如逐次线性化方法，牛顿-拉夫逊迭代法及改进型的牛顿-拉夫逊迭代法等。 通过求解目标函数的极小值点来获得非线性方程组的解：如最速下降法，共轭梯度法等。 逐次线性化方法特点：比较简便，有良好的收敛性，但收敛速度慢，适宜在求解方程的阶数不是特别高时采用； 牛顿-拉夫逊迭代法特点：收敛速度快，但在形成方程组时需要很大的计算量，并要求有很好的初值，否则迭代过程可能不收敛。为了克服这个缺点，但又能发挥其特点，产生了多种改进型的牛顿-拉夫逊迭代法。 1. 牛顿-拉夫逊迭代法 第3章 有限元法基础 对矢量磁位 A，使用直角坐标系， 一般认为 A 是近似解，对于真解必有一个不为零的余矢量 在此方程中，K 和 F 都是 A 的函数，如果 A(0) 是 的近似解，将F(A)用泰勒级数 表 关于同志近三年现实表现材料材料类招标技术评分表图表与交易pdf视力表打印pdf用图表说话 pdf 示，取一次项，则有 若用 J(A(0)) 表示 F(A) 的雅可比矩阵在 A(0) 处的值，则有 令此余矢量对于 A(1) 解时为0，有 第3章 有限元法基础 式中 求余矢量的雅可比矩阵，需要从一个单元的余矢量函数讨论起。 一个单元的余矢量为： 第3章 有限元法基础 单元雅可比矩阵为： 以 为例，看雅可比矩阵第 ii 项元素对单元 e 的贡献。K 是 A 的函数 ， ， 对A的偏导数均不为零，由 (3.66)式可知 因此有： 第3章 有限元法基础 下面求 。由以下各式： 第3章 有限元法基础 可以求出： 式中 所以 同理可求出： 单元雅可比矩阵元素统一表示为： （γ、B 对各单元均不相同） 求出每个单元的雅可比矩阵 J(e) 和余矢量 F (A(k)) 后，将所有单元的 分析 定性数据统计分析pdf销售业绩分析模板建筑结构震害分析销售进度分析表京东商城竞争战略分析 结果进行总体合成，得到： 则牛顿-拉夫逊迭代公式可以归结为下式： 第3章 有限元法基础 第3章 有限元法基础 （3.185) 在上式中，对于J(e) 中各元素，必须在求出 后才能进行计算，因此，由矢量磁位 A 求解时，由每次迭代解 A(K) 计算出单元磁通密度 B，由 B 值查磁化曲线得对应 H 值，再计算 ，磁化曲线的查取归结为数值计算。 对 B = f (H)的数学处理，一般包括插值法和拟合法，比较简单而有效的办法是线性插值和拉格朗日插值等。 当磁场磁化曲线采用逐段线性插值函数逼近，H 值在 HK+1与 HK之间时，有 第3章 有限元法基础 2. 改进型的牛顿-拉夫逊迭代法 牛顿-拉夫逊迭代法的优点： （1）收敛速度快，按平方律收敛，每经一次迭代有效位数基本上增加一倍； （2）自校正功能，即 A(K+1) 仅依赖于 F(A) 与 A(K)，前面迭代的舍入误差不会一步步传递下去。 其缺点： （1）每次迭代都要形成一次 J (K)，而 J (K) 的计算时间往往比增加迭代的次数所需要的时间多； （2）对初值要求较高，选择不当，常会引起振荡。 为了克服所存在的缺点，现已形成多种改进型牛顿-拉夫逊迭代法。 第3章 有限元法基础 (1) 修正的牛顿—拉夫逊迭代法 为了克服每次迭代都需要计算 J (K)，采取全部求解过程都用J (0)，即每次迭代不需要形成新的系数矩阵，保持其斜率不变。这将使迭代次数增加而影响收敛速度，但如初值选取较好，总体上将节省运算时间。 (2) 采用欠松弛因子的方法 欠松弛迭代的方法为 式中ω为收敛因子，0<ω<1；ΔA(K) 为每次迭代近似解的误差。 此外，还有牛顿-拉夫逊迭代法与最速下降法的组合求解法，加速的全牛顿法，采用阻尼因子等，具体方法请参阅文献［2］。 否 读入所需系数及数组 建立右端矢量 R 建立雅可比矩阵及右端剩余矢量 ΔR 计算迭代误差 计算γ的初值 计算各单元 b，c，Δ 是 开始 牛顿－拉夫逊迭代 迭代误差是否小于控制值？ 结束 非线性方程组计算框图： 根据边界条件修正方程组 计算各单元 A，B 解方程组得到各节点 A 计算各单元γ及 3.8 有限元素的自动剖分 采取自动剖分的必要性： 第3章 有限元法基础 在单元剖分中，为了尽可能的压缩存储，减小计算量，提高精度，必须注意以下问题： (1) 三角形各边不要相差太悬殊，避免出现尖锐的三角形； (2) 一个节点周围，不宜集中过多的三角形单元，为压缩存储创造条件； (3) 三角形单元内物理参数(γ或μ)变化连续，即媒质交界面应与单元的边界重合； (4) 精度要求不同的区域，元素的密度应不同； (5) 元素节点编排应规格化。 本节讨论平面内的几种自动剖分方法。 第3章 有限元法基础 3.8.1 直线内插法 适宜于对以直线段为边界的的场域进行自动剖分。只需给出 x, y方向上两端点的坐标，就可算出所有节点的坐标。 1. 确定 x 方向、y方向节点数及总节点数 x方向节点数 y方向节点数 节点总数 2. 确定各节点的坐标 x 和 y 方向节点坐标最小值为x1、y1，最大值为xm、ym，则节点坐标增量分别为 ， 第3章 有限元法基础 可推得第 Ni 列 Nj 行的第 Nk 个节点的编号为 式中 第 Nk 个节点的坐标值为 例3.1 设 x1=0, y1=0, xm=7cm, ym=5cm，取n=7, m=5,可以求出x方向节点数 Nx=7+1=8，y方向节点数Ny =5+1=6，总节点数 N0 =8×6=48。求第2列第3行节点数编号及坐标。 节点编号( Ni=2, Nj=3 )，按(3.192)式，有 Nk=(2-1)6+3=9 其坐标按(3.193)式为 （3.193) 3.确定三角形单元编号及各三角形单元3顶点(节点)按规定顺序的的总体编号 第3章 有限元法基础 对于图(a) 中类型的三角形单元，其单元编号 E 与行、列的关系为 E＝2(Ni-1)(Ny-1)+Nj 式中 Ni＝1，2，3…，(Nx-1); Nj ＝1，2，3…，(Ny-1)。 图(a)中第 E 个三角形单元的3节点，i, j, m的总体编号与行、列的关系分别为 同理，对图(b)类型的三角形单元有 第3章 有限元法基础 对于图(b) 中类型的三角形单元，其单元编号 E 与行、列的关系为 E＝2(Ni-1)(Ny-1)+Nj 式中 Ni＝1，2，3…，(Nx-1); Nj ＝1，2，3…，(Ny-1)。 3节点 i, j, m的总体编号与行、列的关系分别为 例3.2 (1) 求上例中1列3行（b）类型单元的编号 E 和3个节点 i, j, m 的总体节点编号。(2) 求出1列3行（a）类型单元的编号 E 的3个节点 i, j, m的总体节点编号。 当场域形状复杂，要求场域内单元密度分布不均匀，又希望密度是平滑过渡时，这种简单的直线内插法就不适用了，广泛应用的是等势剖分的方法。 第3章 有限元法基础 3.8.2 等势剖分 1962年，W.P. Crowley 1. 微分方程的推导 等势剖分是把网格线看成平面上的两个势函数的等势线，定义两个函数 (x,y) 及ψ(x,y) ，它们在整个计算区域内部满足拉普拉斯方程： (3.198) 第3章 有限元法基础 当边界 和ψs 给定后，可以解出场域内 和ψ的分布，求出等势线。 两族曲线的交点就是剖分所要确定的节点，但是得到的节点坐标是 和ψ 值，我们需要的是等势线交点在(x,y)平面上的坐标！ 若将 x, y看做( ,ψ)平面上的函数 x( ,ψ)及 y( ,ψ)，通过函数的变换，方程式(3.198)可化为所要求解的形式，简单步骤如下： 设 两端对 x 求导 第3章 有限元法基础 可解得 令 可得 (3.204) 第3章 有限元法基础 又因为 两端对 x 求导 将式(3.204)两式展开，与上式比较可得 (3.208) (3.209) (3.210) (3.211) 第3章 有限元法基础 将式(3.208)、 (3.209)代入 ，则有 其中 由式 可得 (3.213) (3.212) (3.214) 将式(3.213)中 展开为 第3章 有限元法基础 将式(3.203)、(3.208)、 (3.209)、 (3.210) 、(3.211) 、(3.213)、 (3.214)及 表达式代入式(3.212) ，化简可得 (3.215) 同理可以导出 (3.216) (3.217) 第3章 有限元法基础 其中 齐次方程式(3.216)、式(3.217)有唯一解的条件为 J≠0，因此有 若在( ,ψ)平面得到式(3.219)的两组数值解x( ,ψ)和y( ,ψ)，则对于定域内每一个确定的 、ψi ， 就有一对 xi 和 yi ，这就是要求的某两条网格线相交处的节点坐标。 (3.218) (3.219) 在( , ψ)平面上网格与节点数与(x,y)平面上是相等的，但在( ,ψ)平面上是等距的网格，在(x,y)平面上却是用相应边界确定的，各单元密度过渡是平滑的。但并不是等距网格。 第3章 有限元法基础 2. 等势剖分计算公式的推导 图3.19 网格与边界的划分 x y平面区域 逻辑网格 （1）边界值的确定 第3章 有限元法基础 （2）三角形网格的逻辑坐标 在空间坐标中用三角形网格剖分时， 在 ( ,ψ)坐标中，等 线和等ψ 线有多种取法，常用的有如图 3.20 所示3种。 图3.20 不同的逻辑网格方式 (a) (b) (c) （3）对等边三角形逻辑网格求(3.219)式的差分方程 第3章 有限元法基础 令等边三角形 , ，推导的关键是对三角形网格求出 ， ， ， ， ， 。 利用积分中值定理有 设函数 格林公式 同理，有 第3章 有限元法基础 对于二次积分，同样应用积分中值定理和格林公式，可得 现把 替换成 x( ,ψ) 或 y ( ,ψ)，可以计算出每个节点( l, m)的 ， ， 与 。 第3章 有限元法基础 图3.21 等边三角形节点编号 如图示，节点( l, m)由6个等边三角形包围，求 的回路积分是沿6个三角形的外边界积分，即 第3章 有限元法基础 同理可得 ， 与上面求法完全相同，即将 xi 换成 yi 以上式子对 y 全部成立。由此可进一步求出 α，β ，γ 。 为求出 ， ， ， ， ， ， 要利用图3.21中虚线网格，即对质心构成包围(l , m)节点的六边形进行积分。 第3章 有限元法基础 首先求出各三角形内平均的 和 。 利用上式 ，求出图3.21式中对应三角形的 。 第3章 有限元法基础 再沿图中虚线做环路积分，求出 对应图示各点 ( m, l )： 1’ [( m+1/3 ), ( l+1/3 )] I [ m, ( l+1/2 )] 2’ [( m-1/3 ), ( l+1/3 )] II [( m-1/2 ), ( l+1/2 )] 3’ [( m-2/3 ), ( l+1/3 )] III [(m-1/2), l ] 4’ [( m-1/3 ), ( l-1/3 )] IV [ m, ( l-1/2 )] 5’ [( m+1/3 ), ( l-2/3 )] V [( m+1/2 ), ( l-1/2 )] 6’ [( m+2/3 ), ( l-1/3 )] VI [( m+1/2 ), l ] 同理可得 第3章 有限元法基础 同理可得到 ， ， 类似的关系。 将它们代入式(3.219)： 整理后可得 上式整理后可得 第3章 有限元法基础 归纳成统一形式为： 式中 wi 是α，β ，和γ 的函数。上面得到的值是对等边三角形的情况。对于不同的逻辑网格，它们取值将不同。 (3.236) 3. 差分方程组的求解 第3章 有限元法基础 方程式(3.236)是非线性代数方程，n个内部节点有2n个坐标，有2n个方程，一般采用高斯-赛德尔选代。计算的步骤为： (1) 边界顶点坐标已知，其它节点赋初值。 (2) 从左下角开始计算 J， α，β ，和γ 。 (3) 计算(3.236)式，求 x, y，经过多次迭代，最后得到各节点坐标值。 网格自动剖分在有限元方法中占有很重要的地位，它直接决定了生成矩阵的特性。目前二维有限元网格剖分算法已经趋于成熟和完善，几乎可以处理任意复杂的场域，已不需要一般性工作。 但网格的自适应细分仍然是一个热门课题，目的是使网格分别自动适应于场域结构或场量分布，使场域中每个单元都能给出几乎相同的计算精度，这样就要求程序本身能自行判断何处的单元需要细分，细分到何种程度。 第3章 有限元法基础 3.8.3 自适应剖分技术 第3章 有限元法基础 自适应技术是靠网格细分与场量计算的循环过程来实现，具体步骤为： (1) 生成开端网格，即采用网格自动细分软件将场域剖成很粗的网格； (2) 求解场量，即形成有限元方程与求解方程； (3) 分析场结果计算误差(对每个单元进行误差分析)； (4) 根据误差分析，确定需要细分的网格单元； (5) 细分局部网格； (6) 返回到(2)，求解场量； (7) 计算最后结果。 第3章 有限元法基础 否 场域信息输入 是 精度满足要求？ 开端网格生成 场量计算 误差分析 网格自动细分 后处理 图3.22 自适应剖分框图 Matlab 提供 PDETOOL 软件包就具有网格自动剖分和自适应细分功能，并且它采用CAD输入方式，既简单又直观，所以可选用它作为有限元计算的前处理程序。 网格初始剖分 采用自适应技术后网格一次细分 网格二次细分 第3章 有限元法基础