ANSYS 有限元分析 后处理 结点解与单元解

Blog Links

• DalNur | 博客总目录

• ANSYS 有限元分析 概述

• ANSYS 有限元分析 坐标系/工作平面

• ANSYS 有限元分析 几何建模

• ANSYS 有限元分析 网格划分

• ANSYS 有限元分析 选择与组件

• ANSYS 有限元分析 修改与编辑

• ANSYS 有限元分析 接触分析

• ANSYS 有限元分析 加载/求解/输出

• ANSYS 有限元分析 后处理 General Postproc

• ANSYS 有限元分析 后处理 结点解与单元解

• ANSYS 有限元分析 命令流 实例

• ANSYS 静力 分析实例 悬臂梁

• nCodeDL 高周疲劳 分析实例 悬臂梁

一、前言

  本文以如下图所示的悬臂梁为例，介绍 ANSYS 后处理中的结点解与单元解的主要区别。

  悬臂梁长度为 60 mm，其横截面尺寸为 H×B = 10 mm × 6 mm。材料为钢材，牌号为 Q235B，其弹性模量为 200 Gpa，泊松比为 0.3。其端部承受集中荷载 P = 100 N，沿梁的长度方向承受均布荷载 q = 1 N/mm² ，如下图所示：

  参考资料： ANSYS工程分析-基础与观念 第三章 ANSYS快速浏览. 李辉煌.

二、前处理

2.1 创建几何

FINISH                                ! 退出当前处理器
/CLEAR,ALL                            ! 清除所有

! Units: N, mm, t

/PREP7                                ! 访问前处理器

BLOCK, 0, 60, -5, 5, -3, 3            ! 创建六面体

/VIEW,, 1, 1, 1                       ! 设置视图
VPLOT                                 ! Plot Volumes

  当用户自定义视图时，即确定观察模型的角度时，需要执行 /VIEW 命令。我们需要确定视线的方向，即从何点看向何点，此时，我们需要给定视线的方向向量，该向量的起点称为视点，ANSYS默认为点 (0，0，1)；向量的终点称为目标点，该点坐标不可更改，始终为全局坐标系的原点即点 (0，0，0)。实际上，/VIEW 命令是用来定义一个新的视点的，第 1 个参数是window 的编号 (graphics window 可以切割成很多个 windows)，默认值为 1，第 2、3、4 个参数分别是视点的 x、y、z 坐标值。此时，我们定义的视线方向是 从点 (1, 1, 1) 看向点 (0, 0, 0) 。

2.2 定义属性

ET, 1, SOLID185                       ! 定义单元类型
MP, EX, 1, 200000                     ! 定义弹性模量
MP, NUXY, 1, 0.3                      ! 定义泊松比

2.3 网格划分

TYPE, 1                               ! 激活单元类型号1
MAT, 1                                ! 激活材料号1

ESIZE, 3                              ! 单元尺寸设置为3
VMESH, ALL                            ! 执行分网

EPLOT                                 ! Plot Elements

FINISH                                ! 退出当前处理器

三、加载与求解

3.1 设置边界条件

/SOLU                                 ! 访问求解器

NSEL, S, LOC, X, 0                    ! 选择x坐标值为0的全部结点
D, ALL, ALL, 0                        ! 指定边界条件为固定约束

ALLSEL,ALL                            ! 全选

3.2 施加均布荷载

NSEL, S, LOC, Y, 5                    ! 选择y坐标值为5的全部结点
SF, ALL, PRES, 1                      ! 施加均布荷载

ALLSEL,ALL                            ! 全选

3.3 施加集中荷载

N1 = NODE(60,-5,-3)                   ! 获取坐标值为(60,-5,-3)的结点编号
N2 = NODE(60,-5,3)                    ! 获取坐标值为(60,-5,3)的结点编号

F, N1, FY, -50                        ! 施加集中荷载
F, N2, FY, -50                        ! 施加集中荷载

ALLSEL,ALL                            ! 全选

EPLOT                                 ! Plot Elements

  梁自由端施加集中荷载时，按理说应该在梁宽度中央结点处直接施加一个 100 N 的荷载，但有时宽度中央不一定存在结点 (本例只是恰好有)，比较保险的方式是把 100 N 分成两个 50 N，分别施加到两个端点上，即在编号为 N1 及 N2 的结点上各施加大小为 50 N 方向为 -y 的集中荷载。其中 NODE 为根据结点坐标值获取对应的结点编号的 ANSYS 内置函数。根据圣维南原理，此种加载方式并不影响远端的计算结果。

3.4 求解

ANTYPE,STATIC                         ! 设置分析类型为静力分析
SOLVE                                 ! 提交求解器进行求解

FINISH                                ! 退出当前处理器

四、后处理

  ANSYS 提供了两个后处理器：通用后处理器 (POST1) 和 时间历程后处理器 (POST26)。通用后处理器 (POST1)：用来观察整个模型在某一时刻的结果。时间历程后处理器 (POST26)：用来观察整个模型在不同时间段或荷载步上的结果，常用于处理瞬态分析和动力分析结果。本算例为静力分析，因此，该模型的后处理主要用到 POST1 处理器。

4.1 显示变形形状

/POST1                                ! 访问POST1后处理器
/VIEW,, 0, 0, 1                       ! 设置视图
PLDISP, 2                             ! 显示结构变形形状且重叠显示结构变形前后的形状。
ANDSCL                                ! 生成变形形状的动画

变形形状 (最大变形为0.17mm)

4.2 显示位移云图

/POST1                                ! 访问POST1后处理器
/VIEW,, 0, 0, 1                       ! 设置视图
PLNSOL, U, Y,                         ! 绘制y方向位移云图
/EDGE,1,0                             ! 仅显示边缘轮廓(对于等值线显示)
/REPLOT 

竖向位移云图 (最大值为0.169mm)

  PLNSOL 为用等值线或云图的方式显示结点处的计算结果；PLESOL 为用等值线或云图的方式显示单元的计算结果。

4.3 显示应力云图

4.3.1 显示连续应力云图

/POST1                                ! 访问POST1后处理器
/VIEW,, 0, 0, 1                       ! 设置视图
PLNSOL, S, X                          ! 绘制x方向连续正应力云图
/EDGE,1,1                             ! 显示所有单元面之间的公共线(对于等值线显示)
/REPLOT                               ! Replot

4.3.2 显示非连续应力云图

/POST1                                ! 访问POST1后处理器
/VIEW,, 0, 0, 1                       ! 设置视图
PLESOL, S, X                          ! 绘制x方向非连续正应力云图
/EDGE,1,1                             ! 显示所有单元面之间的公共线(对于等值线显示)
/REPLOT                               ! Replot

  PLESOL 命令绘制的云图的等应力线 (contour lines) 呈锯齿状，PLNSOL 命令绘制的云图的等应力线平滑连续。

  有限位移元的基本未知量为结点位移，结点位移得到求解后，结构的位移场分布便得到确定，位移场在空间上是连续的，但不一定是平滑的，即其导函数不一定保持连续。事实上，在单元内部，位移场连续且平滑，这是由形函数所决定的，但是在跨越单元边界时，位移场通常是连续但不平滑的。所以，位移场在整个求解域内连续但不平滑，在数学上称之为分片平滑函数。这种分片平滑函数经微分 (应力场基本上是位移场的微分) 后，就变成了分片连续的函数，亦即函数在单元内部保持连续，在跨过单元边界时是不连续的，就如上图 (Element Solution) 所示，等应力线的不连续点发生在单元边界上。

  理论上，只要单元尺寸足够小，这种应力/应变的不连续性也会随之得到改善。事实上，这种不连续性可以作为解答精度的量测基准，即可用于判断网格尺寸是否合适，应力不连续越明显说明单元尺寸过大。当然，可以通过细化网格来改善应力场的不连续性，但是，此种做法没有特别的必要，我们只需要将不连续的应力值做一个简单的数学上的处理，使之连续甚至是平滑，这就是 PLNSOL 命令的功能。所谓数学上的简单处理，就是采用绕结点平均法处理角结点上的应力。

  在单元角点处，若该角点为多个单元的公共结点，则该点处的应力值由绕节点平均法得到。对于高阶单元的边中结点，并不是采用绕结点平均法处理应力，而是由该边的两个端点应力值进行插值得到。角点平均，中点插值。

  由上图可以看出，结点解最大正应力为 144.736 Mpa，单元解最大正应力为 159.299 Mpa，两者间的误差为：

$$err=\frac{159.299-144.736}{144.736}\times 100\%=7.3\%$$

  由于 err > 5%，这说明解答不满足精度要求，网格尺寸需进一步减小，以获得理想的解答。

  采用高阶单元及完全积分模式，可改善 PLESOL 命令绘制的云图的不连续性。理论上，只要假设的形函数能保证应力场在全求解域上连续，那么两种方式得到的应力云图是一致的，但这对形函数的要求就相应提高，不仅仅是连续、一阶可导，甚至得保证至少2阶可导且导函数连续，显然，没有必要。

4.4 提取结点应力/单元结点应力

  有限元模型 = 结点 + 单元

  基本思想：几何离散和分片插值。

  基本步骤：结构离散、单元分析和整体分析。

  离散的含义：用假想的线或面将连续物体分割成由有限个单元组成的集合，且单元之间仅在节点处连接，单元之间的作用仅由节点传递。当单元趋近无限小，节点无限多，则这种离散结构将趋近于实际的连续结构。

  将上述有限元模型的单元类型更改为 SOLID186 (20结点二阶单元) 后，重新分析，并选择 141、142 号单元，进行应力数据的提取，以说明 ANSYS 是如何处理应力数据的，141、142 号单元位置如下图所示：

ESEL,S, , ,141                        ! 选择141号单元
ESEL,A, , ,142                        ! 选择142号单元
/PNUM,NODE,1                          ! 显示结点编号
/REPLOT                               ! Replot


/GRAPHICS,POWER                       ! 开启 PowerGraphics Display
/EFACET,2                             ! 对于二阶单元，输出中间结点信息。
! /EFACET,1                           ! 对于二阶单元，不输出中间结点信息。
/REPLOT                               ! Replot

• 结点解应力结果

PRNSOL,S,COMP                         ! 列出结点上的应力结果，六个分量全部列出。

• 单元解应力结果

PRESOL,S,COMP                         ! 列出单元上的应力数据，六个分量全部列出。

  由上图可以看出，PRESOL 命令并没有列出每个单元的全部结点的应力数据，若全部列出，每个单元应该有 20 个结点，该命令只列出了单元 8 个角结点的应力数据。这是因为，在单元层面上，只有单元角结点的应力数据是由高斯积分点外插得到的，高阶单元边中结点的数据不是通过高斯点外插得到的，而是由所在边的角结点直接平均得到的，ANSYS的逻辑是只显示最重要最直接的结果，但凡能通过简单的运算推得的结果均不显示。

  通过 Main Menu >> General Postproc >> Query Results >> Subgrid Solu ，可查询中间结点的应力值。

  按此方式，可查询 337、669 及 670 三个结点在不同单元上的 x 方向正应力解答，如下图所示：

左侧为141号单元 右侧为142号单元

  由上图可以看出，在 141 号单元上，669 号结点的 Sx = 153.365，该值为 377 号结点与 670 号结点 x 方向正应力的平均值，即 (141.12+165.609) / 2 = 153.3645；同样地，在 142 号单元上，669 号结点的 Sx = 150.868，此值为 148.51 和 153.225 的平均值，即 (148.51+153.225) / 2 = 150.8675。这说明，高阶单元的边中结点的应力值为该边两端点相应应力值的平均值，并不是由高斯点外插得到的。

  那么，从141号单元上，推导得到的 669 号结点的 x 方向正应力值为 153.365，从142 号单元上，推导得到的 669 号结点的 x 方向正应力值为 150.8675，以上两个值为单元解。理论上，从不同单元上，推导得到的同一结点处的应力值完全一致，但这对形函数及单元尺寸的要求格外的高，即要求形函数不但连续且其导函数也得连续，单元尺寸足够细小等。但由此付出的计算成本是难以接受的，也没有必要。只要这种不一致满足一定的条件，我们就认为解答可靠。因此，我们可以通过查看单元解在单元边界上的不连续程度，来考察计算误差是否可接受。

  实际上，应力场在整个求解域上是单值连续函数，即一个位置处有且只有一个数值，因此，各个单元求解得到的同一结点应力值需做数学上的简单处理，即取平均值作为应力在该结点上的实际解答，故 669 号结点最终的表征应力值为 (153.365+150.8675) / 2 = 152.11625。这就是 ANSYS 中，Nodal Solu 中列出的 669 号结点 Sx = 152.12 的整个求解过程。

  以 377 号结点为例，进一步说明应力结果在结点处的处理规则，下表为不同处理方式得到的 6 个应力分量，可以看出，结点解为两个单元解的平均值，因为 377 号结点只是 141 号单元和 142 号单元两个单元的公共结点，所以为这两个单元解的平均值，若某个结点为多个单元共有结点，那结点解则为多个单元解的平均值。

数据处理方式	单元编号	结点编号	Sx	Sy	Sz	Sxy	Syz	Szx
结点解	-	377	144.82	-13.434	-11.516	4.9365	0.82502	6.3010
单元解	141	377	141.12	-19.846	-17.569	10.032	0.22062	10.832
单元解	142	377	148.51	-7.0208	-5.4638	-0.15867	1.4294	1.7695

377号结点的单元解与结点解

  通过，以上各种对比分析，不难发现，所谓的结点解与单元解均指的是应力在同一结点处的解答。比较同一结点的结点解与单元解才有意义，不同点你能比较出来个啥，比较出来又有什么意义，啥问题说明不了。不同的是，对于结点的看待角度有所不同。结点解中结点为整个有限元层面的结点，不同的是，虽然单元解指的是单元结点应力解并不是孤立的，而是具有紧密联系的。结点解与单元解都是指应力在应力 (可以是应变等派生解)

  在有限元分析中，我们最最关注的是各物理量在结点处的取值，这些物理量可以是位移、应力、应变等等。

  为什么单元内任意一点处的位移、应力、应变没有结点处的重要 ？？？

  因为，当单元各结点处的位移、应力、应变确定后，单元内任意一点处的位移、应力、应变便可由结点处各相应物理量插值得到，插值函数通常为单元的形函数，什么时候不通常？高阶单元的应力插值函数为其相对应低阶单元的形函数。

  众所周知，位移在结点上精度最高，应力/应变在积分点上精度最高，按理说，应力结果最具代表性的位置应该是高斯积分点，为什么 ANSYS 中，仍然采用单元结点来表征应力，即在 ANSYS 中列出单元结点应力是方便的，若想列出单元高斯点出的应力有点费劲，需通过间接方式： 把高斯点处的应力结果复制到单元结点上，然后在列出单元结点应力。这是因为 ANSYS 的应力后处理规则保证了采用结点表征应力时，高斯积分点处的应力维持不变，显然，用单元结点更为方便。

  有限元分析的基本未知量为结点位移，因此，有限元分析完成后，最先呈现出来的计算结果就是有限元模型中各结点处的结点位移，而诸如应力、应变等其他物理量均是在结点位移确定后，通过各种数学计算得到的，也没什么高级的计算，基本上就是加减乘除内插外插等。有限元中的各种数学运算并不多高级，真的就只是加减乘除内插外插，只不过是数据量特别庞大，显得运算很高级。

  有限元最精髓的操作就是分片插值和取平均值，各种插值，内插、外插，一顿插，不同插值方式得到的同一位置处的物理量不一致时就取平均值。

  结点位移为有限元分析的原生解答，单元上 (单元结点也在单元上) 任意一点处的位移、应力、应变为派生解答，结点应力、结点应变为该结点全部所属单元在该位置处解答的算数平均值。理论上，结点位移求解后，结合形函数可确定单元位移场 ，通过几何方程可确定单元应变场，通过物理方程可确定单元应力场，但由此单元应力场计算得到的单元应力，在高斯点处精度最高，单元结点处精度最低，在单元结点处计算结果往往具有很大误差甚至错误，因此，需要对单元应力进行数学上的处理，以消除误差，提高精度。

  若有足够的高斯积分点数量，则 ANSYS 采用单元应力修匀法修正单元应力场。对于等参元来说，有了单元内应力局部磨平的方法，可以方便地利用精度较高的高斯积分点 (最忌应力点) 的应力值来改进等参元结点应力的近似性质，这种改进结点应力的方法亦称之为应力插值外推。(单元局部应力磨平法)

  由于位移场连续不一定可导，而应力场由位移场求导得到，不是单纯的就求个导，还有其他代数运算，但肯定是要求个导，所以应力场不一定连续。若应力场不连续，则同一位置处通常是多个单元共享结点处，不同单元推导得到的应力值将无法保持一致，要是都一致，那么应力场在该点处就连续了。显然，在真实世界中，结构在外载作用下，各种物理量是保持连续的，而有限元计算只能保证位移连续，应力、应变不一定保持连续。为了解决这一矛盾，当同一位置处，由不同单元推导出的物理量不一致时，则取需要做数学上的简单处理，即算数平均，以使得该位置处该物理量有唯一取值。在 ANSYS 中，有限元模型各结点上应力/应变等物理量未做平均处理前，即保持各个单元计算出的该点应力/应变值不变，不做任何数学上的平均处理，数据在该点处明显不连续，这个解答称为单元解，显然，同一个结点，从不同的单元推导其应力，可能得到不同的数值。结点上的应力/应变等物理量做了绕结点平均处理后的解答称作结点解，即各物理量在结点上保持连续了，这更符合实际情况，因此，我们往往采用结点解作为有限元最终的解答。而单元解只是结点解的上一个步骤，数据在结点处未做平均处理就是单元解，数据在结点处做了平均处理就是结点解。单元解主要用于查看应力/应变的不连续程度，以表明网格质量和解答精度。如果网格足够细小，有限元的解答足够精确，那么单元解显示的应力云图连续性越好，即由不同单元推导出的同一结点上的应力值几乎是一样的，那这解答肯定是完美的。

  综上，单元解与结点解区分的是在结点处不连续的物理量，如应力、应变等，这些都是派生解答，即都是有结点位移推导得到的。而连续的物理量如位移，没必要区分单元解与结点解，因为，位移在整个求解域上是保持连续的，而结点位移本来就是有限元的原生解。所以，在ANSYS 中，云图绘制时，Nodal Solu 内由位移云图绘制选项，而 Element Solu 内没有。

  以应力解答为例，应变等其他不连续的物理量解答类似，来进一步说明单元解与结点解的区别。通过上述大量的分析可以看出，单元解的全称应该是结点应力单元解，结点指明位置，应力指明物理量，单元解指明在结点上应力这个物理量是由单元推导出来的，其实更清晰点还应该指出是由哪个单元推导出的。显然，结点应力单元解不止一个，若某个结点是 k 单元的公共结点那么该结点应力的单元解就有 k 个，每个单元推导出一个解，共 k 个。那么，结点应力结点解就是结点应力各单元解的平均值。

  例如，在上述悬臂梁分析中，377 号结点是 141 号单元和 142 号单元两个单元的公共结点，那么， 377 号结点处 x 方向上的正应力 Sx 就有两个单元解，一个来自 141 号单元，另一个来自 142 号单元，我们可以按如下方式表述：在 141 号单元上第 377 号结点的 x 方向上的正应力 Sx 的解答为 141.12 Mpa，在 142 号单元上第 377 号结点的 x 方向上的正应力 Sx 的解答为 148.51 Mpa。;由此便可计算得 377 号结点 x 方向上的正应力 Sx 的结点解为：

$$\frac{141.12+148.51}{2}=144.815$$

  这就是 ANSYS 软件中，列出 377 号结点 Sx 结点解 144.82 的由来。;可见，单元解为有限元计算的实际反映，即算出来多少就是多少，原封不动的呈现出来，符合不符合物理实际暂时不管，但是它能反映出解答的好坏，而结点解做了数学上的平均处理，使结果更符合物理实际，解决的是同一位置处应力值必须只有一个的问题。

4.5 查看单元高斯点应力

  - ansys查看积分点(高斯点)上应力的方法

4.6 计算合力

五、单元结点应力的计算

  有限位移元分以结点位移为基本未知量，一旦结点位移确定后，单元内任意一点的位移可由下式求得：
$$\mathbf{d}=\mathbf{N}\cdot {\mathbf{\delta }}_{\mathbf{e}}$$

式中， $\mathbf{N}$ 为形函数矩阵，${\mathbf{\delta }}_{\mathbf{e}}$ 为单元各结点位移列阵。

  根据弹性力学几何方程，单元的应变场由下式得到：

$$\mathbf{\epsilon }={\mathbf{A}}^T\mathbf{d}=\mathbf{B}\cdot {\mathbf{\delta }}_{\mathbf{e}}$$

式中，${\mathbf{A}}^T$ 为微分算子，即求导规则矩阵； $\mathbf{B}$ 为应变矩阵，$\mathbf{B}={\mathbf{A}}^T\cdot \mathbf{N}$。

  根据弹性力学物理方程，单元应力场由下式得到：

$$\mathbf{\sigma }=\mathbf{D}\cdot \mathbf{\epsilon }=\mathbf{D}\cdot \mathbf{B}\cdot {\mathbf{\delta }}_{\mathbf{e}}$$

其中，$\mathbf{D}$ 称为弹性矩阵，它完全取决于弹性材料的弹性模量 $E$ 和泊松比 $\mu$ 。

  对于等参元，高斯积分点为单元内最佳应力点，即高斯积分点上应力精度最高。对于面积和体积单元，通常，应力/应变在高斯点上的计算结果更加准确。然而单元结点和形心位置处的应力/应变结果是人们更为关注的，但在这些位置处，由上式计算得到的应力/应变结果精度较差，尤其是单元结点上。

  基于单元形函数或简化形函数，可以将高斯点处的应力向外插值得到单元结点上的应力值，向内插值得到单元形心处的应力值，由此得到的结果较由物理方程得到的结果在精度上大为改善，内插/外插函数如下表所示：

  Uniform stress cases, like a constant stress triangle, do not require the above processing. 若积分点的数量不足以完成外插工作，如4结点四边形减缩积分单元，只有一个积分点，则无法采用积分点应力外插的方式提高单元内应力的精度，此时，则不做上述处理，单元结点上的应力直接取为物理方程的计算结果。

六、尾声

  以上，便是 ANSYS 结点解与单元解 的简单介绍。

  因篇幅有限，某些内容未做更详细介绍，如有疑问，欢迎邮件交流。

  Email： liyang@alu.hit.edu.cn 。

  仅以此文为我 ANSYS 及有限元相关理论的学习做一个总结。

  与此同时，也希望能够为初学者/有需要的人提供多一点参考。

  本文仅用于个人学习，除此之外，无其他任何用途。

  因个人水平有限，文中难免有所疏漏，还请各位大神不吝批评指正。

  胸藏文墨怀若谷，腹有诗书气自华，希望各位都能在知识的 pāo 子里快乐徜徉。

  本文逻辑清楚，内容详实，引例丰富。

  欢迎大家点赞、评论及转载，转载请注明出处！

  为我打call，不如为我打款！

  最后，祝各位攻城狮们，珍爱生命，保护发际线！

七、参考文献

[1]. 有限单元法. 王勖成 编著. 北京. 清华大学出版社. 2003.

[2]. 有限单元法教程. 王焕定 王伟 编著

[3]. ANSYS工程分析-基础与观念. 李辉煌.