电磁场数值方法电子教案：有限元(FEM)课件

有限元（FEM）有限元（FEM）概述概述v历史历史1943 Courant 最早提出思想最早提出思想20世纪世纪50年代年代 用于飞机设计用于飞机设计1960 Clough在著作中首先提出名称在著作中首先提出名称19641965年间数学家冯康独立地开创有限元年间数学家冯康独立地开创有限元方法并奠定其数学基础方法并奠定其数学基础1965 Winslow首次应用于电气工程问题首次应用于电气工程问题1969 Silvester推广应用于时谐电磁场问题推广应用于时谐电磁场问题v应用范围应用范围广泛地被应用于各种结构工程广泛地被应用于各种结构工程成功地用来解决其他工程领域中的问题成功地用来解决其他工程领域中的问题热传导、渗流、流体力学、空气动力学、土壤力学、热传导、渗流、流体力学、空气动力学、土壤力学、机械零件强度分析、电磁工程问题等等机械零件强度分析、电磁工程问题等等概述历史v电磁工程应用及发展电磁工程应用及发展静态场时变场，闭域开域，线性非静态场时变场，闭域开域，线性非线性，散射，波导、腔体、传输线线性，散射，波导、腔体、传输线 标量有限元发展到矢量有限元标量有限元发展到矢量有限元 高阶矢量有限元高阶矢量有限元单一方法发展到混合方法单一方法发展到混合方法(快速算法快速算法)频域求解发展到时域求解频域求解发展到时域求解(区域分解技术区域分解技术)商用软件：比如商用软件：比如HFSS、ANSYS 电磁工程应用及发展v有限元思想有限元思想1有限元法是有限元法是函数逼近理论、偏微分方程、变分与函数逼近理论、偏微分方程、变分与泛函分析泛函分析的巧妙结合。从数学上分析，有限元法的巧妙结合。从数学上分析，有限元法是是Rayleigh-Ritz-Galerkin法的推广。法的推广。传统的有限元以传统的有限元以变分原理变分原理为基础为基础变分问题就是求变分问题就是求泛函极值泛函极值的问题的问题 直接解法把变分问题化为普通多元函数求极值的问题直接解法把变分问题化为普通多元函数求极值的问题Ritz 寻找一组在全域上解析、而又要在边界上满足强加边界寻找一组在全域上解析、而又要在边界上满足强加边界条件的基函数条件的基函数间接解法间接解法 变分原理变分原理 变分问题与对应的变分问题与对应的边值问题等价边值问题等价 有限元思想1v有限元思想有限元思想2有限元法采取了与变分问题间接解法相反的途径,把所要求的微分方程型数学模型边值问题，首先转化为相应的变分问题，即泛函求极值问题；然后利用剖分插值，离散化变分问题为普通多元函数的极值问题，即最终归结为一组多元的代数方程组，解之即得待求边值问题的数值解。有限元思想2v有限元思想有限元思想3有限元法的核心在于剖分插值，它是将所研究的连续场分割为有限个单元，用比较简单的插值函数来表示每个单元的解，但是它并不要求每个单元的试探解都满足边界条件，而是在全部单元总体合成后再引入边界条件。这样，就有可能对内部和边界上的单元采用同样的插值函数，使方法构造极大地得到简化。有限元思想3v有限元思想有限元思想4由于变分原理的应用，使第二、三类及不同媒质分界面上的边界条件作为自然边界条件在总体合成时将隐含地得到满足，也就是说，自然边界条件将被包含在泛函达到极值的要求之中，不必单独列出，而唯一考虑的仅是强制边界条件（第一类边界条件）的处理，这就进一步简化了方法的构造。有限元思想4v有限元法主要特点有限元法主要特点1离散化过程保持了明显的物理意义。因为变分原理描述了支配物理现象的物理学中的最小作用原理（如力学中的最小势能原理、静电学中的汤姆逊定理等）。因此，基于问题固有的物理特性而予以离散化处理，列出计算公式，当可保证方法的正确性、数值解的存在与稳定性等前提要素。有限元法主要特点1v有限元法主要特点有限元法主要特点2优异的解题能力。与其他数值方法相比较，有限元法在适应场域边界几何形状以及媒质物理性质变异情况的复杂问题求解上，有突出的优点：不受几何形状和媒质分布的复杂程度限制；不同媒质分界面上的边界条件是自动满足的；不必单独处理第二、三类边界条件；离散点配置比较随意，通过控制有限单元剖分密度和单元插值函数的选取，可以充分保证所需的数值计算精度。有限元法主要特点2v有限元法主要特点有限元法主要特点3可方便地编写通用计算程序，使之构成模块化的子程序集合。容易并行。从数学理论意义上讲，有限元作为应用数学的一个分支，它使微分方程的解法与理论面目一新，推动了泛函分析与计算方法的发展。有限元法主要特点3&3.1变分原理与尤拉方程变分原理与尤拉方程v在微积分学形成的初期，以数学物理问题为在微积分学形成的初期，以数学物理问题为背景，与背景，与多元函数的极值问题多元函数的极值问题相对应，就已相对应，就已经在几何、力学上提出了经在几何、力学上提出了若干求解泛函极值若干求解泛函极值的问题。的问题。v例如最速降线问题，即在于研究当质点从定例如最速降线问题，即在于研究当质点从定点点A自由下滑到定点自由下滑到定点B时，为使滑行时，为使滑行时间最短时间最短，试求指点应延着怎样形状的光滑试求指点应延着怎样形状的光滑轨道轨道下滑。下滑。&3.1变分原理与尤拉方程在微积分学形成的初期，以数学物理问dxdsA(x1,y1)B(x2,y2)xyO沿曲线滑行弧线所需时间为沿曲线滑行弧线所需时间为滑行总时间为滑行总时间为dxdsA(x1,y1)B(x2,y2)xyO沿曲线滑行电磁场数值方法电子教案：有限元(FEM)课件泛函的极值（泛函的极值（max或或min）问题就称为变分问题。）问题就称为变分问题。对一般一般问题而言，可而言，可导出下列出下列对应于一个自于一个自变量量 、单个函数个函数及其及其导数数的已知函数的已知函数函数族函数族 仅有一个有一个能使定能使定积分分达到极小达到极小值 泛函的极值（max或min）问题就称为变分问题。对一般问题间接解法间接解法是将变分问题转化为是将变分问题转化为尤拉方程尤拉方程（微分方程）的定解问题，即边值问题来求解。（微分方程）的定解问题，即边值问题来求解。称之称之为的的变分分，它反映了，它反映了整个函数的整个函数的变化量化量 相相应于于变分分的泛函增量的泛函增量为间接解法是将变分问题转化为尤拉方程称之为的变分，它反映了整个电磁场数值方法电子教案：有限元(FEM)课件电磁场数值方法电子教案：有限元(FEM)课件任意给定的微量实参数任意给定的微量实参数 满足足齐次次边界条件的可微函数界条件的可微函数 极值极值任意给定的微量实参数 满足齐次边界条件的可微函数 极值 电磁场数值方法电子教案：有限元(FEM)课件简写为简写为 只差一个数只差一个数值因子因子 极值函数解极值函数解必须满足的必要条件必须满足的必要条件 等同于等同于简写为 只差一个数值因子 极值函数解必须满足的必要条件等同于电磁场数值方法电子教案：有限元(FEM)课件泛函泛函的极值问题的的极值问题的尤拉方程尤拉方程泛函的极值问题的尤拉方程简单函数简单函数简单泛函简单泛函自变量的微分自变量的微分 表示自变量值的微小表示自变量值的微小变化变化函数变分函数变分 表示函数形式的表示函数形式的微小变化，其中微小变化，其中 是正的任意给是正的任意给定的常数，定的常数，为可取函数为可取函数 引起的函数值变化可利用引起的函数值变化可利用Taylor级级数展开数展开函数增量的线性部分函数增量的线性部分函数的一阶微分简称微分函数的一阶微分简称微分函数的函数的n阶微分表示为阶微分表示为 引起的泛函值的变化可展开为引起的泛函值的变化可展开为定义：泛函的一阶变分简称变分，是泛定义：泛函的一阶变分简称变分，是泛函增量的线性主部函增量的线性主部同样有二阶直到同样有二阶直到n阶变分阶变分简单函数简单泛函自变量的微分 表示自变量值的微小变化函数变简单函数简单函数简单泛函简单泛函自变量在自变量在 上变化时，函数有极大上变化时，函数有极大和极小点。和极小点。极大点极大点 取极大值取极大值（在（在 领域）领域）极小点极小点 取极小值取极小值（在（在 领域）领域）取极值条件：一阶微分为零，取极值条件：一阶微分为零，的解的解用二阶微分可以判断该点为极大用二阶微分可以判断该点为极大（），极小（），极小（），），还是拐点还是拐点函数定义空间变化时（曲线簇）使值域函数定义空间变化时（曲线簇）使值域数值为极大和极小数值为极大和极小极大曲线极大曲线 是泛函是泛函 极大极大极小曲线极小曲线 是泛函是泛函 极小极小泛函极值条件为一阶变分为零：泛函极值条件为一阶变分为零：的解的解用泛函二阶变分判断极值点的特性：用泛函二阶变分判断极值点的特性：简单函数简单泛函自变量在 上变化时，函数有v泛函的极值泛函的极值问题就称为问题就称为变分问题变分问题 v变分问题与边值问题变分问题与边值问题等价等价 v有限元正是间接求解变分问题过程的有限元正是间接求解变分问题过程的逆逆过程过程 v泛函取极值的过程中泛函取极值的过程中 第二、第三类边界条件为第二、第三类边界条件为自然边界条件自然边界条件无条件变分问题无条件变分问题第一类边界条件为第一类边界条件为强加边界条件强加边界条件 条件变分条件变分 泛函的极值问题就称为变分问题&3.2与线性问题等价的变分问题与线性问题等价的变分问题v与齐次边值问题等价的变分问题与齐次边值问题等价的变分问题&3.2与线性问题等价的变分问题与齐次边值问题等价的变分问题与泊松方程齐次第三类边值问题等价的变分问题与泊松方程齐次第三类边值问题等价的变分问题与泊松方程齐次第二类边值问题等价的变分问题与泊松方程齐次第二类边值问题等价的变分问题 与泊松方程齐次第三类边值问题等价的变分问题与泊松方程齐次第一类边值问题等价的变分问题与泊松方程齐次第一类边值问题等价的变分问题 混合型边界条件混合型边界条件 与泊松方程齐次第一类边值问题等价的变分问题 v与非齐次边值问题等价的变分问题与非齐次边值问题等价的变分问题与泊松方程非齐次第三类边值问题等价的变分问与泊松方程非齐次第三类边值问题等价的变分问题题 与非齐次边值问题等价的变分问题与泊松方程非齐次第二类边值问题等价的变分问与泊松方程非齐次第二类边值问题等价的变分问题题 与泊松方程非齐次混合型边值问题与泊松方程非齐次混合型边值问题 与泊松方程非齐次第二类边值问题等价的变分问题 v二维问题二维问题 体积分变成面积分，面积分变成线积分体积分变成面积分，面积分变成线积分v轴对称场轴对称场 v分层介质中的变分问题分层介质中的变分问题 变为问题中，由于介质分界面上的边界条件为齐变为问题中，由于介质分界面上的边界条件为齐次自然边界条件，所以泛函取极值时自动满足，次自然边界条件，所以泛函取极值时自动满足，不必另行处理。不必另行处理。二维问题 半线系统半线系统 无遗漏、无多余的覆盖无遗漏、无多余的覆盖&3.3基于变分原理的差分方程基于变分原理的差分方程 半线系统&3.3基于变分原理的差分方程 梯形积分公式梯形积分公式线性线性LagrangeLagrange插值公式插值公式梯形公式梯形公式梯形积分公式线性Lagrange插值公式梯形公式电磁场数值方法电子教案：有限元(FEM)课件电磁场数值方法电子教案：有限元(FEM)课件当（当（i，j）为内点）为内点 右边点而非角点右边点而非角点 上边点而非角点（上边点而非角点（i，N）角点（角点（M，N）当（i，j）为内点 右边点而非角点 上边点而非角点（i，N）&3.4有限元法求解有限元法求解v给出与待求边值问题相应的给出与待求边值问题相应的泛函及其等价变分问泛函及其等价变分问题题v应用应用有限单元剖分场区域，有限单元剖分场区域，并选取并选取相应的插值函相应的插值函数数v将将变分问题离散化为一个多元函数的极值变分问题离散化为一个多元函数的极值问题，问题，导出一组联立的代数方程导出一组联立的代数方程 单元分析单元分析/总体合成总体合成/强加边界条件处理强加边界条件处理 v选择选择适当的代数解法，解有限元方程适当的代数解法，解有限元方程，即得待求，即得待求边值问题的近似解（数值解）边值问题的近似解（数值解）v检验（附加计算）检验（附加计算）&3.4有限元法求解给出与待求边值问题相应的泛函及其等价变分任务任务离散化离散化方式方式场域或物体分为有限个子域，如：三角形、四边场域或物体分为有限个子域，如：三角形、四边形、四面体、六面体等形、四面体、六面体等内容内容单元数量、大小、排列单元数量、大小、排列任务任务选择插值函数选择插值函数方式方式选择插值函数的类型，如多项式，用节点（图形选择插值函数的类型，如多项式，用节点（图形顶点）的场值求取子域中各点的场的近似值。一顶点）的场值求取子域中各点的场的近似值。一般用多项式，其次数与节点数有关般用多项式，其次数与节点数有关内容内容插值函数、形式、次数插值函数、形式、次数任务离散化方式场域或物体分为有限个子域，如：三角形、四边形、任务任务建立单元特征式建立单元特征式（单元分析）（单元分析）方式方式推导单元系数矩阵、取决于插值函数、单元推导单元系数矩阵、取决于插值函数、单元几何形状、单元材质。相应的变分问题几何形状、单元材质。相应的变分问题内容内容找到对应的变分问题，将已知插值函数进行找到对应的变分问题，将已知插值函数进行微分、积分运算。整理出单元形函数，单元微分、积分运算。整理出单元形函数，单元系数矩阵系数矩阵任务任务建立系统有限元建立系统有限元（总体合成）（总体合成）方式方式把单元特征式采用简单处理方法加以合并，把单元特征式采用简单处理方法加以合并，然后表示出整域上的线性方程组，节点互联然后表示出整域上的线性方程组，节点互联处的场值相同，一般此过程由计算机自动完处的场值相同，一般此过程由计算机自动完成成内容内容有限元方程有限元方程任务建立单元特征式（单元分析）方式推导单元系数矩阵、取决于插任务任务求解有限元方程求解有限元方程方式方式考虑边界条件并修改上一步得到的方程，采考虑边界条件并修改上一步得到的方程，采用适当的方法求解线性方程组，求得节点处用适当的方法求解线性方程组，求得节点处未知场的数值。再由插值函数求域中任一点未知场的数值。再由插值函数求域中任一点的值的值内容内容任一点的场值任一点的场值任务任务附加计算附加计算方式方式由场值求取其他关心的重要参数，如电荷分由场值求取其他关心的重要参数，如电荷分布，电流分布、电压分布等。可由相应的物布，电流分布、电压分布等。可由相应的物理规律，经离散化处理后，得到各单元的相理规律，经离散化处理后，得到各单元的相应表达式应表达式内容内容关心的其他物理量的值关心的其他物理量的值任务求解有限元方程方式考虑边界条件并修改上一步得到的方程，采&3.4.1场域剖分场域剖分 v区域离散的方式将区域离散的方式将影响影响计算机计算机内存需求内存需求、计计算时间算时间和数值结果的和数值结果的精确度精确度 一维区域一维区域 短直线段短直线段 二维区域二维区域小三角形或矩形小三角形或矩形 三维区域三维区域四面体、三棱柱或矩形块四面体、三棱柱或矩形块&3.4.1场域剖分 区域离散的方式将影响计算机内存需求、计电磁场数值方法电子教案：有限元(FEM)课件电磁场数值方法电子教案：有限元(FEM)课件v划分单元的基本原则划分单元的基本原则总原则总原则 在满足精度要求的前提下，尽量减少剖分单元数目以在满足精度要求的前提下，尽量减少剖分单元数目以节省存储量和计算时间节省存储量和计算时间 其他其他需要详尽了解的部位要切分得细小，其他部位可以粗需要详尽了解的部位要切分得细小，其他部位可以粗糙一些，几何形状变化剧烈的地方电磁场也变化大，糙一些，几何形状变化剧烈的地方电磁场也变化大，单元要细小一些单元要细小一些所有单元应所有单元应接近等边三角形接近等边三角形 当有当有曲线边界和复杂形状边界曲线边界和复杂形状边界时就需要把复杂形状用时就需要把复杂形状用标准形状来逼近标准形状来逼近 节点、分割线或分割面节点、分割线或分割面应应设置设置在几何形状和介质形状在几何形状和介质形状发生发生突变处突变处 划分单元的基本原则v注意事项注意事项 各单元只能在顶点处相交各单元只能在顶点处相交 不同单元在边界处相连，既不能相互分离又不能不同单元在边界处相连，既不能相互分离又不能相互重叠相互重叠 注意事项 v单元、节点编号单元、节点编号单元编号单元编号、全局节点编号全局节点编号、局部节点编号局部节点编号 单元编号单元编号用一组整数给单元编号用一组整数给单元编号没有特殊要求，只要方便没有特殊要求，只要方便 一般按内部单元、第一类、第二类边界条件单元的顺一般按内部单元、第一类、第二类边界条件单元的顺序进行序进行 节点编号节点编号完整描述应包括它的完整描述应包括它的坐标值、局部编号和全局编码坐标值、局部编号和全局编码局部编号表示它在单元中的位置局部编号表示它在单元中的位置：逆时针方向逆时针方向 用两组整数编号用两组整数编号单元、节点编号引入引入3M的整型数组的整型数组 节点的局部编号节点的局部编号 单元编号单元编号 节点全局编号节点全局编号 1241254233524564引入3M的整型数组 节点的局部编号 单元编号 节点全局编号&3.4.2分片插值与形状函数分片插值与形状函数 v插值插值 用一个简单函数去近似代替真实函数，二者在某用一个简单函数去近似代替真实函数，二者在某些积分点上具有相同的函数值甚至直到某阶导数些积分点上具有相同的函数值甚至直到某阶导数值值 插值函数插值函数:一阶（线性）、二阶（二次）、或高一阶（线性）、二阶（二次）、或高阶多项式阶多项式&3.4.2分片插值与形状函数 插值:待定系数待定系数:P点的插值函数点的插值函数（也称为展开函数或基函数），（也称为展开函数或基函数），通常代表了有限单元上通常代表了有限单元上用来逼近带求场分布的近似规律用来逼近带求场分布的近似规律 它们只有在单元它们只有在单元e内才不为零，而在单元内才不为零，而在单元e外均为零外均为零 v形状函数形状函数基函数与单元的形状尺寸有关基函数与单元的形状尺寸有关 :待定系数:P点的插值函数它们只有在单元e内才不为线性三角形单元线性三角形单元单元剖分得足够小，以致可将其上的场量看作不变单元剖分得足够小，以致可将其上的场量看作不变 三角形单元节点数与插值函数的待定系数的个数相等。三角形单元节点数与插值函数的待定系数的个数相等。线性三角形单元三角形单元节点数与插值函数的待定系数的个数相等三角元三角元e上的线性插值基函数上的线性插值基函数 三角元e上的线性插值基函数 电磁场数值方法电子教案：有限元(FEM)课件的几何意义的几何意义 形状函数形状函数 表示以表示以 为一顶点，为一顶点，为对边的三角形面积与三角形单元面积之比。为对边的三角形面积与三角形单元面积之比。的几何意义 形状函数 表示以 为一顶点，为对边的三角形面积v形状函数性质形状函数性质 观察点观察点位于第位于第个节点的对边上时，个节点的对边上时，为零为零 保证了单元两侧解的连续性保证了单元两侧解的连续性 形状函数性质 观察点位于第个节点的对边上时，为零 保证了单元&3.4.3 有限元方程建立有限元方程建立 v齐次自然边界条件下拉氏方程的有限元方程齐次自然边界条件下拉氏方程的有限元方程静电场的边值问题静电场的边值问题 二维直角坐标系二维直角坐标系&3.4.3 有限元方程建立 齐次自然边界条件下拉氏方程的有用三角形剖分，得到用三角形剖分，得到E个单元，个节点，每个单元上的能量泛函为个单元，个节点，每个单元上的能量泛函为 处理过程与变分处理过程与变分原理导出的差分原理导出的差分法有什么异同呢法有什么异同呢？用三角形剖分，得到E个单元，个节点，每个单元上的能量泛函为单元分析单元分析单元分析电磁场数值方法电子教案：有限元(FEM)课件不是坐标不是坐标x x，y y的函数，的函数，而是单元上的节点位置（已知数值）及其场量的函数而是单元上的节点位置（已知数值）及其场量的函数 单元上的泛函表达式单元上的泛函表达式 单元的系数矩阵，它只与剖分有关，而非单元的系数矩阵，它只与剖分有关，而非x x，y y的函数的函数 不是坐标x，y的函数，单元上的泛函表达式 单元的系数矩阵，它电磁场数值方法电子教案：有限元(FEM)课件总体合成总体合成迭加：交汇于同一节点的单元泛函方程迭加：交汇于同一节点的单元泛函方程 对称矩阵对称矩阵 组合时，将各单元系数矩阵扩展为组合时，将各单元系数矩阵扩展为 阶方阵，直接相加阶方阵，直接相加 总体合成迭加：交汇于同一节点的单元泛函方程 对称矩阵 组合例：如何由单元矩阵组合成总体矩阵例：如何由单元矩阵组合成总体矩阵 例：如何由单元矩阵组合成总体矩阵 所对应的等价变分问题为所对应的等价变分问题为所对应的等价变分问题为节点总体编号节点总体编号123456x002244y020202节点局节点局部编号部编号 ijm132342354564单元单元 号号节点总体节点总体编号编号节点总体编号123456x002244y020202 对于对角元素对于对角元素 有贡献的只有以有贡献的只有以 为共同顶点的单元为共同顶点的单元 对非对角元素对非对角元素 有贡献的只是以有贡献的只是以 为公共边的单元为公共边的单元 零元素由不在一个单元、不相干的节点产生零元素由不在一个单元、不相干的节点产生 对于对角元素 有贡献的只有以 为共同顶点的单元 对非对角元素对称对称 主对角线元素占优，正定主对角线元素占优，正定 稀疏稀疏 带状带状 对称 主对角线元素占优，正定 稀疏 带状 v减小系数矩阵的最大半带宽减小系数矩阵的最大半带宽 取决于场域上单元任二节点取决于场域上单元任二节点总体编号的最大差值总体编号的最大差值D D 减小系数矩阵的最大半带宽 取决于场域上单元任二节点先选与其他节点联系最少的节先选与其他节点联系最少的节点作为起始节点，然后将相邻点作为起始节点，然后将相邻节点编为紧接着的号数，使相节点编为紧接着的号数，使相邻节点的编号数均为相差不多邻节点的编号数均为相差不多的数的数节点全局编号的顺序节点全局编号的顺序先选与其他节点联系最少的节点作为起始节点，然后将相邻节点编为D=2，B=3 D=3，B=4 一般沿场域的窄边（节点数少的）编号一般沿场域的窄边（节点数少的）编号 需占内存数需占内存数6*3=186*4=24D=2，B=3 D=3，B=4 一般沿场域的窄边（节点数少的v非齐次自然边界条件下拉氏方程的有限元方程非齐次自然边界条件下拉氏方程的有限元方程二维直角坐标系二维直角坐标系 非齐次自然边界条件下拉氏方程的有限元方程二维直角坐标系 由非齐次自然边界条件引起由非齐次自然边界条件引起 由非齐次自然边界条件引起 剖分时应尽量使场域边界分段线性化剖分时应尽量使场域边界分段线性化 应先编不含非齐次边界条件的单元，再编含非齐次边界条件的单元应先编不含非齐次边界条件的单元，再编含非齐次边界条件的单元 应先编不含非齐次边界条件的节点应先编不含非齐次边界条件的节点1，2，N，再编含非齐次边界条件节点再编含非齐次边界条件节点N+1，N+2，。剖分时应尽量使场域边界分段线性化 应先编不含非齐次边界条件的将正对着边界的那个将正对着边界的那个三角形单元顶点编为三角形单元顶点编为i jm的长度为的长度为边上任一点到边上任一点到j点的距离为点的距离为，该点上的电位为，该点上的电位为将正对着边界的那个三角形单元顶点编为i jm的长度为边上任一电磁场数值方法电子教案：有限元(FEM)课件写为矩阵形式写为矩阵形式 由由 引起引起由由 引起引起写为矩阵形式 由 引起由 引起的值只和与它相联系的节点有关的值只和与它相联系的节点有关 的值只和与它相联系的节点有关 电磁场数值方法电子教案：有限元(FEM)课件v齐次自然边界条件下泊松方程的有限元方程齐次自然边界条件下泊松方程的有限元方程静电场的边值问题静电场的边值问题 二维直角坐标系二维直角坐标系 齐次自然边界条件下泊松方程的有限元方程静电场的边值问题 二维由密度由密度 引起引起 三角元剖分三角元剖分 由密度 引起 三角元剖分 剖分的三角元面积剖分的三角元面积 较小较小 提出积分号外提出积分号外 近似为常量近似为常量 剖分的三角元面积 较小 提出积分号外 近似为常量 积分积分 x-y平面平面 平面平面 积分 x-y平面 平面 由二重积分的变换式由二重积分的变换式雅克比式，得到面积元素变换式为雅克比式，得到面积元素变换式为 由二重积分的变换式雅克比式，得到面积元素变换式为 电磁场数值方法电子教案：有限元(FEM)课件为为N阶列阵阶列阵 完全由场域内的电荷密度分布引起完全由场域内的电荷密度分布引起 为N阶列阵 完全由场域内的电荷密度分布引起 v非齐次自然边界条件下泊松方程的有限元方程非齐次自然边界条件下泊松方程的有限元方程二维直角坐标系二维直角坐标系 非齐次自然边界条件下泊松方程的有限元方程二维直角坐标系 非齐次边界条件的边界单元才存在非齐次边界条件的边界单元才存在 非齐次边界条件的边界单元才存在 电磁场数值方法电子教案：有限元(FEM)课件电磁场数值方法电子教案：有限元(FEM)课件&3.4.3有限元方程的计算有限元方程的计算 v迭代法迭代法系数矩阵正定系数矩阵正定 迭代法收敛迭代法收敛非齐次边界条件下的拉氏方程的有限元方程式非齐次边界条件下的拉氏方程的有限元方程式&3.4.3有限元方程的计算 迭代法系数矩阵正定 迭代法收敛任一节点任一节点i的第的第n+1次场值的次场值的GS迭代公式迭代公式 非齐次自然边界条件的节点编号为非齐次自然边界条件的节点编号为 其他节点编为其他节点编为 SOR迭代格式迭代格式任一节点i的第n+1次场值的GS迭代公式 非齐次自然边界条件v强加边界条件的处理强加边界条件的处理 设第设第m个节点上具有强加边界条件个节点上具有强加边界条件=1=0强加边界条件的处理 设第m个节点上具有强加边界条件=1=0电磁场数值方法电子教案：有限元(FEM)课件v处理强加边界条件小结处理强加边界条件小结将对角线元素将对角线元素 置置1该行的右端项改为强加电位值该行的右端项改为强加电位值 第第m行与行与m列的其他元素全部置零列的其他元素全部置零除第除第m行外，其他各行右端项为原右端项减去强加电行外，其他各行右端项为原右端项减去强加电位置位置 与对应第与对应第m列未变换前的系数的乘积列未变换前的系数的乘积 若强加边界条件节点共有若强加边界条件节点共有N0N0个，则如法处理个，则如法处理N0N0次次 N阶方阵阶方阵 与强加边界节点电位与强加边界节点电位值及有关量组成的列阵值及有关量组成的列阵 处理强加边界条件小结若强加边界条件节点共有N0个，则如法处理v系数矩阵的存储系数矩阵的存储利用其对称性利用其对称性、稀疏性稀疏性 等带宽存储等带宽存储 系数矩阵的存储等带宽存储 变带宽存储变带宽存储 具体作法具体作法 将将KK矩阵的下三角部分矩阵的下三角部分带内元素带内元素按按行行的的顺序顺序，从各行的第一个非零元素起，从各行的第一个非零元素起，至主对角元素（包括带内的零元素及对角元素）为止的元素至主对角元素（包括带内的零元素及对角元素）为止的元素（共有（共有MD+1MD+1），依次），依次存入一维数组存入一维数组KL(n)KL(n)中中 变带宽存储 具体作法 将K矩阵的下三角部分带内元素按行的说明说明KL(n)数组中各元素在原矩阵数组中各元素在原矩阵K中的行、列位置中的行、列位置 指针指针数组数组L(n)L(n)存储原存储原K阵中主对角元素阵中主对角元素 在KL(n)中的位置 说明KL(n)数组中各元素在原矩阵K中的行、列位置 指针v需要计算的量需要计算的量K的下三角阵各行的半带宽的下三角阵各行的半带宽MD 第第i行的半带宽行的半带宽MD=i-min(j)主对角元素在主对角元素在KL(n)的地址的地址JOJO=L(i)=L(i-1)+MD+1 非主对角元素在非主对角元素在KL(n)的地址的地址JO JO=L(i)-(i-j)非主对角元素的地址为该行主对角元素的地址向前推移非主对角元素的地址为该行主对角元素的地址向前推移(i-j)个位置个位置 所有以节点所有以节点i为顶点的各单元上节点总体编号的最小值为顶点的各单元上节点总体编号的最小值 前一行主对角元的地址前一行主对角元的地址 需要计算的量所有以节点i为顶点的各单元上节点总体编号的最小值电磁场数值方法电子教案：有限元(FEM)课件行号行号iMD在在KL(n)中中位置（即位置（即L(n)值）值）在在KL(n)中位置中位置第一行第一行1第二行第二行13的的JO2第三行第三行26的的JO4 4的的JO5第四行第四行29的的JO7 7的的JO8第五行第五行212的的JO1010的的JO11第六行第六行316的的JO1313的的JO1414的的JO15K阵中各元素在一维数组阵中各元素在一维数组KL(n)中的位置为中的位置为 KL(1)KL(2)KL(3)KL(4)KL(5)KL(6)KL(7)KL(8)KL(9)KL(10)KL(11)KL(12)KL(13)KL(14)KL(15)KL(16)行号iMD在KL(n)中在KL(n)第一行1第二行13的JOv高斯消元法求解线性代数方程组高斯消元法求解线性代数方程组消元消元将原方程的系数矩阵化为对角元素为将原方程的系数矩阵化为对角元素为1的上三角阵的上三角阵 回代回代从最末一个方程（只包含一个未知数）求出一个未知从最末一个方程（只包含一个未知数）求出一个未知数，再代入上一个方程，求出另一个未知数数，再代入上一个方程，求出另一个未知数 高斯消元法求解线性代数方程组消元公式消元公式第一个方程的系数第一个方程的系数及右端项的处理及右端项的处理 回代公式回代公式消元公式第一个方程的系数回代公式v电场强度的计算电场强度的计算一单元中一单元中E E的数值解的数值解 与坐标与坐标x x、y y无关，无关，仅与单元节点的坐标仅与单元节点的坐标及电位有关及电位有关 电场强度的计算一单元中E与坐标x、y无关，单元的节点一般与相邻几个单元相联系单元的节点一般与相邻几个单元相联系 节点场强必须将有关单元对场强综合贡献节点场强必须将有关单元对场强综合贡献 分别总和有关单元的分别总和有关单元的，然后取算术平均值，然后取算术平均值 单元的节点一般与相邻几个单元相联系 节点场强必须将有关单元对分别将有关单元的分别将有关单元的 加权后总和，以单元面积加权后总和，以单元面积 为权因子，然后取其平均值为权因子，然后取其平均值 场强的方向场强的方向 剖分细剖分细 对称对称 分界面处细剖分分界面处细剖分 高精度：高精度：分别将有关单元的 加权后总和，以单元面积 为权因子，然后取其例例1：计算方同轴线间的电位分布：计算方同轴线间的电位分布边值问题边值问题等价变分问题等价变分问题有限元方程有限元方程&3.5 应用举例应用举例例1：计算方同轴线间的电位分布边值问题等价变分问题有限元方程 每个单元的三个顶点编号每个单元的三个顶点编号每个顶点的每个顶点的xyxy坐标坐标内边界节点编号内边界节点编号边界边界外边界节点编号外边界节点编号节点总数、单元总数节点总数、单元总数1、将几何信息读入代码、将几何信息读入代码 每个单元的三个顶点编号每个顶点的xy坐标内边界节点编号边界2、计算总体系数矩阵、计算总体系数矩阵K计算每个三角形单元的矩阵系数（子程序）计算每个三角形单元的矩阵系数（子程序）计算三角形单元总面积计算三角形单元总面积计算中间矩阵计算中间矩阵Se单个三角形单元的矩阵系数合成总体系数矩阵单个三角形单元的矩阵系数合成总体系数矩阵3、处理强加边界条件、处理强加边界条件找到内、外边界节点，为其电位赋值找到内、外边界节点，为其电位赋值 节点编号，重新排序节点编号，重新排序4、求解有限元方程、求解有限元方程2、计算总体系数矩阵K计算每个三角形单元的矩阵系数（子程序）例例2：波导场的有限元解：波导场的有限元解v用有限元法求解的优越用有限元法求解的优越 将截面逐步细分，将导致本征值向极限单调减小将截面逐步细分，将导致本征值向极限单调减小 可保证有一个更快的收敛速度趋于本征值可保证有一个更快的收敛速度趋于本征值 对难处理的边界形状更容易处理，而不导致非对对难处理的边界形状更容易处理，而不导致非对称矩阵称矩阵 波导中的奇点无需特别处理波导中的奇点无需特别处理 例2：波导场的有限元解用有限元法求解的优越 v假设假设 波导壁为完纯导体波导壁为完纯导体 波导内的介质系均匀、线性且各向同性的理想介波导内的介质系均匀、线性且各向同性的理想介质质 波导中无自由电荷和传导电流波导中无自由电荷和传导电流 波导工作在匹配状态，具有均匀截面，所以在分波导工作在匹配状态，具有均匀截面，所以在分析时只考虑入射波，无反射波析时只考虑入射波，无反射波 假设 波导中传播的电磁波可分为波导中传播的电磁波可分为TETE波波或或TMTM波波 求解相应的场纵向分量所描述的定解问题求解相应的场纵向分量所描述的定解问题 波导中传播的电磁波可分为TE波或TM波 求解相应的场纵向分量TE TM TE TM 选取矩形波导选取矩形波导BJ-100 中中TE波的截止波长波的截止波长 的分布问题进行分析的分布问题进行分析 三角剖分、三节点插值三角剖分、三节点插值 选取矩形波导BJ-100 中TE波的截止波长 的分布问题进行电磁场数值方法电子教案：有限元(FEM)课件电磁场数值方法电子教案：有限元(FEM)课件对称阵对称阵 对称正定阵对称正定阵（广义代数特征值问题）（广义代数特征值问题）对称阵对称阵 的特征值问题的特征值问题 对称阵 对称正定阵（广义代数特征值问题）对称阵 的特征值平方根法平方根法 下三角阵下三角阵 平方根法 下三角阵 利用豪斯豪尔法利用豪斯豪尔法：对称三角矩阵：对称三角矩阵 经过多次相似变换经过多次相似变换 满足指定精度的对角阵满足指定精度的对角阵 对角阵的每个元素就是其特征值对角阵的每个元素就是其特征值 反变换反变换 求得的一系列特征值中非负的求得的一系列特征值中非负的最小非零最小非零特征值特征值就给出相应波导中最低型（主模）的截止波长就给出相应波导中最低型（主模）的截止波长 利用豪斯豪尔法：对称三角矩阵 经过多次相似变换 满足指定精截截 止止 波波 长长理论值理论值相对误差相对误差波型波型(cm)(cm)数值解数值解TETE10104.5724.5724.5124.5124.5694.5694.5744.5740.0660.066TETE20202.2862.2862.1932.1932.2622.2622.2752.2751.051.05TETE01012.0322.0321.8471.8471.9791.9792.612.61TETE11111.8571.8571.6351.6351.7901.7903.613.61TETE30301.5241.5241.4101.4101.4881.4882.362.36TETE21211.5191.5191.2501.2501.4401.4405.205.20TETE31311.2191.2191.0211.0211.1411.1416.406.40TETE40401.1431.1430.9430.9431.0931.0934.374.37TETE02021.0161.0160.9220.9220.9230.9239.159.15TETE41410.9960.9960.8590.8590.9160.9168.038.03TETE22220.9280.9280.7950.7950.8910.8913.993.99截 止 波 长理论值相对误差波型(cm)数值解TE104.5