焊接模拟ansys实例(2020年整理).doc
- 格式:doc
- 大小:43.50 KB
- 文档页数:7
《基于ANSYS的焊接温度场和应力的数值模拟研究》篇一一、引言随着制造业和机械工程的不断发展,焊接技术已成为一种重要的连接方法,广泛应用于各个工业领域。
然而,焊接过程中涉及到的高温、高应力等复杂因素往往会对焊接件的性能产生影响。
因此,对焊接过程中的温度场和应力进行准确预测和模拟,对于优化焊接工艺、提高产品质量具有重要意义。
本文基于ANSYS 软件,对焊接过程中的温度场和应力进行数值模拟研究。
二、研究背景及意义焊接过程中,温度场和应力的分布直接影响着焊接件的质量和性能。
传统的焊接工艺往往依赖于经验公式和试验数据,但这些方法往往难以准确预测焊接过程中的温度场和应力分布。
因此,采用数值模拟方法对焊接过程进行模拟,可以为优化焊接工艺、提高产品质量提供有力支持。
ANSYS作为一种功能强大的工程仿真软件,广泛应用于各种工程领域的数值模拟,包括焊接过程的温度场和应力模拟。
三、研究方法本研究采用ANSYS软件,对焊接过程中的温度场和应力进行数值模拟。
首先,建立焊接过程的有限元模型,包括材料属性、热源模型、边界条件等。
其次,利用ANSYS的热分析模块,对焊接过程中的温度场进行模拟,分析温度分布和变化规律。
最后,利用ANSYS的结构分析模块,对焊接过程中的应力进行模拟,分析应力的分布和变化规律。
四、结果与讨论1. 温度场模拟结果通过ANSYS的热分析模块,我们得到了焊接过程中的温度场分布。
结果表明,在焊接过程中,焊缝处的温度最高,随着距离焊缝的远离,温度逐渐降低。
此外,我们还发现在焊接过程中,温度分布呈现出明显的非线性特征,这与实际焊接过程中的情况相符。
2. 应力模拟结果在得到温度场分布的基础上,我们进一步利用ANSYS的结构分析模块对焊接过程中的应力进行模拟。
结果表明,在焊接过程中,焊缝处由于高温产生的热应力较大。
此外,由于焊接过程中材料的热膨胀和收缩不均匀,还会产生残余应力。
这些应力对焊接件的性能和寿命具有重要影响。
第一章前处理第1例 关键点和线的创建实例—正弦曲线FINISH/CLEAR, NOSTART /PREP7K,100,0,0,0CIRCLE,100,1,,,90 CSYS,1KFILL,2,1,4,3,1K,7,1+3.1415926/2,0,0 CSYS,0KFILL,7,1,4,8,1 KGEN,2,7,11,1,,1 LSTR,8,13 LSTR,9,14 LSTR,10,15 LSTR,11,16 LANG,5,6,90,,0 LANG,4,5,90,,0 LANG,3,4,90,,0 LANG,2,3,90,,0BSPLIN,1,17,18,19,20,12 LSEL,U,,,14LDELE,ALL LSEL,ALL KWPAVE,12 CSYS,4 LSYMM,X,14NUMMRG,KP,,,,LOWLCOMB,ALL,,0FINISH/CLEAR, NOSTART /PREP7 PI=3.14159 J=0*DO,I,0,PI,PI/10.0 J=J+1 X=IY=SIN(I) I=I+1 K,J,X,Y *ENDDOBSPLIN,1,2,3,4,5,6 BSPLIN,6,7,8,9,10,11 csys,4 KWPAVE,11LSYMM,y,1,2,,,,0 KWPAVE,11LSYMM,x,3,4,,,,1以上程序有意没算到2 为了使用几个命令第2例工作平面的应用实例—相交圆柱体[本例提示]通过相交圆柱体的创建,本例主要介绍了工作平面的使用方法。
通过本例,读者可以了解并掌握工作平面与所创建体的位置、方向的关系,学习工作平面的设置、偏移、旋转和激活为当前坐标系的方法。
FINISH/CLEAR,NOSTART/PREP7CYLIND,0.015,0,0,0.08,0,360CYLIND,0.03,0,0,0.08,0,360/VIEW,1,1,1,1/PNUM,VOLU,1WPOFF,0,0.05,0.03WPROT,0,60CYLIND,0.012,0,0,0.055,0,360CYLIND,0.006,0,0,0.055,0,360VSEL,S,,,2,3,1CM,VV1,VOLUVSEL,INVECM,VV2,VOLUVSEL,ALLVSBV,VV1,VV2BLOCK,-0.002,0.002,-0.013,-0.009,0,0.008WPSTYLE,,,,,,1CSYS,4VGEN,3,1,,,,120VSBV,5,1VSBV,4,2VSBV,1,3WPROT,0,0,90VSBW,ALLVDELE,1,4,3VADD,ALLVPLOT/REPLOT第3例复杂形状实体的创建实例—螺栓[本例提示]在使用ANSYS软件进行结构分析时,建立实体模型是最复杂最难以掌握的一个过程。
基于ANSYS的某焊接件两焊缝在顺序焊接过程中的分析(生死单元应用案例)焊接几何模型如下图所示,左右两侧90度扇区为焊接材料,其余为钢板材料。
其他更多已知条件请参考命令流,这里不再赘述。
网格单元本实例中顺序焊接分为如下步骤:第一步0-1秒:右侧焊接稳态分析(杀死左焊缝,施加右焊缝温度和焊接件参考温度)第二步1-100秒:相变分析(删除温度载荷,施加对流热传导)第三步100-1000秒:右侧焊缝凝固分析第四步1000-1001秒:激活左侧焊缝单元进行稳态分析(施加左焊缝温度)第五步1001-1100秒:左焊缝相变分析第六步1100-2000秒:左侧焊缝凝固分析第七步:结果后处理ANSYS命令流:FINISH/FILNAME,Exercise ! 定义隐式热分析文件名/PREP7 ! 进入前处理器ET,1,SOLID70 ! 选择8节点实体热分析单元MP,KXX,1,.5e-3MP,C,1,.2MP,DENS,1,.2833MPTEMP,1,0,2643,2750,2875,3000MPDATA,ENTH,1,1,0,128.1,163.8,174.2,184.6 ! 定义右焊缝材料热物理性能MP,KXX,2,.5e-3MP,C,2,.2MP,DENS,2,.2833MP,KXX,3,0.5e-3 ! 定义两块钢板的热物理性能MP,DENS,3,.2833MPTEMP,1,0,2643,2750,2875,3000MPDATA,ENTH,3,1,0,128.1,163.8,174.2,184.6 ! 定义左焊缝材料热物理性能BLOCK,-0.17,0.17,0,0.34,0,1.2BLOCK,0.17,0.34,0,0.34,0,1.2BLOCK,0.34,1,0,0.34,0,1.2BLOCK,-0.17,0.17,0.34,0.51,0,1.2BLOCK,-0.17,0.17,0.51,1.34,0,1.2WPAVE,0.17,0.34,0CYLIND,0.17,0,0,1.2,0,90WPAVE,0,0,0CSYS,0FLST,3,3,6,ORDE,3FITEM,3,2FITEM,3,-3FITEM,3,6VSYMM,X,P51X, , , ,0,0 ! 建立焊接件的几何模型VGLUE,ALL ! 粘接各体VSEL,S,,,10VATT,1,1,1 ! 附于右焊缝的材料属性VSEL,S,,,1VSEL,A,,,12,17,1VATT,2,1,1 ! 附于两块钢板的材料属性VSEL,S,,,11VATT,3,1,1 ! 附于左焊缝的材料属性ALLSEL,ALLESIZE,0.05 ! 定义单元划分尺寸VSWEEP,ALL ! 划分单元ESEL,S,MAT,,3TOFFST,460 ! 定义温度偏移量!第一步:稳态分析EKILL,ALL ! 杀死左焊缝单元ALLSEL,ALL/SOLUANTYPE,TRANS ! 定义瞬态分析类型TIMINT,OFF ! 关闭时间积分ESEL,S,MAT,,1NSLE ! 选择右焊缝节点D,ALL,TEMP,3000 ! 施加右焊缝初始温度载荷NSEL,INVE ! 选择其它节点D,ALL,TEMP,70 ! 施加初始温度载荷TIME,1 ! 定义求解时间KBC,0 ! 设置为斜坡载荷SOLVE ! 求解!第二步:右侧焊缝相变分析(1到100秒)DDELE,ALL,TEMP ! 删除温度载荷TIMINT,ON ! 打开时间积分TINTP,,,,1 ! 定义瞬态积分参数TIME,100 ! 定义求解时间DELTIME,1,.5,10 ! 定义时间子步AUTOTS,ON ! 打开自动时间开关KBC,1 ! 设置为阶越载荷OUTRES,ERASEOUTRES,ALL,ALL ! 设置结果输出ASEL,S,EXTASEL,U,LOC,Y,0SFA,ALL,,CONV,5E-5,70 ! 施加对流换热载荷ALLSEL,ALLSOLVE ! 求解!第三步:右侧焊缝凝固分析(100到1000秒)TIME,1000 ! 定义求解时间DELTIME,50,10,100 ! 定义时间子步AUTOTS,ON ! 打开自动时间开关SOLVE ! 求解!第四步:激活左侧焊缝单元进行分析(1000到1001秒) EALIVE,ALL ! 激活左侧焊缝单元ALLSEL,allESEL,S,MAT,,3NSLE ! 选择左焊缝节点D,ALL,TEMP,3000 ! 施加左焊缝初始温度载荷TIME,1001 ! 定义求解时间DELTIME,1,1,1 ! 定义时间子步ALLSEL,ALLSOLVE ! 求解!第五步:左侧焊缝相变分析(1001到1100秒) DDELE,ALL,TEMP ! 删除温度载荷TIME,1100 ! 定义求解时间DELTIME,1,.5,10 ! 定义时间子步SOLVE ! 求解!第六步:左侧焊缝凝固分析(1100到2000秒) TIME,2000 ! 定义求解时间DELTIME,100,10,200 ! 定义时间子步SOLVE ! 求解!第七步:后处理/POST1 ! 进入通用后处理器SET,,,,,1, , ! 读取1秒中分析结果PLNSOL, TEMP,, 0 ! 显示1秒钟后焊接件的温度分布SET,,,,,100, , ! 读取100秒中分析结果PLNSOL, TEMP,, 0 ! 显示100秒钟后焊接件的温度分布SET,,,,,1000, , ! 读取1000秒中分析结果PLNSOL, TEMP,, 0 ! 显示1000秒钟后焊接件的温度分布SET,,,,,1001, , ! 读取1001秒中分析结果PLNSOL, TEMP,, 0 ! 显示1001秒钟后焊接件的温度分布SET,,,,,1100, , ! 读取1100秒中分析结果PLNSOL, TEMP,, 0 ! 显示1100秒钟后焊接件的温度分布SET,,,,,2000, , ! 读取2000秒中分析结果PLNSOL, TEMP,, 0 ! 显示2000秒钟后焊接件的温度分布/POST26 ! 进入时间历程后处理器NSOL,2,4727,TEMP,, TEMP_2STORE,MERGENSOL,3,4752,TEMP,,TEMP_3STORE,MERGENSOL,4,4808,TEMP,,TEMP_4STORE,MERGENSOL,5,4833,TEMP,,TEMP_5STORE,MERGENSOL,6,4883,TEMP,,TEMP_6STORE,MERGENSOL,7,4908,TEMP,,TEMP_7STORE,MERGENSOL,8,5088,TEMP,,TEMP_8STORE,MERGENSOL,9,5308,TEMP,, TEMP_9STORE,MERGE ! 定义焊接件某些位置8个节点的时间温度变量/AXLAB,X,TIME/AXLAB,Y,TEMPERATURE ! 更改坐标轴标识/XRANGE,0,2000 ! 设定横坐标轴范围PLVAR,2,3,4,5,6,7,8,9, ! 绘制8节点温度随时间的变化曲线温度结果显示1秒时:100秒时:1000秒时:1001秒时:1100秒时:2000秒时:关键点温度变化曲线关键点位置分布:关键点与曲线对应情况:NSOL,2,4727,TEMP,,TEMP_2NSOL,3,4752,TEMP,,TEMP_3NSOL,4,4808,TEMP,,TEMP_4NSOL,5,4833,TEMP,,TEMP_5NSOL,6,4883,TEMP,,TEMP_6NSOL,7,4908,TEMP,,TEMP_7NSOL,8,5088,TEMP,,TEMP_8NSOL,9,5308,TEMP,,TEMP_9关键点温度变化曲线:从该图中可以明显看到:离焊接热影响区域距离较远的关键点温升较小,距离较近的关键点温升较大。
ANSYS案例——20例ANSYS经典实例】针对【典型例题】3.3.7(1)的模型,即如图3-19所示的框架结构,其顶端受均布力作用,用有限元方法分析该结构的位移。
结构中各个截面的参数都为:113.010PaE=,746.510mI-=,426.810mA-=,相应的有限元分析模型见图3-20。
在ANSYS平台上,完成相应的力学分析。
图3-19框架结构受一均布力作用图3-20单元划分、节点位移及节点上的外载解答对该问题进行有限元分析的过程如下。
1.基于图形界面的交互式操作(tepbytep)(1)进入ANSYS(设定工作目录和工作文件)程序→ANSYS→ANSYSInteractive→Workingdirectory(设置工作目录)→Initialjobname(设置工作文件名):beam3→Run→OK(2)设置计算类型(3)选择单元类型(4)定义材料参数ANSYSMainMenu:Preproceor→MaterialProp→MaterialModel→Struc tural→Linear→Elatic→Iotropic:E某:3e11(弹性模量)→OK→鼠标点击该窗口右上角的“”来关闭该窗口(5)定义实常数以确定平面问题的厚度ANSYSMainMenu:Preproceor→RealContant…→Add/Edit/Delete→Add→Type1Beam3→OK→RealContantSetNo:1(第1号实常数),Cro-ectionalarea:6.8e-4(梁的横截面积)→OK→Cloe(6)生成几何模型生成节点ANSYSMainMenu:Preproceor→Modeling→Creat→Node→InActiveCS→Nodenumber1→某:0,Y:0.96,Z:0→Apply→Nodenumber2→某:1.44,Y:0.96,Z:0→Apply→Nodenumber3→某:0,Y:0,Z:0→Apply→Nodenumber4→某:1.44,Y:0,Z:0→OK生成单元ANSYSMainMenu:Preproceor→Modeling→Create→Element→AutoNum bered→ThruNode→选择节点1,2(生成单元1)→apply→选择节点1,3(生成单元2)→apply→选择节点2,4(生成单元3)→OK(7)模型施加约束和外载左边加某方向的受力ANSYSMainMenu:Solution→DefineLoad→Apply→Structural→Force/Moment→OnNode→选择节点1→apply→Directionofforce:F某→VALUE:3000→OK→上方施加Y方向的均布载荷ANSYSMainMenu:Solution→DefineLoad→Apply→Structural→Preure→OnBeam→选取单元1(节点1和节点2之间)→apply→VALI:4167→VALJ:4167→OK左、右下角节点加约束(8)分析计算(9)结果显示(10)退出系统(11)计算结果的验证与MATLAB支反力计算结果一致。
ansys模拟焊接/FILNAME,weldingstress,0/PREP7!定义焊接参数L=1E-1 !焊件的长度W=1E-1 !焊件的宽度H=6E-3 !焊件的高度U=20 !焊接电压I=160 !焊接电流V=0.01 !焊接速度YITA=0.7 !焊接热效率R=0.007 !电弧有效加热半径Q=U*I*YITA !电弧热功率Qm=3/3.1415/R**2*Q !加热斑点中心最大热流密度!*ET,1,PLANE55ET,2,SOLID70MPTEMP,,,,,,,,MPTEMP,1,20MPTEMP,2,200MPTEMP,3,500 MPTEMP,4,750 MPTEMP,5,1000 MPTEMP,6,1500 MPTEMP,7,1700 MPTEMP,8,2500MPDATA,KXX,1,,50 MPDATA,KXX,1,,47 MPDATA,KXX,1,,40 MPDATA,KXX,1,,27 MPDATA,KXX,1,,30 MPDATA,KXX,1,,35 MPDATA,KXX,1,,40 MPDATA,KXX,1,,55MPDATA,DENS,1,,7820 MPDATA,DENS,1,,7700 MPDATA,DENS,1,,7610 MPDATA,DENS,1,,7550 MPDATA,DENS,1,,7490 MPDATA,DENS,1,,7350MPDATA,DENS,1,,7300 MPDATA,DENS,1,,7090MPDATA,C,1,,460 MPDATA,C,1,,480 MPDATA,C,1,,530MPDATA,C,1,,675 MPDATA,C,1,,670 MPDATA,C,1,,660 MPDATA,C,1,,780 MPDATA,C,1,,820MPDATA,EX,1,,2.05E11 MPDATA,EX,1,,1.87E11 MPDATA,EX,1,,1.5E11 MPDATA,EX,1,,0.7E11 MPDATA,EX,1,,0.2E11 MPDATA,EX,1,,0.19E2 MPDATA,EX,1,,0.18E2 MPDATA,EX,1,,0.12e2MPDATA,PRXY,1,,0.28MPDATA,PRXY,1,,0.29 MPDATA,PRXY,1,,0.31 MPDATA,PRXY,1,,0.35 MPDATA,PRXY,1,,0.4 MPDATA,PRXY,1,,0.45 MPDATA,PRXY,1,,0.48 MPDATA,PRXY,1,,0.5UIMP,1,REFT,,,20 MPDATA,ALPX,1,,1.1e-5 MPDATA,ALPX,1,,1.22e-5 MPDATA,ALPX,1,,1.39e-5 MPDATA,ALPX,1,,1.48e-5 MPDATA,ALPX,1,,1.34e-5 MPDATA,ALPX,1,,1.33e-5 MPDATA,ALPX,1,,1.32e-5 MPDATA,ALPX,1,,1.31e-5TB,BISO,1,6,2, TBTEMP,20 TBDATA,,220e6,0,,,,TBTEMP,250 TBDATA,,175e6,0,,,, TBTEMP,500 TBDATA,,80e6,0,,,, TBTEMP,750 TBDATA,,40E6,0,,,, TBTEMP,1000 TBDATA,,10E6,0,,,, TBTEMP,1500 TBDATA,,1E-5,0,,,,K,1,0,0,0K,2,0,L,0K,3,-W/2*0.15,L,0 K,4,-W/2*0.3,L,0 K,5,-W/2*0.5,L,0 K,6,-W/2,L,0K,7,-W/2,0,0K,8,-W/2*0.5,0,0 K,9,-W/2*0.3,0,0K,10,-W/2*0.15,0,0 K,11,0,0,HA,1,2,3,10A,10,3,4,9A,9,4,5,8A,8,5,6,7ESIZE,0.0012 AMESH,1ESIZE,0.0025 AMESH,2ESIZE,0.005 AMESH,3ESIZE,0.0065 AMESH,4TYPE, 2 EXTOPT,ESIZE,2,0, EXTOPT,ACLEAR,1!*EXTOPT,ATTR,1,0,0 REAL,_Z4ESYS,0!*VOFFST,1,H, , VOFFST,2,H, , VOFFST,3,H, , VOFFST,4,H, , EPLOT NUMMRG,ALL, , , ,LOW/SOL!*ANTYPE,4!*TRNOPT,FULLLUMPM,0*DEL,_FNCNAME*DEL,_FNCMTID*DEL,_FNC_C1*DEL,_FNC_C2*DEL,_FNC_C3*DEL,_FNCCSYS*SET,_FNCNAME,'GAOSI'*DIM,_FNC_C1,,1*DIM,_FNC_C2,,1*DIM,_FNC_C3,,1*SET,_FNC_C1(1),QM*SET,_FNC_C2(1),V*SET,_FNC_C3(1),R*SET,_FNCCSYS,0! /INPUT,HANJIE.func,,,1*DIM,%_FNCNAME%,TABLE,6,19,1,,,,%_FNCCSYS%!! Begin of equation: Qm*exp(-3*({X}^2+({Y}-V*{TIME})^2)/R^2) *SET,%_FNCNAME%(0,0,1), 0.0, -999 *SET,%_FNCNAME%(2,0,1), 0.0*SET,%_FNCNAME%(3,0,1), %_FNC_C1(1)%*SET,%_FNCNAME%(4, 0,1), %_FNC_C2(1)%*SET,%_FNCNAME%(5,0,1), %_FNC_C3(1)%*SET,%_FNCNAME%(6,0,1), 0.0*SET,%_FNCNAME%(0,1,1), 1.0, -1, 0, 0, 0, 0, 0*SET,%_FNCNAME%(0,2,1), 0.0, -2, 0, 1, 0, 0, -1*SET,%_FNCNAME%(0,3,1), 0, -3, 0, 1, -1, 2, -2 *SET,%_FNCNAME%(0,4,1), 0.0, -1, 0, 3, 0, 0, -3*SET,%_FNCNAME%(0,5,1), 0.0, -2, 0, 1, -3, 3, -1 *SET,%_FNCNAME%(0,6,1), 0.0, -1, 0, 2, 0, 0, 2*SET,%_FNCNAME%(0,7,1), 0.0, -3, 0, 1, 2, 17, -1 *SET,%_FNCNAME%(0,8,1), 0.0, -1, 0, 1, 18, 3, 1*SET,%_FNCNAME%(0,9,1), 0.0, -4, 0, 1, 3, 2, -1*SET,%_FNCNAME%(0,10,1), 0.0, -1, 0, 2, 0, 0, -4 *SET,%_FNCNAME%(0,11,1), 0.0, -5, 0, 1, -4, 17, -1 *SET,%_FNCNAME%(0,12,1), 0.0, -1, 0, 1, -3, 1, -5 *SET,%_FNCNAME%(0,13,1), 0.0, -3, 0, 1, -2, 3, -1 *SET,%_FNCNAME%(0,14,1), 0.0, -1, 0, 2, 0, 0, 19 *SET,%_FNCNAME%(0,15,1), 0.0, -2, 0, 1, 19, 17, -1 *SET,%_FNCNAME%(0,16,1), 0.0, -1, 0, 1, -3, 4, -2 *SET,%_FNCNAME%(0,17,1), 0.0, -1, 7, 1, -1, 0, 0 *SET,%_FNCNAME%(0,18,1), 0.0, -2, 0, 1, 17, 3, -1 *SET,%_FNCNAME%(0,19,1), 0.0, 99, 0, 1, -2, 0, 0! End of equation: Qm*exp(-3*({X}^2+({Y}-V*{TIME})^2)/R^2) TUNIF,20, !定义初始温度!定义对流换热边界SFA,15,1,CONV,30,20SFA,20,1,CONV,30,20SFA,9,1,CONV,30,20SFA,14,1,CONV,30,20SFA,19,1,CONV,30,20SFA,24,1,CONV,30,20SFA,23,1,CONV,30,20SFA,7,1,CONV,30,20SFA,12,1,CONV,30,20SFA,17,1,CONV,30,20SFA,22,1,CONV,30,20SFA,1,1,CONV,30,20SFA,2,1,CONV,30,20SFA,3,1,CONV,30,20SFA,4,1,CONV,30,20!施加高斯热源SFA,5,1,HFLUX, %GAOSI% SFA,10,1,HFLUX, %GAOSI%OUTRES,ALL,ALL,TIME,L/V !设置求解时间AUTOTS,-1NSUBST,50,50,50KBC,0TSRES,ERASELSWRITE,1, !写入载荷文件为1 !*TIME,20AUTOTS,1NSUBST,20,20,20KBC,0!*TSRES,ERASELSWRITE,2,TIME,50AUTOTS,1NSUBST,30,30,30KBC,0!*TSRES,ERASELSWRITE,3,TIME,1100AUTOTS,1NSUBST,105,105,105KBC,0!*TSRES,ERASELSWRITE,4,LSSOLVE,1,4,1, !开始求解!进入热应力求解/PREP7ETCHG,TTS !热到结构分析转换/SOL ANTYPE,4TRNOPT,FULLLUMPM,0NLGEOM,1NROPT,FULL, ,OFF!定义热应力计算参考温度TREF,20,DA,6,SYMMDL,13, ,UY,DL,10, ,UY,DL,7, ,UY,DL,4, ,UY,DL,1, ,UZ,*DO,I,1,50LDREAD,TEMP,,,0.2*I, ,'weldingstress','rth',' ' !读入热分析的计算结果OUTRES,ALL,ALL,TIME,0.2*IDELTIM,0.2,0.075,0.2,1SOLVE*ENDDO*DO,I,1,20LDREAD,TEMP,,,10+I*0.5, ,'weldingstress','rth',' ' OUTRES,ALL,ALL,TIME,10+I*0.5DELTIM,0.5,0.5,1,1SOLVE*ENDDO*DO,I,1,30LDREAD,TEMP,,,20+I, ,'weldingstress','rth',' ' OUTRES,ALL,ALL, TIME,20+IDELTIM,1,1,1,1SOLVE*ENDDO*DO,I,1,105LDREAD,TEMP,,,50+10*I, ,'weldingstress','rth',' ' OUTRES,ALL,ALL,TIME,50+10*IDELTIM,10,10,10,1 SOLVE*ENDDO。
ANSYS Mechanical在焊接仿真中的应用详细解析1 前言焊接作为现代制造业必不可少的工艺,在材料加工领域一直占有重要地位。
焊接是一个涉及到电弧物理、传热、冶金和力学等各学科的复杂过程,其涉及到的传热过程、金属的融化和凝固、冷却时的相变、焊接应力和变形等是企业制造部门和设计人员关心的重点问题。
焊接过程中产生的焊接应力和变形,不仅影响焊接结构的制造过程,而且还影响焊接结构的使用性能。
这些缺陷的产生主要是焊接时不合理的热过程引起的。
由于高能量的集中的瞬时热输入,在焊接过程中和焊后将产生相当大的残余应力和变形,影响结构的加工精度和尺寸的稳定性。
因此对于焊接温度场合应力场的定量分析、预测有重要意义。
传统的焊接温度场和应力测试依赖于设计人员的经验或基于统计基础的半经验公式,但此类方法带有明显的局限性,对于新工艺无法做到前瞻性的预测,从而导致实验成本急剧增加,因此针对焊接采用数值模拟的方式体现出了巨大优势。
ANSYS作为世界知名的通用结构分析软件,提供了完整的分析功能,完备的材料本构关系,为焊接仿真提供了技术保障。
文中以ANSYS为平台,阐述了焊接温度场仿真和热变形、应力仿真的基本理论和仿真流程,为企业设计人员提供了一定的参考。
2 焊接数值模拟理论基础焊接问题中的温度场和应力变形等最终可以归结为求解微分方程组,对于该类方程求解的方式通常为两大类:解析法和数值法。
由于只有在做了大量简化假设,并且问题较为简单的情况下,才可能用解析法得到方程解,因此对于焊接问题的模拟通常采用数值方法。
在焊接分析中,常用的数值方法包括:差分法、有限元法、数值积分法、蒙特卡洛法。
差分法:差分法通过把微分方程转换为差分方程来进行求解。
对于规则的几何特性和均匀的材料特性问题,编程简单,收敛性好。
但该方法往往仅局限于规则的差分网格(正方形、矩形、三角形等),同时差分法只考虑节点的作用,而不考虑节点间单元的贡献,常常用来进行焊接热传导、氢扩散等问题的研究。
第一章前处理第1例 关键点和线的创建实例—正弦曲线FINISH/CLEAR, NOSTART /PREP7K,100,0,0,0CIRCLE,100,1,,,90 CSYS,1KFILL,2,1,4,3,1K,7,1+3.1415926/2,0,0 CSYS,0KFILL,7,1,4,8,1 KGEN,2,7,11,1,,1 LSTR,8,13 LSTR,9,14 LSTR,10,15 LSTR,11,16 LANG,5,6,90,,0 LANG,4,5,90,,0 LANG,3,4,90,,0 LANG,2,3,90,,0BSPLIN,1,17,18,19,20,12 LSEL,U,,,14LDELE,ALL LSEL,ALL KWPAVE,12 CSYS,4LSYMM,X,14NUMMRG,KP,,,,LOWLCOMB,ALL,,0FINISH/CLEAR, NOSTART /PREP7 PI=3.14159 J=0*DO,I,0,PI,PI/10.0 J=J+1 X=IY=SIN(I) I=I+1 K,J,X,Y *ENDDOBSPLIN,1,2,3,4,5,6 BSPLIN,6,7,8,9,10,11 csys,4 KWPAVE,11LSYMM,y,1,2,,,,0 KWPAVE,11LSYMM,x,3,4,,,,1以上程序有意没算到2 为了使用几个命令第2例工作平面的应用实例—相交圆柱体[本例提示]通过相交圆柱体的创建,本例主要介绍了工作平面的使用方法。
通过本例,读者可以了解并掌握工作平面与所创建体的位置、方向的关系,学习工作平面的设置、偏移、旋转和激活为当前坐标系的方法。
FINISH/CLEAR,NOSTART/PREP7CYLIND,0.015,0,0,0.08,0,360CYLIND,0.03,0,0,0.08,0,360/VIEW,1,1,1,1/PNUM,VOLU,1WPOFF,0,0.05,0.03WPROT,0,60CYLIND,0.012,0,0,0.055,0,360CYLIND,0.006,0,0,0.055,0,360VSEL,S,,,2,3,1CM,VV1,VOLUVSEL,INVECM,VV2,VOLUVSEL,ALLVSBV,VV1,VV2BLOCK,-0.002,0.002,-0.013,-0.009,0,0.008WPSTYLE,,,,,,1CSYS,4VGEN,3,1,,,,120VSBV,5,1VSBV,4,2VSBV,1,3WPROT,0,0,90VSBW,ALLVDELE,1,4,3VADD,ALLVPLOT/REPLOT第3例复杂形状实体的创建实例—螺栓[本例提示]在使用ANSYS软件进行结构分析时,建立实体模型是最复杂最难以掌握的一个过程。
ansys模拟焊接/FILNAME,weldingstress,0 /PREP7!定义焊接参数L=1E-1!焊件的长度W=1E-1!焊件的宽度H=6E-3!焊件的高度U=20!焊接电压I=160!焊接电流V=0.01!焊接速度YITA=0.7!焊接热效率R=0.007!电弧有效加热半径Q=U*I*YITA!电弧热功率Qm=3/3.1415/R**2*Q!加热斑点中心最大热流密度!*ET,1,PLANE55ET,2,SOLID70MPTEMP,,,,,,,,MPTEMP,1,20MPTEMP,2,200MPTEMP,3,500MPTEMP,4,750MPTEMP,5,1000MPTEMP,6,1500MPTEMP,7,1700MPTEMP,8,2500MPDATA,KXX,1,,50 MPDATA,KXX,1,,47 MPDATA,KXX,1,,40 MPDATA,KXX,1,,27 MPDATA,KXX,1,,30 MPDATA,KXX,1,,35 MPDATA,KXX,1,,40 MPDATA,KXX,1,,55 MPDATA,DENS,1,,7820 MPDATA,DENS,1,,7700 MPDATA,DENS,1,,7610 MPDATA,DENS,1,,7550 MPDATA,DENS,1,,7490 MPDATA,DENS,1,,7350 MPDATA,DENS,1,,7300 MPDATA,DENS,1,,7090 MPDATA,C,1,,460 MPDATA,C,1,,480 MPDATA,C,1,,530 MPDATA,C,1,,675 MPDATA,C,1,,670 MPDATA,C,1,,660 MPDATA,C,1,,780 MPDATA,C,1,,820 MPDATA,EX,1,,2.05E11 MPDATA,EX,1,,1.87E11 MPDATA,EX,1,,1.5E11 MPDATA,EX,1,,0.7E11 MPDATA,EX,1,,0.2E11 MPDATA,EX,1,,0.19E2 MPDATA,EX,1,,0.18E2 MPDATA,EX,1,,0.12e2 MPDATA,PRXY,1,,0.28 MPDATA,PRXY,1,,0.29 MPDATA,PRXY,1,,0.31 MPDATA,PRXY,1,,0.35 MPDATA,PRXY,1,,0.4 MPDATA,PRXY,1,,0.45 MPDATA,PRXY,1,,0.48MPDATA,PRXY,1,,0.5 UIMP,1,REFT,,,20 MPDATA,ALPX,1,,1.1e-5 MPDATA,ALPX,1,,1.22e-5 MPDATA,ALPX,1,,1.39e-5 MPDATA,ALPX,1,,1.48e-5 MPDATA,ALPX,1,,1.34e-5 MPDATA,ALPX,1,,1.33e-5 MPDATA,ALPX,1,,1.32e-5 MPDATA,ALPX,1,,1.31e-5 TB,BISO,1,6,2, TBTEMP,20 TBDATA,,220e6,0,,,, TBTEMP,250 TBDATA,,175e6,0,,,, TBTEMP,500 TBDATA,,80e6,0,,,, TBTEMP,750 TBDATA,,40E6,0,,,, TBTEMP,1000 TBDATA,,10E6,0,,,, TBTEMP,1500 TBDATA,,1E-5,0,,,,K,1,0,0,0K,2,0,L,0K,3,-W/2*0.15,L,0K,4,-W/2*0.3,L,0K,5,-W/2*0.5,L,0K,6,-W/2,L,0K,7,-W/2,0,0K,8,-W/2*0.5,0,0K,9,-W/2*0.3,0,0K,10,-W/2*0.15,0,0K,11,0,0,HA,1,2,3,10A,10,3,4,9A,9,4,5,8A,8,5,6,7ESIZE,0.0012 AMESH,1ESIZE,0.0025 AMESH,2ESIZE,0.005 AMESH,3ESIZE,0.0065 AMESH,4TYPE,2 EXTOPT,ESIZE,2,0, EXTOPT,ACLEAR,1!*EXTOPT,ATTR,1,0,0 REAL,_Z4ESYS,0!*VOFFST,1,H,, VOFFST,2,H,, VOFFST,3,H,, VOFFST,4,H,,EPLOTNUMMRG,ALL,,,,LOW /SOL!*ANTYPE,4!*TRNOPT,FULL LUMPM,0*DEL,_FNCNAME*DEL,_FNCMTID*DEL,_FNC_C1*DEL,_FNC_C2*DEL,_FNC_C3*DEL,_FNCCSYS*SET,_FNCNAME,'GAOSI'*DIM,_FNC_C1,,1*DIM,_FNC_C2,,1*DIM,_FNC_C3,,1*SET,_FNC_C1(1),QM*SET,_FNC_C2(1),V*SET,_FNC_C3(1),R*SET,_FNCCSYS,0!/INPUT,HANJIE.func,,,1*DIM,%_FNCNAME%,TABLE,6,19,1,,,,%_FNCCSYS%!!Begin of equation:Qm*exp(-3*({X}^2+({Y}-V*{TIME})^2)/R^2) *SET,%_FNCNAME%(0,0,1),0.0,-999*SET,%_FNCNAME%(2,0,1),0.0*SET,%_FNCNAME%(3,0,1),%_FNC_C1(1)%*SET,%_FNCNAME%(4,0,1),%_FNC_C2(1)%*SET,%_FNCNAME%(5,0,1),%_FNC_C3(1)%*SET,%_FNCNAME%(6,0,1),0.0*SET,%_FNCNAME%(0,1,1),1.0,-1,0,0,0,0,0*SET,%_FNCNAME%(0,2,1),0.0,-2,0,1,0,0,-1*SET,%_FNCNAME%(0,3,1),0,-3,0,1,-1,2,-2*SET,%_FNCNAME%(0,4,1),0.0,-1,0,3,0,0,-3*SET,%_FNCNAME%(0,5,1),0.0,-2,0,1,-3,3,-1*SET,%_FNCNAME%(0,6,1),0.0,-1,0,2,0,0,2*SET,%_FNCNAME%(0,7,1),0.0,-3,0,1,2,17,-1*SET,%_FNCNAME%(0,8,1),0.0,-1,0,1,18,3,1*SET,%_FNCNAME%(0,9,1),0.0,-4,0,1,3,2,-1*SET,%_FNCNAME%(0,10,1),0.0,-1,0,2,0,0,-4*SET,%_FNCNAME%(0,11,1),0.0,-5,0,1,-4,17,-1*SET,%_FNCNAME%(0,12,1),0.0,-1,0,1,-3,1,-5*SET,%_FNCNAME%(0,13,1),0.0,-3,0,1,-2,3,-1*SET,%_FNCNAME%(0,14,1),0.0,-1,0,2,0,0,19*SET,%_FNCNAME%(0,15,1),0.0,-2,0,1,19,17,-1*SET,%_FNCNAME%(0,16,1),0.0,-1,0,1,-3,4,-2*SET,%_FNCNAME%(0,17,1),0.0,-1,7,1,-1,0,0*SET,%_FNCNAME%(0,18,1),0.0,-2,0,1,17,3,-1*SET,%_FNCNAME%(0,19,1),0.0,99,0,1,-2,0,0!End of equation:Qm*exp(-3*({X}^2+({Y}-V*{TIME})^2)/R^2) TUNIF,20,!定义初始温度!定义对流换热边界SFA,15,1,CONV,30,20SFA,20,1,CONV,30,20SFA,9,1,CONV,30,20SFA,14,1,CONV,30,20SFA,19,1,CONV,30,20SFA,24,1,CONV,30,20SFA,23,1,CONV,30,20SFA,7,1,CONV,30,20SFA,12,1,CONV,30,20SFA,17,1,CONV,30,20SFA,22,1,CONV,30,20SFA,1,1,CONV,30,20SFA,2,1,CONV,30,20SFA,3,1,CONV,30,20SFA,4,1,CONV,30,20!施加高斯热源SFA,5,1,HFLUX,%GAOSI% SFA,10,1,HFLUX,%GAOSI% OUTRES,ALL,ALL,TIME,L/V!设置求解时间AUTOTS,-1NSUBST,50,50,50KBC,0TSRES,ERASELSWRITE,1,!写入载荷文件为1 !*TIME,20AUTOTS,1NSUBST,20,20,20KBC,0!*TSRES,ERASELSWRITE,2,TIME,50AUTOTS,1NSUBST,30,30,30KBC,0!*TSRES,ERASELSWRITE,3,TIME,1100AUTOTS,1NSUBST,105,105,105 KBC,0!*TSRES,ERASELSWRITE,4,LSSOLVE,1,4,1,!开始求解!进入热应力求解/PREP7ETCHG,TTS!热到结构分析转换/SOLANTYPE,4TRNOPT,FULLLUMPM,0NLGEOM,1NROPT,FULL,,OFF!定义热应力计算参考温度TREF,20,DA,6,SYMMDL,13,,UY,DL,10,,UY,DL,7,,UY,DL,4,,UY,DL,1,,UZ,*DO,I,1,50LDREAD,TEMP,,,0.2*I,,'weldingstress','rth', ''!读入热分析的计算结果OUTRES,ALL,ALL,TIME,0.2*IDELTIM,0.2,0.075,0.2,1SOLVE*ENDDO*DO,I,1,20LDREAD,TEMP,,,10+I*0.5,,'weldingstress','rth 9;,''OUTRES,ALL,ALL,TIME,10+I*0.5DELTIM,0.5,0.5,1,1SOLVE*ENDDO*DO,I,1,30LDREAD,TEMP,,,20+I,,'weldingstress','rth',& #39;'OUTRES,ALL,ALL,TIME,20+IDELTIM,1,1,1,1SOLVE*ENDDO*DO,I,1,105LDREAD,TEMP,,,50+10*I,,'weldingstress','rth 9;,''OUTRES,ALL,ALL,TIME,50+10*IDELTIM,10,10,10,1SOLVE*ENDDO。
ANSYS的生死单元模拟焊接过程长安CAE2016年10月11日1388ANSYS的生死单元模拟焊接过程1 概述焊接模拟计算在CAE仿真是比较大的一块内容,也是比较复杂的一个过程,几个比较关键的问题是热源函数的描述、单元的融覆、热源的移动等等,通过单纯的GUI操作,无论使ANSYS还是Abaqus都不大可能完成这个过程,通常需要借助软件的内置语言。
本次主要介绍单元生死的应用,单元生死主要用于单元缺失的场合,比如凝固溶解过程,断裂过程,焊接过程等等,这些过程都是非线性或者时间历程过程,计算需要很多子步和迭代,为了在此过程中避免一遍一遍修改单元,便引入生死单元的概念,通俗的讲就是通过一些方法让单元失效,具体的改变是单元的弹性模量的改变,当单元死时,修改其弹性模量为非常小的值,让其在求解过程中不起作用。
详细地说,激活单元死这个状态时,ANSYS程序将单元刚度矩阵乘以很小的因子,程序默认值为1E-6,死单元的单元载荷为0,从而不对载荷向量生效,同样的,死单元的质量、阻尼、比热等等参数也设置为0,单元的应力应变也因此为0。
2 前处理前处理包括单元定义、材料定义和建模,单元定义是需要注意单元属性,此次定义13号二维耦合单元,具有温度和位移自由度。
材料属性包括结构参数和热参数,具体包含弹性模量,泊松比,屈服强度,塑性属性,材料密度,热膨胀系数,热传导系数,比热容。
焊接时温度较高,定义材料通常需要定义多个温度下的值。
例如定义各材料在各温度点下的屈服应力和屈服后的弹性模量:tb,bkin,1,5tbtemp,20,1tbdata,1,1200e6,0.193e11tbtemp,500,2tbdata,1, 933e6,0.150e11tbtemp,1000,3tbdata,1, 435e6,0.070e11tbtemp,1500,4tbdata,1, 70e6,0.010e11tbtemp,2000,5tbdata,1, 7e6,0.001e11建立的二维模型如图1所示,中间三角区域为焊接区域。
焊接模拟ansys实例!下面的命令流进行的是一个简单的二维焊接分析, 利用ANSYS单元生死和热-结构耦合分析功能进!行焊接过程仿真, 计算焊接过程中的温度分布和应力分布以及冷却后的焊缝残余应力。
finish/clear/filnam,1-2D element birth and death/title,Weld Analysis by "Element Birth and Death"/prep7/unit,si !采用国际单位制!******************************************************et,1,13,4 !13号二维耦合单元, 同时具有温度和位移自由度et,2,13,4!1号材料是钢!2号材料是铝!3号材料是铜!铝是本次分析中的焊料, 它将钢结构部分和铜结构部分焊接起来!下面是在几个温度点下, 各材料的弹性模量mptemp,1,20,500,1000,1500,2000mpdata,ex,1,1,1.93e11,1.50e11,0.70e11,0.10e11,0.01e11mpdata,ex,2,1,1.02e11,0.50e11,0.08e11,0.001e11,0.0001e11mpdata,ex,3,1,1.17e11,0.90e11,0.30e11,0.05e11,0.005e11!假设各材料都是双线性随动硬化弹塑性本构关系!下面是各材料在各温度点下的屈服应力和屈服后的弹性模量tb,bkin,1,5tbtemp,20,1tbdata,1,1200e6,0.193e11tbtemp,500,2tbdata,1, 933e6,0.150e11tbtemp,1000,3tbdata,1, 435e6,0.070e11tbtemp,1500,4tbdata,1, 70e6,0.010e11tbtemp,2000,5tbdata,1, 7e6,0.001e11!tb,bkin,2,5tbtemp,20,1tbdata,1,800e6,0.102e11tbtemp,500,2tbdata,1,400e6,0.050e11tbtemp,1000,3tbdata,1, 70e6,0.008e11tbdata,1, 1e6,0.0001e11tbtemp,2000,5tbdata,1,0.1e6,0.00001e11!tb,bkin,3,5tbtemp,20,1tbdata,1,900e6,0.117e11tbtemp,500,2tbdata,1,700e6,0.090e11tbtemp,1000,3tbdata,1,230e6,0.030e11tbtemp,1500,4tbdata,1, 40e6,0.005e11tbtemp,2000,5tbdata,1, 4e6,0.0005e11!!材料密度(假设为常值)mp,dens,1,8030mp,dens,2,4850mp,dens,3,8900! 热膨胀系数(假设为常值)mp,alpx,1,1.78e-5mp,alpx,2,9.36e-6mp,alpx,3,1.66e-5!泊松比(假设为常值)mp,nuxy,1,0.29mp,nuxy,2,0.30mp,nuxy,3,0.30!热传导系数(假设为常值)mp,kxx,1,16.3mp,kxx,2,7.44mp,kxx,3,393!比热(假设为常值)mp,c,1,502mp,c,2,544mp,c,3,385!热膨胀系数(假设为常值)!由于该13号单元还有磁自由度, 此处假设一磁特性, 但并不影响我们所关心的结果mp,murx,1,1mp,murx,2,1mp,murx,3,1!假设焊料(铝)焊上去后的初始温度是1500℃mp,reft,2,1500mp,reft,3,20!******************************************************!下面建立几何模型csys,0k,1,0,0,0k,2,0.5,0,0k,3,1,0,0 !长1米k,4,0,0.3,0 !厚度0.3米(二维中叫做宽度)k,5,0.35,0.3,0k,6,0.65,0.3,0k,7,1,0.3,0a,1,2,5,4a,2,6,5a,2,3,7,6!!!!!!!!!!!!!!!!!!!!!!!!!!划分网格!!!!!!!!!!!!!!!!!esize,0.025type,2mat,2amesh,2!esize,0.05 !网格划分出现问题type,1mat,1amesh,1!mat,3amesh,3eplot!/soluantype,4 ! 瞬态分析trnopt,full!!!!!!!!!!!!!!!!!!!!!!!!在模型的左边界加位移约束!!!!!!!!!!!!!!!!!!!!!!!!!!!nsel,all*get,minx,node,,mnloc,xnsel,s,loc,x,minxd,all,ux,0*get,miny,node,,mnloc,ynsel,r,loc,y,minyd,all,uy,0!*****假设模型的左右边界处温度始终保持在20摄氏度左右*****!其他边界条件如对流和辐射等均可施加,此处因为只是示意而已,故只施加恒温边界条件nsel,all*get,minx,node,,mnloc,xnsel,s,loc,x,minxd,all,temp,20nsel,all*get,maxx,node,,mxloc,xnsel,s,loc,x,maxxd,all,temp,20!!!!!!!!!!!!!!!由于第2个面是焊接所在区域,因此首先将该区域的单元“杀死”!!!!!!!!!!!!!!!!nna=2esel,all*get,emax,elem,,num,maxasel,s,area,,nnaesla*get,nse,elem,,count*dim,ne,,nse*dim,nex,,nse*dim,ney,,nse*dim,neorder,,nsemine=0!**********************************************!下面的do循环用于将焊料区的单元按其形心y坐标排序!以便后面模拟焊料由下向上逐步“生长”过程*do,i1,1,nseesel,u,elem,,mine*get,nse1,elem,,countii=0*do,i,1,emax*if,esel(i),eq,1,thenii=ii+1ne(ii)=i*endif*enddo*do,i,1,nse1*get,ney(i),elem,ne(i),cent,y*get,nex(i),elem,ne(i),cent,x*enddominy=1e20minx=1e20*do,i,1,nse1*if,ney(i),lt,miny,thenminy=ney(i)minx=nex(i)mine=ne(i)*else*if,ney(i),eq,miny,then*if,nex(i),lt,minx,thenminy=ney(i)minx=nex(i)mine=ne(i)*endif*endif*endif*enddoneorder(i1)=mine*enddo!************************************************************** max_tem=1500 !按照前面假设,焊料的初始温度为1500℃dt1=1e-3 !用于建立初始条件的一个很小的时间段dt=5 !焊接一个单元所需的时间t=0 !起始时间esel,alleplot/auto,1/replot*do,i,1,nseekill,neorder(i)esel,s,liveeplot*enddoallsel,alloutres,all,allic,all,temp,20kbc,1timint,0,structtimint,1,thermtimint,0,magtintp,0.005,,,1,0.5,0.2!nsub1=2nsub2=40!**************************************************do,i,1,nseealive,neorder(i)esel,s,liveeplotesel,all!******下面的求解用于建立温度的初始条件******t=t+dt1time,tnsubst,1*do,j,1,4d,nelem(neorder(i),j),temp,max_tem*enddosolve!****下面的求解用于保证初始的升温速度为零****t=t+dt1time,tsolve!*********下面的步骤用于求解温度分布***********do,j,1,4ddele,nelem(neorder(i),j),temp*enddot=t+dt-2*dt1time,tnsubst,nsub1solve*enddot=t+50000 !*********下面的步骤用于冷却过程求解***** time,tnsubst,nsub2solvesavefinish!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!后处理过程!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!/post1!**************下面的一系列命令用于生成应力的动画文件******* /seg,dele/cont,1,15,0,1200e6/16,1200e6/dscale,1,1.0avprin,0,0avres,1/seg,multi,stress1,0.1esel,all*do,i,1,nseesel,u,elem,,neorder(i)*enddo*do,i,1,nseesel,a,elem,,neorder(i)set,(i-1)*3+1,1plnsol,s,eqv*do,j,1,nsub1set,(i-1)*3+3,jplnsol,s,eqv*enddo*enddo*do,i,1,nsub2set,(nse-1)*3+4,iplnsol,s,eqv*enddo/seg,off,stress1,0.1/anfile,save,stress1,avi!**********下面的一系列命令用于生成温度的动画文件************ /seg,dele/cont,1,15,0,1500/16,1500/dscale,1,1.0avprin,0,0avres,1/seg,multi,temp1,0.1esel,all*do,i,1,nseesel,u,elem,,neorder(i)*enddo*do,i,1,nseesel,a,elem,,neorder(i)set,(i-1)*3+1,1plnsol,temp*do,j,1,nsub1set,(i-1)*3+3,jplnsol,temp*enddo*enddo*do,i,1,nsub2set,(nse-1)*3+4,iplnsol,temp*enddo/seg,off,temp1,0.1/anfile,save,temp1,avifinish。