实例6 框架结构弹性时程分析
1)问题描述:
本例仍采用实例4(实例5)的框架结构,为了方便对比,改采用弹性截面,主要进行弹性时程分析,材料为弹性,时程分析即动力分析。结构荷载情况与实例4相同(侧向力荷载不需要施加,实际侧向力为地面加速度)。计算结构在地震作用下的响应(主要提取位移结果)。(重力荷载代表值组合为1.0×DEAD+0.5×LIVE)。
注意:上述实例讲到了质量矩阵(质量源)的定义,刚度矩阵通过结构几何与材料属性得到,那么接下来只需要定义了结构阻尼,就可以进行结构动力分析,即时程分析。
2)ETABS模型建模
(1) 结构模型与实例5相同,相关建模细节请看实例5。
图 ETABS建立框架的几何模型
(2) 为了对比OPENSEES弹性时程的分析结果,在ETABS模型同样进行弹性时程分析。
为了进行弹性时程分析,需要输入地震波,(本例只进行单向地震分析)。
(3) 地震波数据导入:(实例的地震波放在光盘“/EXAM06/ETABS/”目录)本算例采用的为单向地震波,地震波文件为:GM1X.txt,通过EXCEL图表画出整个地震波时程曲线如下图所示。(该时程的时间间隔为0.02 sec)
图时程曲线GM1X
从时程曲线可知,曲线最大值为3621(该值为无单位数)。
【Define】→【Time History Functions】,如下窗口所示。选取【Function from File】点击【Add New Functions】,弹出以下窗口。
图 GM1X时程曲线定义窗口
(4) 弹性时程分析工况定义,点击【Define】→【Time History Cases】,如下窗口所示。
按窗口的内容填写数据,点击【Modal Damping】右边的按钮【Modify/Show】,可以看到阻尼比的填写框,由于是混凝土结构,阻尼比为0.05。
图时程分析工况内容
窗口中输入的内容简介如下:
荷载步数为2000步,每步代表0.02s,即总分析时间长度为40s,主轴1方向的地
震波时程数据为 GM1X,其放大倍数为0.138(单位为mm/s2),也就是说整个地
震波的最大地面加速度为3621×0.138=500 mm/s2,即50gal,属于小震量级。(5) 完成上述步骤后建立完ETABS模型。
注意:实例的ETABS模型存放在光盘“/EXAM06/ETABS/”目录。
3)ETABS弹性时程分析结果
(1) 完成ETABS模型后,运行分析。分析完成后,点击菜单中的【Display】→【Show Time History Trace】,选取顶层的4号节点,显示其水平位移的时程,如下图所示。最大位移为7.217mm。
图选取结点位移窗口
图顶部结点(4号结点)位移时程曲线
4)OPENSEES建模
(1) 打开ETABS模型,导出S2K文件。打开ETO程序,导入S2K文件,得到转化的OPENSEES模型,如下图所示。再打开转化TCL按扭,将模型转化成OPENSEES 代码,如下图所示。将代码另存为“Exam06.tcl”。
ETO导入ETABS模型
(2) 在ETO程序中,点击按钮,可以设置结构分析工况。本实例选择OPENSEES 的分析类型为【Time History Analysis】,即单工况的时程分析。
分析设置窗口
(3) 点击按钮,可设置OPENSEES的输出命令(Recorder),勾选如下图所示。
ETO结果输出定义窗口
(4) 点击按钮生成OPENSEES命令流。
(5) 以下将对OPENSEES命令流进行解释并修改,最后提交运算。
5)OPENSEES命令流解读
(1) 从ETO程序中生成的OPENSEES的命令流如下所示。
wipe
puts "System"
model basic -ndm 3 -ndf 6
puts "restraint"
node 1 4.500E+003 5.000E+003 1.050E+004
………………
node 28 9.000E+003 5.000E+003 0.000E+000
puts "rigidDiaphragm"
puts "mass"
mass 1 8.604E+000 8.604E+000 8.604E+000 0.000E+000 0.000E+000 0.000E+000 ………………
mass 22 4.302E+000 4.302E+000 4.302E+000 0.000E+000 0.000E+000 0.000E+000
puts "node"
fix 23 1 1 1 1 1 1;
………………
fix 28 1 1 1 1 1 1;
puts "material"
uniaxialMaterial Elastic 1 1.999E+005
uniaxialMaterial Elastic 2 2.680E+004
uniaxialMaterial Elastic 3 1.999E+005
puts "transformation"
geomTransf Linear 1 1.000 0.000 0.000
………………
geomTransf Linear 47 0.000 0.000 1.000
puts "element"
element elasticBeamColumn 1 1 2 1.600E+005 2.680E+004 1.117E+004 3.605E+009 2.133E+009 2.133E+009 1
………………
element elasticBeamColumn 47 19 20 1.800E+005 2.680E+004 1.117E+004 3.708E+009 5.400E+009 1.350E+009 47
puts "recorder"
recorder Node -file node0.out -time -nodeRange 1 28 -dof 1 2 3 disp
recorder Element -file ele0.out -time -eleRange 1 47 localForce
set xDamp 0.05;
set nEigenI 1;
set nEigenJ 2;
set lambdaN [eigen [expr $nEigenJ]];
set lambdaI [lindex $lambdaN [expr $nEigenI-1]];
set lambdaJ [lindex $lambdaN [expr $nEigenJ-1]];
set omegaI [expr pow($lambdaI,0.5)];
set omegaJ [expr pow($lambdaJ,0.5)];
set alphaM [expr $xDamp*(2*$omegaI*$omegaJ)/($omegaI+$omegaJ)];
set betaKcurr [expr 2.*$xDamp/($omegaI+$omegaJ)];
rayleigh $alphaM $betaKcurr 0 0
set IDloadTag 1001;
set iGMfile "GMX.txt";
set iGMdirection "1";
set iGMfact "0.1";
set dt 0.02;
foreach GMdirection $iGMdirection GMfile $iGMfile GMfact $iGMfact {
incr IDloadTag;
set GMfatt [expr 1*$GMfact];
set AccelSeries "Series -dt $dt -filePath $iGMfile -factor $GMfatt";
pattern UniformExcitation $IDloadTag $GMdirection -accel $AccelSeries;
}
constraints Transformation;
numberer Plain;
system UmfPack;
test EnergyIncr 1.0e-4 200;
algorithm Newton
integrator Newmark 0.5 0.25 analysis Transient analyze 1000 0.02
(2) 实例6与实例5的结构模型基本上一样,采用了弹性模型,大部分命令流一样,其
它方面的命令流看上述的实例。
(3) 弹性时程分析需要定义结构阻尼比,在ETABS 中阻尼比定义比较简单,只需要输入
一个参数,在OPENSEES 中,结构采用瑞利(Rayleigh )阻尼,即阻尼矩阵的大小
与结构的质量矩阵,刚度矩阵都相关,瑞利阻尼的计算公式如下,阻尼与刚度质量的关系如下图所示。
瑞利阻尼与频率、质量、刚度的关系
01[][][]c a m a k =+
0121
m n
m n a a ωωξωω????=????+??
?? 式中,ξ为阻尼比,0a 为质量相关系数,1a 为刚度相关系数,[]c 为阻尼矩阵,[]m 为质量矩阵,[]k 为刚度矩阵,,m n ωω为结构两个主振型的圆频率,
由于OPENSEES 能够直接求解振型的特征值,那么特征值与圆频率的关系如下式所示。
ω=
命令流的解读如下:
set xDamp 0.05;————设置阻尼比为0.05 set nEigenI 1; ————主振型1为第1振型 set nEigenJ 2; ————主振型2为第2振型
set lambdaN [eigen [expr $nEigenJ]]; ————求解两阶振型即可
set lambdaI [lindex $lambdaN [expr $nEigenI-1]]; ————提取第1阶特征值 set lambdaJ [lindex $lambdaN [expr $nEigenJ-1]]; ————提取第2阶特征值
set omegaI [expr pow($lambdaI,0.5)]; ————从特征值求圆频率 set omegaJ [expr pow($lambdaJ,0.5)]; ————从特征值求圆频率
set alphaM [expr $xDamp*(2*$omegaI*$omegaJ)/($omegaI+$omegaJ)]; ————alphaM 为0a ,即质量相关系数
set betaKcurr [expr 2.*$xDamp/($omegaI+$omegaJ)]; ————betaKcurr 为1a ,即刚度相关系数 rayleigh $alphaM $betaKcurr 0 0
————定义瑞利阻尼,只需要填写0a ,1a ,其它值为0
(4) 单向地震波数据设置的命令流如下:
set IDloadTag 1001; ————地震波工况号为1001
set iGMfile "GM1X.txt"; ————地震波数据文件名为GM1X.txt ,需要与TCL 文件放在同一个目录下。
set iGMdirection "1"; ————地震波方向为x 方向,即系统1方向 set iGMfact "0.138"; ————地震波峰值放大系数为0.138mm/s2 set dt 0.02; ————地震波时间间隔为0.02s
foreach GMdirection $iGMdirection GMfile $iGMfile GMfact $iGMfact { incr IDloadTag;
set GMfatt [expr 1*$GMfact];
set AccelSeries "Series -dt $dt -filePath $iGMfile -factor $GMfatt"; pattern UniformExcitation $IDloadTag $GMdirection -accel $AccelSeries; } ————以上为多维地震波的输入标准样式,可用于三向地震波,也可以用于单向,数组如果只有一个值,如算例所示,只作为单向波。 pattern UniformExcitation $IDloadTag $GMdirection -accel $AccelSeries; 为荷载工况的命令流,用于输入地震波数据,$IDloadTag 为地震波工况号,$GMdirection 为地震波的方向,-accel $AccelSeries ,为地震波的其它参数,包括文件名,时间间隔,峰值放大倍数等。
(5) 弹性时程分析命令流如下:
constraints Transformation; ————用于动力时程分析 numberer Plain; ————普通的排位方法
system UmfPack; ————采用Umfpack 自由度排列 test EnergyIncr 1.0e-4 200; ————采用能量收敛准则 algorithm Newton ————采用牛顿迭代法
integrator Newmark 0.5 0.25 ————采用Newmark 法对时间进行离散 analysis Transient ————采用时程分析
analyze 1000 0.02————分析步数为1000,时间间隔为0.02s ,即分析20s 。
(6) 为了查看8号结点的水平位移时间数据,记录的命令流如下:
recorder Node -file node8.out -time -node 8 -dof 1 disp
(7) 综上所述,主要修改命令流中的地震波输入的命令流、阻尼比的命令流、结果记录命令流,修改后,可以提交进行分析,修改后的文件可查看
“Exam06\OpenSEES\Exam06.tcl”。
6)OPENSEES分析及分析结果
(1) 打开OpenSEES程序,输入命令:
source Exam06.tcl
运行完成后,程序自动保存结果。
(2) 打开node8.out,提取8号结点的水平位移与ETABS提取对应位置结点的水平位移进行比较,得到以下对比结果,表明结构分析与ETABS分析结果一致。
(3) 打开OPENSEES前后处理程序ETO,点击按钮,显示结构变形。弹出窗口如下图所示
点击【Load Node Deform Data】,选取Exam06.tcl文件,窗口显示结构变形。
【Scaling Factor】需要调整合适,如100000,可以显示合理的变形形状如下图所示。
【load step】为254,即254步,时间为5.08s。
T=5.08s时,结构的整体变形
7)知识点回顾:
(1) ETABS中地震波数据的输入
(2) ETABS中弹性时程分析的设置
(3) OPENSEES中瑞利阻尼的计算及设置
(4) OPENSEES中地震波的输入
(5) OPENSEES中弹性时程分析工况的输入
1、利用零长单元模拟阻尼, uniaxialMaterial Elastic 1 6.8098e6; uniaxialMaterial Viscous 2 3.24e5 1; uniaxialMaterial Parallel 3 3 5; element zeroLength 1 $iNode $jNode -mat 3 -dir 1; 通常有两种方式: (1)truss element and viscous material.(桁架单元和阻尼材料) (2)force-based beam-column element and Maxwell material(基于力的梁柱单元和Maxwell 材料)。 -、如何运行OpenSEES 有三种方法可以执行OpenSees/Tcl命令: 1、interactive交互式 直接将命令输入Prompt。 2、执行文件输入 这种方法是最常用的一种,以source inputfile.tcl方式执行已写好的外部命令文件。 3、Batch模式 即以Opensees inputFile.tcl方式在MS-DOS/Unix promt中运行。 二、定义单位和常数 在编写一个较大的Opensees命令时。最好先定义好单位及常数。在Opensees中,编译器不能自行转换单位。所以一开始就要先定义好。 单位定义包括两部分:首先定义基本单位;再定义合成单位。其中基本单位要相互独立。同时,在定义单位时,既可以按国际公制单位,也可以按私制单位。因些在单位定义文件中可能是混合的。我个人建议,还是采用国际公制单位较好。像国外常用英制单位。很不习惯。对于一些常数,如 和g等常数要事先定义好。 在定义这些单位时所用的命令是“set”。
OPENSEES opensees中的单元问题 梁柱单元 1. Nonlinear BeamColumn 基于有限单元柔度法理论。允许刚度沿杆长变化,通过确定单元控制截面各自的截面抗力和截面刚度矩阵,按照Gauss-Lobatto积分方法沿杆长积分计算出整个单元的抗力与切线刚度矩阵。NonlinearBeamColumn单元对于截面软化行为,构件反应由单元积分点数控制,为保证不同积分点数下构件反应的一致性,可以通过修正材料的应力-应变关系来实现,但同时会造成截面层次反应的不一致,因此需要在截面层次进行二次修正。一根构件不需要单元划分,使用1个单元即可,建议单元内使用4个截面积分点,截面上使用6*6的纤维积分点。 [5] 2. Displacement – Based BeamColumn 基于有限单元刚度法理论。允许刚度沿杆长变化,按照Gauss -Legendre积分方法沿杆长积分计算出整个单元的抗力与切线刚度矩阵。 Displacement - BasedBeam- Column单元对于截面软化行为,构件反应由遭受软化行为的单元长度控制,为保证计算结果的精确性,一般需要将构件离散为更多的单元,而截面层次的反应与构件的单元离散数无关,可以较为准确地反应截面的软化行为。 建议一根构件划分为5个单元,单元内使用4个截面积分点,截面上使用6*6的纤维积分点。[5] 3. Beam With Hinges 基于有限单元柔度法理论。假定单元的非弹性变形集中在构件的两端,在杆件端部设置2个积分控制截面,并设定恰当的塑性铰长度,按照Gauss - Radau积分方法沿塑性铰长度积分来模拟构件和整体结构的非线性反应特点,而杆件中部的区段仍保持弹性。
wipe # Openseesdandun # #Units:kN, m, sec # ----------------- # Start of model generation # ----------------- # CreateModeBulider (with two-dimensions and 3 DOF/node) model basic -ndm 2 -ndf 3 # ----------------- # tag X Y node 1 0.0 0.0 node 2 0.0 0.0 node 3 0.0 2.0 node 4 0.0 4.0 node 5 0.0 6.0 node 6 0.0 8.0 node 7 0.0 10.0 node 8 0.0 12.0 node 9 0.0 14.0 node 10 0.0 16.0 node 11 0.0 18.0 node 12 0.0 20.0 # ----------------- # Fix supports at base of columns # tag DX DY RZ fix 1 1 1 1 # ---------------- # Concrete tag fc ec0 fcuecu # Core concrete (confined) uniaxialMaterial Concrete01 1 -25600.0 -0.00219 -17780.0 -0.01 #Cover concrete (unconfined) uniaxialMaterial Concrete01 2 -23400.0 -0.002 -0.0 -0.006 # STEEL # Reinforcing steel setfy 400000.0; #Yield stress set E 200000000.0;# Young's modulus # tag fy E0 b uniaxialMaterial Steel02 3 $fy $E 0.01 18.5 0.925 0.15 uniaxialMaterial Elastic 11 29043600 uniaxialMaterial Elastic 12 12326600 uniaxialMaterial Elastic 13 587247596 #Define cross-section for nonlinear columns # ---------------------
1、定义梁柱单元局部坐标轴的命令流为: geomTransf Linear $transfTag $vecxzX $vecxzY $vecxzZ 其中,$transfTag 代表局部坐标轴矢量的编号,$vecxzX $vecxzY $vecxzZ 表示局部坐标轴的方向矢量值。 2、OPENSEES 的刚性隔板假定命令流格式为: rigidDiaphragm $perpDirn $masterNodeTag $slaveNodeTag1 $slaveNodeTag2 ... 其中,$perpDirn 表示刚性隔板的方法,如实例中楼板的刚性隔板的平移方向为U1(X 方向)与U2(Y 方向),即1-2 平面,该值应为3。$masterNodeTag 为主结点,$slaveNodeTag1 为从结点。主结点一般为刚性隔板刚心。 实例中:rigidDiaphragm 3 35 2,表示刚性隔板平动方向为1-2 平面,刚心主节点为35 点,2号结点为从结点。 3、弹性梁柱单元的命令流: element elasticBeamColumn $eleTag $iNode $jNode $A $E $G $J $Iy $Iz $transfTag 需要提供截面的截面积A、截面Y 轴惯性矩Iy,截面Z 轴惯性矩Iz,截面扭转矩,截面材料的弹性模量E 及剪切模量G。其中:$transfTag 与$eleTag 是一致的,表示一个单元有自已特定的坐标轴向量,为了编程的方便。 陈:例题三 4、非线性材料模型的定义 (1)uniaxialMaterial Steel01 1 335 200000 0.00001 表示,钢筋的屈服强度为335MPa,弹性模量为200000MPa,硬化系数为0.00001,即屈服平台基本上为水平段。 将混凝土材料本构C40 改为非线性混凝土本构【Concrete01】,命令流如下: (2)uniaxialMaterial Concrete01 2 -26.8 -0.002 -10 -0.0033材料参数意见参 考图所示。 注意:混凝土本构Concrete01 是最简单的混凝土本构,注意数值是负数, 即表示受压段。该本构没有受拉段,即受拉强度为0,表示结构一分析即进 入弹塑性。 5、采用纤维单元,需要定义纤维截面,纤维截面的定义如下面代码所示: section Fiber 1 { fiber -1.125E+002 -2.700E+002 4.500E+003 2 ……… fiber 1.150E+002 -2.650E+002 4.900E+002 1 } 以上命令流表示,纤维截面编号为1,{}内部为子命令流,表示每一个纤维的信息,每一个纤维的定义格式如下: fiber $Y $Z $Area $Mat 命令中,$Y 表示每个纤维的截面Y 坐标(截面中心为原点0);$Z 表示每个纤维的截面Z 坐标(截面中心为原点0);$Area 表示每个纤维的贡献面积;$ Mat 表示每个纤维使用的非线性材料本构的编号。 注意:纤维的坐标与材料切线模量可以组装成截面的刚度,而纤维的坐标与材料的应力可以组装成截面的内力(抗力),那么每个纤维的应变可以通过截面的变形与坐标求出。采用纤维截面的单元,即为基于平截面假定。截面变形求解应变是基于平截面假定的。 6、采用的单元为非线性梁柱单元,即基于柔度法的纤维单元(Nonlinear BeamColumn Element or Force Beam Column Element),需要输入命令流如下: element nonlinearBeamColumn $eleTag $iNode $jNode $numIntgrPts $secTag $transfTag 其中,$eleTag 为单元编号;$iNode 为开始结点;$jNode 为结束结点;$numIntgrPts 为积分点数量;$secTag 为纤维截面编号,$transfTag 为局部坐标轴编号。积分点数量,也就纤维单元的计算截面数量,纤维单元的刚度与抗力是由截面刚度与抗力沿杆件长度积分所得,显然,不能将全部截面积分,只能采用
o p e n s e e s总结 -标准化文件发布号:(9556-EUATWK-MWUB-WUNN-INNUL-DDQTY-KII
1、定义梁柱单元局部坐标轴的命令流为: geomTransf Linear $transfTag $vecxzX $vecxzY $vecxzZ 其中,$transfTag 代表局部坐标轴矢量的编号,$vecxzX $vecxzY $vecxzZ 表示局部坐标轴的方向矢量值。 2、OPENSEES 的刚性隔板假定命令流格式为: rigidDiaphragm $perpDirn $masterNodeTag $slaveNodeTag1 $slaveNodeTag2 ... 其中,$perpDirn 表示刚性隔板的方法,如实例中楼板的刚性隔板的平移方向为U1(X 方向)与U2(Y 方向),即1-2 平面,该值应为3。$masterNodeTag 为主结点,$slaveNodeTag1 为从结点。主结点一般为刚性隔板刚心。 实例中:rigidDiaphragm 3 35 2,表示刚性隔板平动方向为1-2 平面,刚心主节点为35 点,2号结点为从结点。 3、弹性梁柱单元的命令流: element elasticBeamColumn $eleTag $iNode $jNode $A $E $G $J $Iy $Iz $transfTag 需要提供截面的截面积A、截面Y 轴惯性矩Iy,截面Z 轴惯性矩Iz,截面扭转矩,截面材料的弹性模量E 及剪切模量G。其中:$transfTag 与$eleTag 是一致的,表示一个单元有自已特定的坐标轴向量,为了编程的方便。 陈:例题三 4、非线性材料模型的定义 (1)uniaxialMaterial Steel01 1 335 200000 0.00001 表示,钢筋的屈服强度为335MPa,弹性模量为200000MPa,硬化系数为0.00001,即屈服平台基本上为水平段。 将混凝土材料本构C40 改为非线性混凝土本构【Concrete01】, 命令流如下: (2)uniaxialMaterial Concrete01 2 -26.8 -0.002 -10 -0.0033材料参 数意见参考图所示。 注意:混凝土本构Concrete01 是最简单的混凝土本构,注意数值是负数,即表示受压段。该本构没有受拉段,即受拉强度为0,表示结构一分析即进入弹塑性。 5、采用纤维单元,需要定义纤维截面,纤维截面的定义如下面代码所示:
OpenSEES解题一般规律、技巧总结 单位 OpenSEES中是可以用公制单位(N,m)的(而并不是像某些文章中说的“OpenSees默认为英制单位”)。实际上我认为OpenSEES中并没有什么默认单位,只要编程者自己保持单位一致就行;这点类似于SAP2000的风格。 建模顺序 做事要讲究顺序,OpenSEES建模亦如是:必须先定义材料才能离散截面(因为离散截面时要对所划分的截面指定材料属性)。 与之类似的,必须先定义(离散)截面,才能定义非线性梁柱单元(因为定义非线性梁柱单元时要指定单元截面)。 关于BandSPD求解方式 官网关于BandSPD方程形式的评价: "This is a good choice for most small size models. " 并且后面紧跟了一句: "The equations have to be numbered so the widely used RCM (Reverse Cuthill-McKee) numberer is used. " 可见numberer 类型不是随便选,而是要根据方程类型来决定的!
(不过直到作业做完,我对numberer, system, test, algorithm, analysis(还包括geomTransf, constraints)等求解控制命令还是一知半解!我觉得要想弄明白这些命令——得先回头好好翻翻有限元和数值分析的书了!) OpenSEES中默认的计算精度比较高! “0.1000000000000001≠0.1”:(自行总结,未找到官方说明)这是一个真实的故事:我曾在程序中自以为是的将一连串相邻均只有0.1左右的数的差强行赋值为0.1,而没有采用循环命令将两数作差并将结果赋给新变量——其中即有这样的强行截断!我以为小数点后都n位了,即使我带着它最后也会被系统截断,还不如我直接预处理来得清爽!没想到这样做直接导致计算不收敛!真是失之毫厘谬以千里!可见在OpenSEES中默认的计算精度比较高! 后来我还在老师给的一份范例程序(Silvia Mazzoni & Frank McKenna, 2006)中发现了这么一段: …… set Ubig 1.e10; # a really large number set Usmall [expr 1/$Ubig]; # a really small number ……
OpenSEES前后处理程序ETO 及实例教程 WSP HONG KONG LTD Principle Engineer Dr. Chen Xuewei
主要内容(Main Content) ?ETO 简介及开发思路?ETO 主要功能介绍? OpenSEES 实例教程 15分钟简单介绍
ETO 简介及开发思路 ETO 是一款具有与ETABS 交互接口的OpenSees 前后处理程序。 节点信息质量源信息截面信息 单元几何信息约束信息荷载信息线性材料信息单元定义 (如:纤维截面、Transformation ) 变形显示 分析参数设置 (如:分析类型、控制参数) 记录设置 (如:节点位移、单元内力、截面变形、模态) Tcl 脚本文件 针对实际问题进行适当修改即可提交运算 读取.OUT 结果文件
当前操作状态 菜单栏 快捷按钮工具栏 三维可视化界面导入ETABS生成的S2K文件ETO主要功能介绍
ETO主要功能介绍 单元定义 单元类型包括Elastic BeamColumn, Nolinear BeamColumn, Disp BeamColumn, Beam with Hinge, Truss。 界面类型包括工字型界面、矩形截 面。 材料选择。 GeoTransf包括Linear, P-Delta, Corotational。 截面配筋。 纤维划分定义。
分析类型包括Single Load Control, Single Displacement Control, Gravity+Pushover, Modal Analysis, Time History Analysis, D+L Time Hist Analysis。 控制工况包括荷载控制、位移控制。 材料选择。 非线性设置包括梁柱单元积分点数、 钢筋材料序号。 截面组装。 模态数量。 记录设置 所有节点位移。 所有框架单元内力。 非线性梁柱单元截面变形。 非线性梁柱单元截面应力-应变关系。 振型特征值。
转自https://www.wendangku.net/doc/b99076480.html,/zpz_87/blog/static/2746458420133149589288/ OpenSees中坐标转化的说明geomTransf 建立OpneSees三维模型时,有人往往在转换坐标geomTransf的设置上出错或犯晕,下面结合OpenSees官网说明如何进行局部到整体坐标的转化。 OpenSees中坐标变化本质上为了确定截面section在整体坐标系中的位置或方向,一般截面局部坐标如下图,即我们所说的强轴是y,弱轴为z(补充注明:此处强弱轴的说法是广义的,强轴不是根据截面惯性矩确定的,而是根据受弯方向确定),局部x沿着杆件的轴向。 I. 2维问题 命令格式:geomTransf Linear $transfTag <-jntOffset $dXi $dYi $dXj $dYj> 说明:在平面问题中,局部z始终指向平面外,当局部x沿着轴向时,y的方向根据右手法则自然就确定了。 II. 3维问题 命令格式:geomTransf Linear $transfTag $vecxzX $vecxzY $vecxzZ <-jntOffset $dXi $dYi $dZi $dXj $dYj $dZj> 说明:任何情况下局部坐标x始终沿着杆件的轴向,geomTransf命令中$vecxzX $vecxzY $vecxzZ这3个参数主要确定局部坐标z与整体坐标X、Y、Z的关系,理解为局部坐标zx平面中的向量在整体坐标中的投影,如下图所示。
所以由局部z的投影与XYZ轴夹角的余弦值即可确定上述命令中的三个参数,举例如下: 若局部z与整体坐标X方向相同,则$vecxzX $vecxzY $vecxzZ分别为1,0,0(即cos(0)=1); 若局部z与整体坐标Z方向相反,则$vecxzX $vecxzY $vecxzZ分别为0,0,-1(即cos(180)=-1); 依次类推,即可确定任何情况下geomTransf命令的参数,基于上述解释就更加容易理解官网的解释,如下: 1.Element 1 : tag 1 : vecxZ = zaxis geomTransf Linear 1 0 0 -1 1.Element 2 : tag 2 : vecxZ = y axis geomTransf Linear 2 0 1 0
wipe # Opensees dandun # #Units:kN, m, sec # # Start of model generation # --------------- # Create ModeBulider (with two -dimensions and 3 DOF/node) model basic -ndm 2 -ndf 3 # --------------- # tag X Y node 1 0.0 0.0 node 2 0.0 0.0 node 3 0.0 2.0 node 4 0.0 4.0 node 5 0.0 6.0 node 6 0.0 8.0 node 7 0.0 10.0 node 8 0.0 12.0 node 9 0.0 14.0 node 10 0.0 16.0 node 11 0.0 18.0 node 12 0.0 20.0 # --------------- # Fix supports at base of columns # tag DX DY RZ fix 1 1 1 1 # --------------- # Concrete tag fc ec0 fcu ecu # Core concrete (confined) uniaxialMaterial Concrete01 1 -25600.0 -0.00219 -17780.0 -0.01 #Cover concrete (unconfined) uniaxialMaterial Concrete01 2 -23400.0 -0.002 -0.0 -0.006 # STEEL # Reinforcing steel set fy 400000.0; #Yield stress set E 200000000.0;# Young's modulus # tag fy E0 b uniaxialMaterial Steel02 3 $fy $E 0.01 18.5 0.925 0.15 uniaxialMaterial Elastic 11 29043600 uniaxialMaterial Elastic 12 12326600 uniaxialMaterial Elastic 13 587247596 #Define cross-section for nonlinear columns # ------------------- # set some parameters set colWidth 8.18 set colDepth 4.28 set cover 0.05 set As 0.00049 # some variables derived from the parameters set y1 [expr $colDepth/2.0] set z1 [expr $colWidth/2.0] section Fiber 1 { # Create the concrete core fibers patch rect 1 20 30 [expr $cover -$y1] [expr $cover -$z1] [expr $y1 -$cover] [expr $z1 -$cover] # Create the concrete cover fibers (top, bottom, left, right) patch rect 2 20 5 [expr -$y1] [expr $z1 -$cover] $y1 $z1 patch rect 2 20 5 [expr -$y1] [expr -$z1] $y1 [expr $cover -$z1] patch rect 2 5 10 [expr -$y1] [expr $cover -$z1] [expr $cover -$y1] [expr $z1 -$cover] patch rect 2 5 10 [expr $y1 -$cover] [expr $cover -$z1] $y1 [expr $z1 -$cover] # Create the reinforcingfibers (left, middle, right) layer straight 3 175 $As [expr $y1 -$cover] [expr $z1 -$cover] [expr $y1 -$cover] [expr $cover -$z1] layer straight 3 175 $As [expr $cover -$y1] [expr $z1 -$cover] [expr $cover -$y1] [expr $cover -$z1] layer straight 3 115 $As [expr $y1 -$cover] [expr $z1 -$cover] [expr $cover -$y1] [expr $z1 -$cover] layer straight 3 115 $As [expr $y1 -$cover] [expr $cover -$z1] [expr $cover -$y1] [expr $cover -$z1] }
OpenSEES 解题一般规律、技巧总结单位 OpenSEES 中是可以用公制单位(N,m )的(而并不是像某些文章中说的“ OpenSees 默认为英制单位”)。实际上我认为OpenSEES 中并没有什么默认单位,只要编程者自己保持单位一致就行;这点类似于SAP2000 的风格。 建模顺序 做事要讲究顺序,OpenSEES 建模亦如是:必须先定义材料才能离散截面(因为离散截面时要对所划分的截面指定材料属性)。 与之类似的,必须先定义(离散)截面,才能定义非线性梁柱单元(因为定义非线性梁柱单元时要指定单元截面)。 关于BandSPD 求解方式 官网关于BandSPD 方程形式的评价: "This is a good choice for most small size models. " 并且后面紧跟了一句: "The equations have to be numbered so the widely used RCM (Reverse Cuthill-McKee) numberer is used. " 可见numberer 类型不是随便选,而是要根据方程类型来决定的!
(不过直至 M 乍业做完,我对 nu mberer, system, test, algorithm, an alysis geomTra nsf, con strai nts )等求解控制命令还是一知半解!我觉得要想弄明白这些命令一 OpenSEES 中默认的计算精度比较高! “0.1000000000000001 M 0.1 ” :(自行总结,未找到官方说明)这是一个真实的故事: 我曾在程序中自以为是的将一连串相邻均只有 0.1左右的数的差强行赋值为 0.1,而没有采 用循环命令将两数作差并将结果赋给新变量一一其中即有这样的强行截断!我以为小数点 后都n 位了,即使我带 着它最后也会被系统截断,还不如我直接预处理来得清爽!没想到 这样做直接导致计算不收敛!真是失之毫厘谬 以千里!可见在 OpenSEES 中默认的计算精 度比较高! 后来我还在老师给的一份范例程序( Silvia Mazzo ni & Frank McKe nna, 2006 )中发现了 这么一段: set Ubig 1.e10; # a really large nu mber set Usmall [expr 1/$Ubig]; # a really small nu mber (还包括 —得先回头好好翻翻有限元和数值分析的书了!
1、利用零长单元模拟阻尼, uniaxialMaterial Elastic 1 6、8098e6; uniaxialMaterial Viscous 2 3、24e5 1; uniaxialMaterial Parallel 3 3 5; element zeroLength 1 $iNode $jNode -mat 3 -dir 1; 通常有两种方式: (1)truss element and viscous material、(桁架单元与阻尼材料) (2)force-based beam-column element and Maxwell material(基于力得梁柱单元与Maxwell材料)。 -、如何运行OpenSEES 有三种方法可以执行OpenSees/Tcl命令: 1、interactive交互式 直接将命令输入Prompt。 2、执行文件输入 这种方法就是最常用得一种,以source inputfile、tcl方式执行已写好得外部命令文件。 3、Batch模式 即以Opensees inputFile、tcl方式在MS-DOS/Unix promt中运行。 二、定义单位与常数 在编写一个较大得Opensees命令时。最好先定义好单位及常数。在Opensees中,编译器不能自行转换单位。所以一开始就要先定义好。 单位定义包括两部分:首先定义基本单位;再定义合成单位。其中基本单位要相互独立。同时,在定义单位时,既可以按国际公制单位,也可以按私制单位。因些在单位定义文件中可能就是混合得。我个人建议,还就是采用国际公制单位较好。像国外常用英制单位。很不习惯。 对于一些常数,如与g等常数要事先定义好。 在定义这些单位时所用得命令就是“set”。
【转载】[教程]OpenSEES超简单易懂的入门第一课 转载一下华南理工陈学伟博士关于OpenSEES的教程,呵呵。 陈博士还搞了OpenSEES的中国区论坛,https://www.wendangku.net/doc/b99076480.html,/index.php,不错的。 现在很多研究生学习OpenSEES这个程序,非常地火。OpenSEES是美国伯克利大学又一个成功的程序。正确的简写是OpenSEES,"pen" 是小写的,呵呵。OpenSEES我也学得不好,但有网友要求我写一个简单的操作过程,我觉得应没多大问题,但是有个要求,请高手不要拍砖!献丑了。 一般在我的博客可以下载到OpenSEES,网址是: https://www.wendangku.net/doc/b99076480.html,/article.asp?id=21 下载下来就以下几个文件,有范例,程序及说明书,这几个不到30MB的程序就可以完成非线性分析,是不是很神奇呢? (1)安装文件: OpenSEES要安装的文件只有一个,就是ActiveTcl8.4.6.1-win32-ix86-99631.exe。双击它进行安装,不断地按NEXT就可以装完了,如图所示。 (2)打开OpenSEES: OpenSEES.exe文件是不需要安装的,它是一个基于tcl的Dos窗口程序,双击打开,如图所示。不像现在的商业软件,一打开就是图形界面,这一点很多研究生接受不了,但对于早期接触有限元或电脑的人一点都不陌生。如以前的3DS,AutoCAD的第一个版本就是Dos 程序,Abaqus 1.0与Sap4都是Dos程序,就是通过如下窗口输入命令流的。当然,学习OpenSEES的人不需要一个个字打命令流进去,这个容易出错,我们可以学习Ansys做法,把ADPL写到一个文本文件然后Load进去。正是这样,以下的方法就是:先写tcl文件,再Load进OpenSEES做计算。
最新文件---------------- 仅供参考--------------------已改成-----------word文本 --------------------- 方便更改 1、利用零长单元模拟阻尼, uniaxialMaterial Elastic 1 6.8098e6; uniaxialMaterial Viscous 2 3.24e5 1; uniaxialMaterial Parallel 3 3 5; element zeroLength 1 $iNode $jNode -mat 3 -dir 1; 通常有两种方式: (1)truss element and viscous material.(桁架单元和阻尼材料) (2)force-based beam-column element and Maxwell material(基于力的梁柱单元和Maxwell 材料)。 -、如何运行OpenSEES 有三种方法可以执行OpenSees/Tcl命令: 1、interactive交互式 直接将命令输入Prompt。 2、执行文件输入 这种方法是最常用的一种,以source inputfile.tcl方式执行已写好的外部命令文件。 3、Batch模式 即以Opensees inputFile.tcl方式在MS-DOS/Unix promt中运行。 二、定义单位和常数 在编写一个较大的Opensees命令时。最好先定义好单位及常数。在Opensees中,编译器不能自行转换单位。所以一开始就要先定义好。 单位定义包括两部分:首先定义基本单位;再定义合成单位。其中基本单位要相互独立。同时,在定义单位时,既可以按国际公制单位,也可以按私制单位。因些在单位定义文件中可能是混合的。我个人建议,还是采用国际公制单位较好。像国外常用英制单位。很不习惯。 对于一些常数,如 和g等常数要事先定义好。 在定义这些单位时所用的命令是“set”。
OPENSEES Opensees模型 OpenSEES中有限元对象被划分成更多的子对象,其中包括节点对象、材料对象、截面对象、单元对象、荷载对象和约束对象等,并且为其子对象提供了多种不同的选择,包括不同的材料类型,截面形式,荷载模式以及约束方式等,再由它们组合成为有限元模型对象。在程序中建立子对象的命令主要有:Node、Mass、Material、Section、Element、LoadPattern、TimeSeries、Transformation、Block和Constraint等等。通过上述命令,我们可以分别确定对象中各节点的位置、节点集中质量、材料本构关系、截面恢复力模型、单元类型、外加荷载模式、几何坐标转换类型和约束形式等。这些命令构建了有限元模型相应的子对象,由这些子对象组合构成有限元模型对象ModelBuilder。 纤维模型 纤维模型是指将纤维截面赋予梁柱构件(即定义构件的每一截面为纤维截面),纤维截面是将构件截面划分成很多小纤维(包括钢筋纤维和混凝土纤维)对每一根纤维只考虑它的轴向本构关系,且各个纤维可以定义不同的本构关系。纤维模型假定构件的截面在变形过程中始终保持为平面,这样只要知道构件截面的弯曲应变和轴向应变就可以得到截面每一根纤维的应变,从而可以计算得到截面的刚度。纤维模型能很好的模拟构件的弯曲变形和轴向变形,但不能模拟构件的剪切非线性和扭曲非线性。 构件 零长度构件 可以赋予零长度构件BARSLIPMaterial(这种材料的本构关系可以精确模拟循环加载时在构件节点处由于钢筋的滑移和混凝土的开裂所引起的构件的刚度退化和强度退化现象)来模拟构件节点处的变形,另外用Bond-SP01Material可以模拟节点处钢筋的应力渗透现象(节点处钢筋还没有整体滑移)所引起的构件的强度和刚度变化。
实例6 框架结构弹性时程分析 1)问题描述: 本例仍采用实例4(实例5)的框架结构,为了方便对比,改采用弹性截面,主要进行弹性时程分析,材料为弹性,时程分析即动力分析。结构荷载情况与实例4相同(侧向力荷载不需要施加,实际侧向力为地面加速度)。计算结构在地震作用下的响应(主要提取位移结果)。(重力荷载代表值组合为1.0×DEAD+0.5×LIVE)。 注意:上述实例讲到了质量矩阵(质量源)的定义,刚度矩阵通过结构几何与材料属性得到,那么接下来只需要定义了结构阻尼,就可以进行结构动力分析,即时程分析。 2)ETABS模型建模 (1) 结构模型与实例5相同,相关建模细节请看实例5。
图 ETABS建立框架的几何模型 (2) 为了对比OPENSEES弹性时程的分析结果,在ETABS模型同样进行弹性时程分析。 为了进行弹性时程分析,需要输入地震波,(本例只进行单向地震分析)。 (3) 地震波数据导入:(实例的地震波放在光盘“/EXAM06/ETABS/”目录)本算例采用的为单向地震波,地震波文件为:GM1X.txt,通过EXCEL图表画出整个地震波时程曲线如下图所示。(该时程的时间间隔为0.02 sec) 图时程曲线GM1X 从时程曲线可知,曲线最大值为3621(该值为无单位数)。 【Define】→【Time History Functions】,如下窗口所示。选取【Function from File】点击【Add New Functions】,弹出以下窗口。
图 GM1X时程曲线定义窗口 (4) 弹性时程分析工况定义,点击【Define】→【Time History Cases】,如下窗口所示。 按窗口的内容填写数据,点击【Modal Damping】右边的按钮【Modify/Show】,可以看到阻尼比的填写框,由于是混凝土结构,阻尼比为0.05。 图时程分析工况内容 窗口中输入的内容简介如下: 荷载步数为2000步,每步代表0.02s,即总分析时间长度为40s,主轴1方向的地 震波时程数据为 GM1X,其放大倍数为0.138(单位为mm/s2),也就是说整个地 震波的最大地面加速度为3621×0.138=500 mm/s2,即50gal,属于小震量级。(5) 完成上述步骤后建立完ETABS模型。
1、瑞利阻尼 在OPENSEES 中,结构采用瑞利(Rayleigh )阻尼,即阻尼矩阵的大小与结构的质量矩阵,刚度矩阵都相关,瑞利阻尼的计算公式如下,阻尼与刚度质量的关系如下图所示。 [][][]01c a m a k =+ 式中,ξ为阻尼比,a 0为质量相关系数,a 1为刚度相关系数,[c]为阻尼矩阵,[m]为质量矩阵,[k]为刚度矩阵;ωm 、ωn 为结构两个主振型的圆频率,由于 OPENSEES 能够直接求解振型的特征值,那么特征值与圆频率的关系:ω。 命令流的解读如下: set xDamp 0.05 ;————设置阻尼比为0.05 set nEigenI 1;————主振型1为第1振型 set nEigenJ 2;————主振型2为第2振型 set lambdaN [eigen [expr $nEigenJ]];————求解两阶振型即可 set lambdaI [lindex $lambdaN [expr $nEigenI-1]];————提取第1阶特征值 set lambdaJ [lindex $lambdaN [expr $nEigenJ-1]];————提取第2阶特征值 set omegaI [expr pow ($lambdaI,0.5)];——— —从特征值求圆频率 set omegaJ [expr pow($lambdaJ,0.5)];————从特征值求圆频率 set alphaM [expr $xDamp*(2*$omegaI*$omegaJ)/($omegaI+$omegaJ)]; ————alphaM 为a 0,即质量相关系数;
set betaKcurr [expr 2.*$xDam p/($omegaI+$omegaJ)]; ———betaKcurr 为a1,即刚度相关系数; rayleigh $alphaM $betaKcurr 0 0———定义瑞利阻尼,只需要填写a0、a1,其它值为0。 2、OPENSEES 质量源 OPENSEES基本上采用结点质量源的形式,与大部分有限元程序一样,也就是每个结点有6 个自由度,6 个自由度上都有广义质量,平动方向UX,UY,UZ 称为质量,而转动方向RX,RY,RZ称为转动惯量(除了刚体,基本上结点没有转动惯量)。普通结点只有UX,UY,UZ的平动质量,且三个值是相等的。puts "mass" mass $NUM $MUX $MUY $MUZ $MRX $MRY $MRZ 其中,$NUM为结点编号,$MUX $MUY $MUZ $MRX $MRY $MRZ代表各个自由度的质量,以单位制(N,mm)的规定,质量的单位应该为ton(吨)。普通结点的质量定义如下: mass $NUM $M $M $M 0 0 0 3、模态分析 模态分析的命令流与普通静力分析的命令流最大的区别在于记录与分析设置。记录命令如下: puts "recorder" recorder Node -file eigen1_node0.out -time -nodeRange 1 28 -dof 1 2 3 "eigen 1" 该命令流用于输出振型位移,即振型形状。 其中,eigen1_node0.out 为输出振型位移的文件名,-dof 1 2 3 代输出的自由度;"eigen1"代表输出的为“第1 振型”。 模态分析设置的命令流: set numModes 12 set lambda [eigen $numModes] set period "Periods.txt" set Periods [open $period "w"] puts $Periods " $lambda"
# SET UP ---------------------------------------------------------------------------- wipe; # clear opensees model model basic -ndm 3 -ndf 6; # 2 dimensions, 3 dof per node file mkdir Data; # create data directory logFile Data/errorFile.txt; node 1 0 0 0; #两个节点定义一竖向悬臂梁 node 2 0 0 10; # Single point constraints -- Boundary Conditions fix 1 1 1 1 1 1 1; #固定1节点 mass 2 100 100 100 0. 0. 0.; # Define ELEMENTS ------------------------------------------------------------- #坐标转换 geomTransf Linear 1 0 -1 0; element elastic BeamColumn 1 1 2 0.01 1.0e8 1.0e20 2.0e-5 1.0e-5 1.0e-5 1; #element elastic BeamColumn $eleTag $iNode $jNode $A $E $G $J $Iy $Iz $transfTag puts "[eigen 1]"; set eigenvalues [eigen 1]; puts "eigenvalues ... [lindex $eigenvalues 0]"; set lambda1 [expr [lindex $eigenvalues 0]]; set omega1 [expr pow($lambda1,0.5)]; set T1 [expr 2*3.1416/$omega1]; puts "period ... $T1"; # Define RECORDERS --------------------------------------------------------- recorder Node -file Data/node2eig.txt -time -node 2 -dof 1 2 "eigen 1"; # define GRA VITY ----------------------------------------------------------- pattern Plain 1 Linear { load 2 0. 0. -1000. 0. 0. 0. }; constraints Plain; numberer Plain; system BandGeneral; test NormDispIncr 1.0e-8 10; algorithm Newton; integrator LoadControl 0.1; analysis Static analyze 1100; loadConst -time 0.0; 这个例子我从在仿真论坛里发过,但那个论坛需要邀请码才能注册,新手一般是没有这个论坛的账号的。这个例子对新手有用,所以也发到这里来吧。类似的例子还有很多,不过最近比较忙,等做完论文了,我会整理一下发上来。 现在就上面的例子做一些说明: