数值分析小论文

基于ABAQUS软件的混凝土柱的有限元分析摘要：有限元法是工程分析中广泛应用的数值计算方法，由于它的通用性和有效性，受到工程技术界的高度重视。ABAQUS软件是国际上公认的最好的CAE大型通用分析软件之一。本文对有限单元法进行简单介绍并采用ABAQUS软件分析一混凝土柱的受力问题。关键词：ABAQUS，混凝土柱，有限元分析1有限元理论概述1.1有限元法基本思想有限元法的基本思想是将连续的求解区域离散为一组有限个、且按一定方式相互联结在一起的单元组合体。由于单元能按不同的联结方式进行组合，且单元本身可以有不同形状，因此可以模型化几何形状复杂的求解区域。有限元法作为数值分析方法的一个重要特点是利用在每一个单元内假设的近似函数，分片地表示全求解域上待求的未知场函数，单元内的近似函数通常由未知场函数或其导数在单元的各个节点的数值和其插值函数表达。这样，一个问题的有限元分析中，未知场函数或其导数在各个节点上的数值就成为新的未知量（即自由度），从而使一个连续的无限自由度问题变成离散的有限自由度问题。一经求解出这些未知量，就可通过插值函数计算出各个单元内场函数的近似值，从而得到整个求解域上的近似解。显然，随着单元数目的增加，即单元尺寸的缩小，或者随着单元自由度的增加及插值函数精度的提高，解的近似程度将不断改进，如果单元是满足收敛要求的，近似解最后将收敛于精确解。1.2有限元法分类1.2.1线弹性有限元法线弹性有限元法以理想弹性体为研究对象，所考虑的变形建立在小变形假设的基础上。在这类问题中，材料的应力与应变呈线性关系，满足广义胡克定律；应变与位移也是线性关系。线弹性有限元问题归结为求解线性方程组问题，所以只需要较少的计算时间。如果采用高效的代数方程组求解方法，也有助于降低有限元分析的时间。线弹性有限元一般包括线弹性静力分析与线弹性动力分析两个主要内容。学习这些内容需具备材料力学、弹性力学、结构力学、数值方法、矩阵代数、算法语言、振动力学、弹性动力学等方面的知识。1.2.2非线性有限元法非线性有限元问题与线弹性有限元问题有很大不同，主要表现在如下三个方面：①非线性问题的方程是非线性的，因此一般需要迭代求解；②非线性问题不能采用叠加原理；③非线性问题不总有一致解，有时甚至没有解。以上三方面的因素使非线性问题的求解过程比线弹性问题更加复杂、费用更高和更具有不可预知性。有限元法所求解的非线性问题可以分为如下三类：（1）材料非线性问题材料的应力与应变是非线性关系，但当应变与位移很微小时，可以认为应变与位移呈线性关系，这类问题属于材料非线性问题。由于从理论上还不能提供能普遍接受的本构关系，所以，一般来说，材料的应力与应变之间的非线性关系要基于试验数据，有时非线性特性可用数学模型进行模拟，尽管这些模型总是有它们的局限性。在工程实际中较为重要的材料非线性问题有：非线性（包括分段线弹性）、弹塑性、粘塑性及蠕变等。（2）几何非线性问题几何非线性是由于位移之间存在非线性关系引起的。当物体的位移较大时，应变与位移的关系是非线性关系，这意味着结构本身会产生大位移或大转动，而单元中的应变却可大可小。研究这类问题时一般都假定材料的应力与应变呈线性关系。这类问题包括大位移大应变问题及大位移小应变问题。如结构的弹性屈曲问题属于大位移小应变问题，橡胶部件形成过程为大应变问题。（3）非线性边界（接触问题）在加工、密封、撞击等问题中，接触和摩擦的作用不可忽视，接触边界属于高度非线性边界。平时遇到一些接触问题，如齿轮传动、冲压成型、轧制成型、橡胶减振器、紧配合装配等，当一个结构与另一个结构或外部边界相接触时通常要考虑非线性边界条件。实际的非线性可能同时出现上述两种或三种非线性问题。1.3有限元分析系统的组成有限元法的通用性使得它可以把固体力学、流体力学、动力学与控制等不同分支中课题的求解统一在一个框架，组织在一个分析系统中。基于数理模型，有限元分析系统一般由前处理器、模型求解器、后处理器三个部分组成。其中前后处理器是算法与空间模型的接口，进行相应数据的前期准备与后期整理，完成算式表达和结果显示。模型求解部分实现数理方程的解算。对线性化模型，目前算法己近于成熟；当前数理方法的主要研究方向是非线性问题和多体系统建模。1.4有限元法的求解过程对于不同物理性质和数学模型的问题，有限元求解法的基本步骤是相同的，只是具体公式推导和运算求解不同。有限元请求解问题的基本步骤通常为：（1）问题及求解域定义。根据实际问题近似确定求解域的物理性质和几何区域。（2）求解域离散化。将求解域近似为具有不同大小和形状且彼此相连的有限个单元组成的离散域，习惯上称为有限元网络划分。显然单元越小、网络越细、则离散域的近似程度越好，计算结果也越精确。但是计算量及误差都将增大，因此求解域的离散化是有限元法的核心技术之一。（3）确定状态变量及控制方法。一个具体的物理问题通常可以用一组包含问题状态变量边界条件的微分方程式表示。为适合有限元求解，通常将微分方程化为等价的泛函数形式。（4）单元推导。对单元构造一个适合的近似解，即推导有限单元的列式。其中包括选择合理的单元坐标系，建立单元式函数，以某种方法给出单元各状态变量的离散关系，从而形成单元矩阵（结构力学中称刚度阵或柔度阵）。为保证问题求解的收敛性，单元推导有许多原则要遵循，对工程应用而言，重要的是应注意每一种单元的解题性能与约束。例如，单元形状应以规则为好，畸形时不仅精度低，而且有缺值的危险，将导致无法求解。（5）总装求解。将单元总装形成离散域的总矩阵方程（联合方程组），反映对近似求解域的离散域的要求，即单元函数的连续性要满足一定的连续条件。总装是在相邻单元结点进行，状态变量及其导数（可能的话）连续性建立在结点处。（6）联立方程组求解和结果解释。有限元法最终导致联立方程组。联立方程组的求解可用直接法、选代法、和随机法，求解结果是单元结点处状态变量的近似值。对于计算结果的质量，将通过与设计准则提供的允许值比较来评价并确定是否需要重复计算。2问题描述ABAQUS是世界上最先进的大型通用有限元分析软件之一，它具有丰富的材料本构模型和平易近人的开发平台。本文利用ABAQUS软件进行柱的有限元分析，模拟柱在位移荷载下的最终破坏形态，进行应力分析。利用ABAQUS有限元软件分析一个高2m，截面尺寸为的对称配筋的钢筋混凝土短柱的受力。混凝土柱配置4Φ22受压钢筋，箍筋为Φ8@200，保护层厚度为30mm，其立面尺寸如图2.1所示，截面尺寸如图2.2所示。图2.1柱立面尺寸图2.2柱截面尺寸材料特性：混凝土：抗压强度，抗拉强度钢筋：弹性模量为，泊松比，屈服强度3建模过程3.1创建部件在ABAQUS/CAE窗口顶部的环境栏中，进入part模块。创建混凝土部件使用此工具，进行混凝土部件创建。在Name后输入column，将ModelingSpace设置为3D，Shape设为solid，Type选择Extrusion，Approximatesize项后输入1，其他参数取默认值，点击continue进入二维绘图界面。将柱的二维几何模型绘出之后并在弹出的EditBaseExtrusion对话框，在Depth项后输入2，完成三维模型的建立。如图3.1所示。图3.1column部件的三维模型创建箍筋部件同样使用此工具，进行箍筋部件创建。在Name后输入gujin，将ModelingSpace设置为3D，Shape设为wire，其他参数取默认值，点击continue进入二维绘图界面。将柱的二维几何模型绘出之后点击鼠标中键完成三维模型的建立。如图3.2所示。创建受压筋部件与创建箍筋方法相似，将Name定义为shouyajin，生成三维模型如图3.3所示。图3.2gujin部件三维模型图3.3受压筋部件三维模型3.2创建材料和截面属性在环境栏的Module列表中选择Property功能模块。1、创建材料混凝土的本构关系点击左侧工具区的，弹出EditMaterial对话框，在Name后输入concrete。在对话框中选择General的Density选项，在MassDensity后输入2400，表明混凝土的密度为2400kg/m3。点击Mechanical→Elasticity→Elastic，在数据表中的Young’sModulus为2.9e10，Poisson’sRatio为0.2。点击Plasticity→ConcreteDamagedPlasticity，填入数据如图3.4所示。在CompressiveBehavior中填入表3.1的数据；在TensileBehavior填入表3.2的数据。完成混凝土本构模型的建立。YieldStressCrackingStrain2400000019600000.000115000000.000312900000.000410800000.00057230000.00084840000.0012170000.002984000.003539000.005YieldStressInelasticStrain182000000221000000.0004240000000.0008245000000.0012230000000.002216000000.0024166000000.0036113000000.00522300000.01表3.1压缩特性数据表3.1压缩特性数据表3.2拉伸特性数据图3.4混凝土塑性损伤塑性模型数据钢筋本构模型钢筋本构模型的建立过程与混凝土本构模型的建立过程相似，但注意，钢筋的质量密度为7800.弹性模量为1.9e11，泊松比为0.3，塑性屈服应力为2.1e8，泊松应变为0。2、创建截面属性点击左侧工具区的，创建Section-1的截面属性，Category设为solid。Type设为Homogeneous，保持其余参数不变继续下一步，在Material对话框中选择concrete，保持其余参数不变，完成Section-1的截面属性创建。点击左侧工具区的，创建箍筋的截面属性，在Name后面输入gujin，Category设为beam。Type设为Truss，保持其余参数不变继续下一步，在Material对话框中选择gangjin，由于箍筋的直径为8mm，因此在Cross-sectional后填写50.3e-6，保持其余参数不变，完成箍筋的截面属性创建。以同样的方式创建受压筋的截面属性，Name为shouyajin，由于受压筋的直径为22mm，因此Cross-sectional后填入3.8e-4。给部件赋予截面属性使用左侧工具区的，分别对混凝土、箍筋和受压筋赋予截面属性。赋予截面属性后，各个部件的颜色变成青色。3.3定义装配件在窗口环境栏的Module列表中选择Assembly功能模块。利用旋转、阵列和移动工具的操作方法，将箍筋以及受压筋组合成一个钢筋骨架，并将column部件消隐，如图3.5所示。同时为了方便操作，使用左侧工具区的，将整个钢筋骨架合成一个整体，并命名为gujia。最终将骨架沿Z轴正方向移动0.03，是钢筋与混凝土有确定的相对位置，最终结果如图3.6所示。图3.5混凝土柱的钢筋骨架图3.6装配件三维视图3.4设置分析步在窗口环境栏的Module列表中选择Step功能模块。点击左侧工具区的，在Proceduretype项中，选择General，下拉菜单中选择Static，General，继续下一步，将Basic中的Timeperiod改为50，Incrementation标签页的数据如图3.7所示。图3.7EditStep对话框3.5定义约束在窗口环境栏的Module列表中选择Interaction功能模块进行模型之间的约束关系。点击左侧区的，在CreteConstraint界面的name后输入gangjin，Type选择EmbededRegion，将gujia作为嵌入的部分，整体模型作为主区域的方法从而完成钢筋骨架与混凝土约束关系的定义，如图3.8所示。点击左侧工具区的，提示区显示“Selectastartpointtoactasreferencepoint-orenterX,Y,Z”,输入（0,0，-0.05）建立第一个参考点，输入（0,0,2.05）建立第二个参考点，如图3.9所示。图3.8定义约束后的模型图3.9建立参考点后的模型点击左侧区的，在CreteConstraint界面的name后输入RP-1，Type选择Coupling，将RP-1与柱的下底面耦合在一起，用同样的方法将RP-2与柱的上帝么耦合在一起，如图3.10所示。图3.10定义参考点耦合后的模型4计算条件计算条件可分为边界条件约束、施加荷载以及网格划分。4.1定义边界条件在窗口环境栏的Module列表中选择Load功能模块。点击左侧区的，在弹出来的定义边界条件的对话框中，Step选为initial，Typesforselectedsteps选择为Displacement/Rotation，将RP-1作为边界条件的区域，在定义边界条件是选中所有的选项，即添加固定约束作为边界条件，从而完成模型边界条件的定义，如图4-1所示。图4-1创建约束后的模型4.2施加荷载本算例采用位移加载的方式施加荷载，利用Amplitudes建立加载的规律，进行加载。建立加载规律选择主菜单Tools→Amplitudes→Create，在弹出来的对话框中将Name定义为Amp-1，Type选择Tabular,即以表格的形式输入振幅变化规律。继续进行振幅的定义，如图4-2所示，将4-3所示的数据整理成Excel表格文件，复制粘贴到4-2所示的表格中即可，完成加载规律的定义。施加荷载点击左侧区的，在弹出来的定义边界条件的对话框中，Step选为step-1，Typesforselectedsteps仍选择为Displacement/Rotation，将RP-2作为边界条件的区域，EditBoundaryCondition对话框填入数据如图4-4所示，即在RP-2处施加位移荷载，从而完成模型荷载的定义，如图4-5所示。位移时间位移时间图4-2EditAmplitude对话框图4-3位移加载曲线图4-4EditBoundaryCondition对话框图4-5定义边界条件和荷载后的模型4.3网格划分在窗口环境栏的Module列表中选择Mesh功能模块进行网格划分。（1）布置边上种子在进行网格划分时要确保整个模型为independent，否则，网格划分进行不了。点击左侧工具区的，在GlobalSeeds对话框中的Approximateglobalsize中输入0.05，表明每隔0.05m布置一个种子。布置完种子的模型如图4-6所示。划分网格点击左侧工具区中，将模型按照网格种子自动划分网格，如图4-7所示。图4-6模型网格种子分布图图4-7划分网格后的装配件查看单元类型由于杆单元默认为beam单元，因此需要修改gujia部件单元类型，否则模型计算时会认为没有对骨架部件进行材料属性定义，导致错误产生。点击左侧工具区，将混凝土部件隐藏，选中骨架部件，将元素类型改为Truss单元，如图4-8所示，完成单元类型的修改。图4-8ElementType对话框4.4提交分析作业在窗口环境栏的Module列表中选择Job功能模块进行作业提交。创建分析作业利用左侧工具中的进行工作创建。提交分析选择主菜单Job→Manager，在JobManager对话框中，点击Submit进行作业提交。可以看到对话框中的状态提示依次变为Submitted、Running和Completed，由此计算完成。5提取结果在Visualization模块中，进行结果提取。5.1变形图图5-1表明模型在位移荷载的作用下的变形，图5-2显示了钢筋骨架变形图，图5-3显示了混凝土的变形图。图5-1变形后的模型图5-2钢筋骨架变形图图5-3混凝土变形图由图5-1至5-3可以看出，此柱子在底部为固定约束，顶部没有约束但受随着时间呈线性增长位移荷载的作用下，变形沿着柱中截面呈马鞍形状的对称分布，最大变形主要出现在柱子的中部。但在柱中间变形不是最大，相比在柱中截面上下5cm处变形最大。由变形图可以看出在柱子的中部是易于破坏区，需要进行加固处理。但在实际工程中，柱的节点处钢筋需要加密或者通常加密，这是由于实际工程中最先考虑“强柱弱梁，强节点弱构件”的概念，且梁、板、柱的节点为易于破坏区，所以整体分析与单个柱子分析不尽相同。5.2云纹图图5-4变形后模型的云纹图图5-5钢筋骨架变形后的云纹图图5-6混凝土变形后的云纹图从图5-4至5-6可以看出，该混凝土柱的Mises屈服应力沿着柱的中截面处程对称分布。整体模型的最大应力为2.1×108pa，最小应力为3.62×104pa；钢筋骨架的最大应力为2.1×108pa，最大应力出现在骨架中部，最小应力为2.445×106pa，最小应力出现在上下底面最外层箍筋的中部；混凝土的最大应力为1.7×107pa，最大应力在混凝土柱的上下底面的四个角节点处，最小应力为3.62×104pa，最小应力出现在柱中截面附近。钢筋骨架部分变形越大，应力越大，而混凝土部分变形最大的地方，应力却最小。5.3混凝土部分切片图通过绘制混凝土切片下的云纹图，可以更清楚的看到混凝土在各个部位的Mises应力分布情况。由图5-7至5-9可以看出沿3个主轴方向的截面应力分布较小。由于荷载施加在混凝土柱的顶部，从图5-7和5-8可以清晰的看出顶部有微裂纹的出现。图5-7混凝土部分沿X轴所得切片图5-8混凝土部分沿Y轴所得切片图5-9混凝土部分沿Y轴所得切片6参数分析及其影响规律6.1混凝土部分应力分析在Visualization模块中，可以进行参数的提取及其分析。以下应力、应变分析中所选的混凝土部分中央节点区域如图6-1所示。图6-2为混凝土部分中央节点区域的Mises应力随时间的变化图。从图中可以看出，在6s之前，应力随时间先呈线性增长趋势，表明此时在混凝土中央节点区域处于弹性阶段。随着时间的推移，该区域进入塑性区域，最大屈服应力达到24.85MPa。随后应力随着时间的增加呈下降趋势，表明该节点部分进入破坏阶段。图6-1混凝土部分中央节点区域所在位置图6-2混凝土中央节点Mises应力与时间的关系曲线图6-3混凝土中央节点应力S33与时间的关系曲线图6-2表明混凝土中央节点应力S33与时间的关系曲线。应力S33指正应力的大小，此应力为负值，表明为压应力。由图6-2可以看出有着与Mises屈服应力相似的应力增长过程。一开始处于线性增长阶段，也处于弹性阶段，随着时间的推移，进入塑性阶段直至混凝土产生裂缝发生破坏。从图中以及ABAQUS的表格中可以得到最大值为-25.5MPa，大于Mises屈服应力。图6-4混凝土中央节点应力S11与时间的关系曲线图6-2表明混凝土中央节点应力S11与时间的关系曲线。应力S11指正应力的大小，此应力为负值，表明为压应力。在22.3s之前，随着时间的增长也逐渐增大，最大应力值为-1.0624MPa。随后又随着时间的增加而逐渐减小至应力为0并逐渐变为拉应力，最大拉应力出现在36.3s，拉应力的值为24267.8Pa。随后拉应力减小为0再次变为压应力，压应力逐渐增大。由图可以看出的变化有明确的变化规律，应力变化值较小，但是出于压、拉应力的变换过程。在ABAQUS中，可以采用同样的方法，将、、、、、以及Tresca屈服应力、静水压力等等应力通过图表的形式表达出。此处就不一一列举。6.2钢筋骨架部分应力分析以下钢筋骨架部分应力、应变分析中，箍筋所选分析部位如图6-5所示，受压筋所选分析部位如图6-6所示。图6-5箍筋中央节点区域所在位置图6-6受压筋分析区域所在位置6.2.1箍筋部分应力分析在Visualization模块中，同样可以进行钢筋骨架部分参数的提取及其分析。图6-7为柱中央区域箍筋Mises屈服应力与试件的关系曲线。从该曲线中可以看出，在20s之前，箍筋应力随着时间的增加而增长，但并没有呈现出线性增长。当20s之后，箍筋应力达到极限，并保持不变，最大应力值为21MPa。6.2.2受压筋部分应力分析图6-8为柱中央区域受压筋Mises屈服应力与试件的关系曲线。从该曲线中可以看出，在10.5s之前，受压筋应力随着时间的增加而增长，并且呈现出线性增长的趋势。当10.5s之后，受压筋应力达到极限，并基本趋于稳定，最大应力值为21MPa。与箍筋的Mises屈服应力相比，受压筋进入塑性阶段时间较早，但最终屈服应力相同，这是由于箍筋和受压筋均采用HRB300的钢筋的缘故。图6-7混凝土柱中央区域箍筋Mises屈服应力与试件的关系曲线图6-8混凝土柱中央区域箍筋Mises屈服应力与试件的关系曲线6.3混凝土部分应变分析图6-9混凝土中央节点应变与时间的关系曲线由图6-9可以看出混凝土中央节点的应变在加载的50s内一直处于增长状态，在前30s应变增长较缓和，后20s增长较迅速，在完成最终加载后最大应变为0.01