1、三维应力有限元分析 及传热分析方法,杨建国浙江工业大学 化工机械设计研究所,第一部分 三维应力分析有限元法,第一部分三维应力分析有限元法,等应变三角形单元中的坐标为两个方向，而在三维状态下节点坐标变成了3个。因而相应的节点位移也就采用三个分量来表示，即：u、v和w。此时存在6个应变向量：,L为63的偏微分算子矩阵。,第一部分三维应力分析有限元法,此关系与平面应力时的关系类似，只不过此时应力与应变分别包含六个参量，而平面应力分析时分别仅为3个分量。,第一部分三维应力分析有限元法,对于有M个节点的弹性应力单元来说，单元内任意一点的位移可以通过插值函数的方式获得，相对于三个方向，获得如下三个关系式：
2、 ui、vi和wi分别为节点i的三个方向的位移；Ni为插值函数,第一部分三维应力分析有限元法,单元节点的位移向量（列向量） 位移向量可以记为如下矩阵形式：,1M,33M,第一部分三维应力分析有限元法,63M,第一部分三维应力分析有限元法,应用虚功原理可建立如下的节点位移与节点力之间的关系 定义如下的单元刚度矩阵,3M 3M对称矩阵,第一部分三维应力分析有限元法,单元的刚度矩阵之后，就可以用前面经常提到的组装的概念，将单元的刚度矩阵组装成为整体刚度矩阵，其组装过程与平面应力的组装过程一样。最终我们获得了整个结构的矩阵方程形式为：,第一部分三维应力分析有限元法,通过已知的力边界条件以及位移边界条件
3、，可以最终解以上的方程获得节点位移及相应的节点力。获得以上结果后，则可以通过位移应变关系，应变应力关系计算出结构的应变与应力。,第二部分热过程的数值模拟,第二部分热过程的数值模拟,早期的关于焊接热过程的解析计算结果表明计算结果与实际测量结果之间存在较大的误差，这主要与这些解析计算的假设条件有关，通常的假设条件包含以下几条： a）应用线热源或者点热源来简化焊接电弧热作用； b）二维传热（忽略深度方向的热传导）； c）材料导热系数不随温度变化； d）传热过程为准稳态； e）忽略对流与辐射对于温度场的影响； f）忽略相变过程中的热效应。,第二部分热过程的数值模拟,但是目前仍然无法完全对复杂结构的传热
4、过程进行解析计算，这样就要借助于其他的一些手段来分析焊接温度场，目前比较常用的数值计算方法有两种：有限单元法和有限差分法。但是相对于有限差分法，有限元方法具有如下的优势： a）可以模拟任意形状复杂的结构； b）比较容易的实现在温度场结果基础上进行应力场的计算； c）可以处理诸如接触在内的复杂边界条件的传热问题。,第二部分热过程的数值模拟 均质棒材的一维稳态导热分析,将一个棒材除了端面外的其余部分均设置为绝热状态，棒材的横截面面积为A，棒材长度为L。对于此问题的有限元建模可以简化为两个节点1和2连接的线性单元。,周围绝热棒材的一维导热分析示意图及有限元模型,第二部分热过程的数值模拟 均质棒材的一
5、维稳态导热分析,根据傅立叶定律，对于一维传热问题，热流可以表达为如下的公式Q为单位时间通过横截面A的热量，T为温度，x为热流传播方向的坐标，为导热材料的导热系数。,第二部分热过程的数值模拟 均质棒材的一维稳态导热分析,稳态传热,第二部分热过程的数值模拟 均质棒材的一维稳态导热分析,第二部分热过程的数值模拟 复合棒材中的热传导有限元分析,第二部分热过程的数值模拟 复合棒材中的热传导有限元分析,第二部分热过程的数值模拟 复合棒材中的热传导有限元分析,第二部分热过程的数值模拟 复合棒材中的热传导有限元分析,第二部分热过程的数值模拟 加权余量法,以上的结果是采用直接法进行的有限元分析，但是直接法只能用
6、来推导比较简单的有限元方程。例如假设温度、位移是线性变化的，同时在单元边界上热流、应力、表面力是常数，容易转化成等效的节点热流和节点力。直接法形式上把连续区域化分为有限元网格，对每个有限元用直接法分析得到单元刚度矩阵再组合成总体刚度矩阵。然而对于比较复杂的系统，以上的方法实现起来就比较困难。 事实上，大部分工程实际问题都受控于该问题的微分方程。由于复杂的几何形状以及载荷状态，很多工程实际问题的微分方程都无法直接获得其精确解。因而在工程领域中就会使用很多的近似方法来求解微分方程。有限元方法就是其中一种近似方法，但实际上有限元方法也要基于一些其他的更基础的近似方法，其中加权余量法是一种应用十分普遍
7、的近似方法。,第二部分热过程的数值模拟 加权余量法,这里以一维的问题入手，介绍加权余量法的基本思想及基本分析过程，实际上基于这样的思想，也同样可以实现二维或者三维问题的分析。讨论如下的微分形式,第二部分热过程的数值模拟 加权余量法,而是通过将y*带入微分方程所获得其残余误差（即余量）,第二部分热过程的数值模拟 加权余量法,所以加权余量法是一种解特定边界条件的近似方法，它应用特定的试探函数使之满足相应的边界条件，同时使余量的加权积分为零来求解微分方程的近似解。 这里的权函数是可以任意选择的，目前依据选择的权函数不同，存在多种比较成型的加权余量法。这里仅对有限元中经常采用的伽辽金法进行介绍。,第二
8、部分热过程的数值模拟 加权余量法,伽辽金方法中直接采用试探函数作为权函数。这样原来微分方程的加权余值的积分变为以下的形式,第二部分热过程的数值模拟 加权余量法,获得了n个方程，这样可以解出近似解中的n个未知量ci。从而获得微分方程的近似解。,第二部分热过程的数值模拟 加权余量法,问题：应用伽辽金方法，解以下微分方程的近似解：,边界条件为y(0)=y(1)=0，即为协调边界条件。,第二部分热过程的数值模拟 加权余量法,由于微分方程中存在二次项，所以这里采用抛物线关系的试探函数即能满足要求，一般来说，如果在区间a b上的同时为了满足边界条件，可以选择以下的多项式形式,第二部分热过程的数值模拟 加权
9、余量法,应用以上的试探函数，则近似解可以表达为：将这个试探函数带入微分方程中，则余量的表达式为：,第二部分热过程的数值模拟 加权余量法,解出c14。所以近似解为,第二部分热过程的数值模拟 加权余量法,精确解：带入边界条件,第二部分热过程的数值模拟 加权余量法,第二部分热过程的数值模拟 加权余量法,应用两项伽辽金近似解以上的微分方程，试探函数分别为：两项近似解的形式为,第二部分热过程的数值模拟 加权余量法,解得：,第二部分热过程的数值模拟 加权余量法,以上微分方程的两项近似解为,第二部分热过程的数值模拟 伽辽金法,使用伽辽金法进行整个问题求解区域内的近似解计算通常会带来一些问题： 一是当问题比较
10、复杂的时候，需要使用项数非常多的近似解，才能获得比较接近的计算结果，这个计算过程非常复杂，计算周期非常长，尤其是在进行二维或者三维分析的时候还涉及到偏微分方程的问题，那时候获得合适的试探函数将是一件非常困难的事情； 另外近似解计算的结果可能满足所有边界条件，但是会在其他的位置出现近似解与精确解的结果相震荡等问题。为了解决这个问题，可以采用将整个求解区域划分成有限元网格的方法，此时每个单元的积分区域都会变得比较小，这样可以避免以上的两个问题。 这就提出了有限元的伽辽金方法。,第二部分热过程的数值模拟 伽辽金法,问题：应用伽辽金方法，解以下微分方程的近似解：,边界条件为y(a)=ya， y(b)=
11、yb即为协调边界条件。 这里不对整个的分析区域即进行直接求解，而是将这个范围划分成M个单元进行分析,第二部分热过程的数值模拟 伽辽金法,第二部分热过程的数值模拟 伽辽金法,值得注意的是此时积分的上下限已经不再是整个积分范围，而是每个单元的区域范围。而被积分的对象也转变为与单元相应的积分对象。从上面的公式可以获得M+1个代数方程，从而可以解出M+1个节点的函数值yi。这些代数方程可以写成矩阵形式：,第二部分热过程的数值模拟 伽辽金法,将精确解、2个单元计算结果以及4个单元得计算结果绘于图中。两个单元得结果与精确解之间差别较大，仅仅在节点的函数值方面与精确解接近，而其在节点处导数的连续性比较差。四
12、个单元的计算结果与精确解之间的差别很小，更加接近精确解，但是节点处的导数不连续现象依旧比较明显，可以预期如果进一步增加单元的数目，计算解果将更加逼近真实解，同时节点的导数不连续现象能够得到进一步缓解。,第二部分热过程的数值模拟 一维稳态传热问题的伽辽金有限元解法,本节应用以上提到的伽辽金有限元法来进行传热问题的有限元分析。为了便于分析，这里就一维稳态传热问题的有限元分析过程进行介绍。,热量仅沿着x轴方向传导，而其他的与x轴垂直的方向假定为完全的绝热边界条件，即这些面上没有热量损失。,第二部分热过程的数值模拟 一维稳态传热问题的伽辽金有限元解法,在dt时间内，以上公式对应的方程可写为：,第二部分
13、热过程的数值模拟 一维稳态传热问题的伽辽金有限元解法,假设材料的导热系数为常数 对于稳态传热进行分析，即,第二部分热过程的数值模拟 一维稳态传热问题的伽辽金有限元解法,应用加辽金方法有限单元的方法来解以上的微分方程。这里的单元类型为两节点线性插值单元，则通过上一节分析，其插值函数的表达式为： 其中T1和T2分别代表单元上的两个节点的温度值。令微分方程的加权余量的积分为零，即：,i=1,2,第二部分热过程的数值模拟 一维稳态传热问题的伽辽金有限元解法,对于积分中的第二项进行分部积分，则以上公式转变为将插值函数及相应的Ni带入上式，与上一节的分析类似，可以获得如下两个公式：,i=1,2,第二部分热
14、过程的数值模拟 一维稳态传热问题的伽辽金有限元解法,经计算,第二部分热过程的数值模拟 一维稳态传热问题的伽辽金有限元解法,圆棒的直径为60mm，长1m，圆棒周向绝热，其左半部为铝，导热系数为200W/(m)，右半部分为铜，其导热系数为389W/(m)。铝棒左端的热流为4000W/m2；铜棒右端温度一直保持在80。应用4个等长度单元分析稳态条件下棒材的温度分布。这个例子我们可以简单的采用本部分刚刚开始的时提到的直接法进行求解，但是这里主要为了说明伽辽金法的工作原理，所以采用伽辽金法进行求解。以期读者掌握这种方法热分析的基本流程。,第二部分热过程的数值模拟 一维稳态传热问题的伽辽金有限元解法,对于
15、铝相关的单元1和单元2，其传热矩阵为：,W/,第二部分热过程的数值模拟 一维稳态传热问题的伽辽金有限元解法,铜相关的单元3和单元4，其传热矩阵为：将边界条件q14000W/m2，T580，带入到整体矩阵方程中去，则,W/,第二部分热过程的数值模拟 一维稳态传热问题的伽辽金有限元解法,将T5带入到前四个方程中，则前四个方程变成如下的矩阵方程 这里四个方程，四个未知数，最终解得：T195.15、T290.14、T385.15、T482.57。将T4和T5带入原矩阵方程的第五个方程，解得q54038.6W/m2。,第二部分热过程的数值模拟 一维稳态传热问题的伽辽金有限元解法,从以上的结果可以看出传入
16、圆棒的热流q14000W/m2与传出圆棒的热流q54038.6W/m2之间存在一定的差距，而实际上，他们应该相等，这主要是由于这里采用了手工计算的形式，中间数据的有效数字不够，如果采用高精度的计算，则可以近似获得q54000W/m2。,第三部分 三维焊接结构应力变形的有限元法,第三部分 三维焊接结构应力变形的有限元法,焊接过程中的应力与变形的有限元分析的实质是热应力的分析。目前关于热应力的分析一般有两种解法，其中一种是解耦算法，另外一种是热机耦合算法。所谓解耦算法就是仅考虑温度对于结构应变、应力的作用，而不考虑由于结构变形而产生的热效应；而热机耦合算法则在计算中同时考虑这种相互作用。实际上对于
17、焊接过程而言，解耦算法就能满足计算精度要求。所以本节主要介绍解耦算法的基本思想。,第三部分 三维焊接结构应力变形的有限元法,解耦算法的基本思想是先进行结构的温度场的计算，然后以温度场为边界条件计算应力场。,第三部分 三维焊接结构应力变形的有限元法,当物体各部分的温度发生变化时，物体由于热胀冷缩（通常的材料符合这个规律，在此这样说，是希望容易理解）而产生膨胀或收缩，从而产生线应变：当物体的热变形不受任何条件的约束，则物体不会产生热应力，但是当物体受到外界条件的约束、或者物体受热不均匀各部分存在一定的温度差时，各部分的热应变之间相互作用而产生应力。,第三部分 三维焊接结构应力变形的有限元法,物体由
18、于热膨胀的作用只产生线性应变，而剪应变为零。 这种由热应变产生的应变可以看作物体的自由应变,第三部分 三维焊接结构应力变形的有限元法,其中的T可以通过前面提到的温度场的计算中的方法求得。在获得了自由应变之后，可以根据外观变形、自由变形和内部变形之间的关系，确定内部应变的表达式为：外观应变与具体的坐标的位置相关，就是我们前面进行应力应变分析中的应变，其为节点坐标的函数。即：,第三部分 三维焊接结构应力变形的有限元法,而结构的应力则与结构的内部应变相对应，存在如下的关系 通过虚功原理，我们可以建立起来节点位移与节点力之间的关系，即 而ff是体积力及表面力等作用所产生的节点力。,第三部分 三维焊接结构应力变形的有限元法,由以上的公式可以看出，结构的热应力问题和无热载荷的应力分析相比，除了增加一项初始温度应变所产生的温度载荷以外，其他的完全相同，这样对于稳态温度场的应力问题，可采用关于三维应力分析的有限元计算方法进行计算分析。而对于瞬态温度所产生的应力场，可以在每一步瞬态温度场计算后进行相应的应力场计算。亦可在所有瞬态温度场计算完成后，对于所有时间步或者特定的时间步的应力场进行计算。 针对于焊接应力场的计算，如果在获得温度场的计算的情况下，是可以利用这里提到的分析方法进行相应的应力场计算的。对于焊接变形问题，由于我们获得了所有节点的位移情况，所以可以方便的获得结构的变形情况。,
