2024年1月21日发(作者:)
if { [file exists output] == 0 } {
file mkdir output
}
建立活页夹
model basic -ndm 2 -ndf 3
node 5 0.0 162.0 -mass $m1 $m1 0.0
约束
品质
# tag fy E0 b
uniaxialMaterial Steel01 3 $fy $E 0.01
Steel02 Material tag fy E0 b $R0 $CR1 $CR2 $R0=between 10 and 20, $cR1=0.925,
$cR2=0.15
CONCRETE tag f'c ec0 f'cu ecu
# Core concrete (confined)
uniaxialMaterial Concrete01 1 -6.0 -0.004 -5.0 -0.014
set fc1U $fc; # UNCONFINED concrete (todeschini parabolic
model), maximum stress
set eps1U -0.003; # strain at maximum strength of unconfined
concrete
set fc2U [expr 0.2*$fc1U]; # ultimate stress
set eps2U -0.01; # strain at ultimate stress
set lambda 0.1; # ratio between unloading slope at $eps2 and
initial slope $Ec
# tensile-strength properties
set ftU [expr -0.14*$fc1U]; # tensile strength +tension
set Ets [expr $ftU/0.002]; # tension softening stiffness
uniaxialMaterial Concrete02 $IDconcU $fc1U $eps1U $fc2U $eps2U $lambda
$ftU $Ets
section Fiber 1 {
# Create the concrete core fibers
patch rect 1 10 1 [expr $cover-$y1] [expr $cover-$z1] [expr $y1-$cover]
[expr $z1-$cover]
# Create the concrete cover fibers (top, bottom, left, right)
patch rect 2 10 1 [expr -$y1] [expr $z1-$cover] $y1 $z1
patch rect 2 10 1 [expr -$y1] [expr -$z1] $y1 [expr $cover-$z1]
patch rect 2 2 1 [expr -$y1] [expr $cover-$z1] [expr $cover-$y1]
[expr $z1-$cover]
patch rect 2 2 1 [expr $y1-$cover] [expr $cover-$z1] $y1 [expr
$z1-$cover]
# Create the reinforcing fibers (left, middle, right)
layer straight 3 3 $As [expr $y1-$cover] [expr $z1-$cover] [expr
$y1-$cover] [expr $cover-$z1]
layer straight 3 2 $As 0.0 [expr $z1-$cover] 0.0 [expr $cover-$z1]
layer straight 3 3 $As [expr $cover-$y1] [expr $z1-$cover] [expr
$cover-$y1] [expr $cover-$z1]
} 定义截面
tag ndI ndJ A E Iz transfTag
element elasticBeamColumn 3 3 4 360 4030 8640 2 弹性梁柱、
定义(二维)
element elasticBeamColumn $eleTag $iNode $jNode $A $E $G $J $Iy $Iz $transfTag(三维)
# tag ndI ndJ nsecs secID transfTag
element nonlinearBeamColumn 1 1 3 $np 1 1 (非线性的梁柱单元)
pattern Plain 1 1 {
eleLoad -ele 13 14 15(单元编号) -type –beamUniform(均布线荷载) -$w1(沿着y轴的线荷载大小)
eleLoad -ele 16 17 18 -type -beamUniform -$w2
eleLoad -ele 19 20 21 -type -beamUniform -$w3
}单元加载
pattern Plain 1 "Linear" {
# Create nodal loads at nodes 3 & 4
# nd FX FY MZ
load 3 0.0 [expr -$P] 0.0
load 4 0.0 [expr -$P] 0.0
}节点载入
geomTransf Linear 1; # beams
geomTransf PDelta 2; # columns
(二维坐标下的坐标转换)
recorder Node -file Data/ -time -node 2 -dof 1 2 3 disp;
{节点记录}
recorder Element -file -ele 1 4 7 10 forces
(单元的记录)
# Define RECORDERS
-------------------------------------------------------------
recorder Node -file Data/ -time -node 2 -dof 1 2 3 disp;
# displacements of free nodes
recorder Node -file Data/ -time -node 1 -dof 1 2 3 disp;
# displacements of support nodes
recorder Node -file Data/ -time -node 1 -dof 1 2 3 reaction; # support reaction
recorder Drift -file Data/ -time -iNode 1 -jNode 2 -dof 1
-perpDirn 2 ; # lateral drift
recorder Element -file Data/ -time -ele 1 globalForce;
# element forces -- column
recorder Element -file Data/ -time -ele 1 deformations;
# element deformations -- column
set lambda [eigen $numModes];(求解$numModes阶振型)
set eigenValue [lindex $eigenValues 0]提取第一节特征值
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 为a0,即质量相关系数;
。
set betaKcurr [expr 2.*$xDam p/($omegaI+$omegaJ)]; — ——betaKcurr 为a1,即刚度相关系数;
rayleigh $alphaM $betaKcurr 0 0———定义瑞利阻尼,只需要填写a0、a1,其它值为0。
set alphaM 0.; # M-prop. damping; D = alphaM*M
set betaKcurr 0.; # K-proportional damping;
+beatKcurr*KCurrent
set betaKcomm [expr 2.*$xDamp/($omega)]; # K-prop. damping
parameter; +betaKcomm*KlastCommitt
set betaKinit 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_ -time -nodeRange 1 28 -dof 1 2 3 "eigen 1"
该命令流用于输出振型位移,即振型形状。
其中,eigen1_ 为输出振型位移的文件名,-dof 1 2 3 代输出的自由度;"eigen 1"代表输出的为“第1 振型”。
for { set k 1 } { $k <= $numModes } { incr k } {
recorder Node -file [format "modes/mode%" $k] -nodeRange 1 6 -dof
1 2 3 "eigen $k"
}
模态分析设置的命令流:
set numModes 12
set lambda [eigen $numModes]
set period ""
set Periods [open $period "w"]
puts $Periods " $lambda"
close $Periods
record
其中,set lambda [eigen $numModes] ,代表计算n 阶振型,将特征值计算结果存为lambda 数组;
set period "" ,定义输出文本文件的名字为“ ”;用于存计算的特征值。
set Periods [open $period "w"],代表打开文本文件进行记录。
puts $Periods " $lambda",代表记录特征值数据至文本当中。
close $Periods,代表记录完成,关闭文本文件。
record 代表记录命令。
4、单向地震波资料设置的命令流如下:
set IDloadTag 1001; ————地震波工况号为 1001
set iGMfile ""; ————地震波数据文件名为 ,需要与 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 -file Path $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; BandGeneral
————采用Umfpack 自由度排列
test EnergyIncr 1.0e-4 200; NormDispIncr 1.0e-8 100
————采用能量收敛准则
algorithm Newton————采用牛顿迭代法
integrator Newmark 0.5 0.25————采用 Newmark法对时间进行离散
integrator DisplacementControl 2 1 0.1(二节点第一个自由度位移每次增加0.1)
analysis Transient————采用时程分析
analyze 1000 0.02————分析步数为 1000,时间间隔为0.02s,即分析20s。
set maxU 6.0; # Max displacement
set numSteps [expr int($maxU/$dU)]
# Perform the analysis
analyze $numSteps
recorder Node -file Data/ -time -node 2 -dof 1 2 3 disp;
# displacements of free nodes
recorder Node -file Data/ -time -node 1 -dof 1 2 3 disp;
# displacements of support nodes
recorder Node -file Data/ -time -node 1 -dof 1 2 3 reaction; # support reaction
recorder Drift -file Data/ -time -iNode 1 -jNode 2 -dof 1
-perpDirn 2 ; # lateral drift
recorder Element -file Data/ -time -ele 1 globalForce;
# element forces -- column
recorder Element -file Data/ -time -ele 1 deformations;
set maxU 10.0; # Max displacement
set controlDisp 0.0;
set ok 0;
while {$controlDisp < $maxU && $ok == 0} {
set ok [analyze 1]
set controlDisp [nodeDisp 3 1]
if {$ok != 0} {
puts "... trying an initial tangent iteration with Newton"
test NormDispIncr 1.0e-8 4000 0
algorithm ModifiedNewton -initial
set ok [analyze 1]
test NormDispIncr 1.0e-8 10
algorithm Newton
}
}
if {$ok != 0} {
puts "Pushover analysis FAILED"
} else {
puts "Pushover analysis completed SUCCESSFULLY"
}
由位移控制换算法
set du [expr 0.1*$mm]
integrator DisplacementControl 2 1 $du
set maxU [expr 200*$mm]
set controlDisp 0.0;
set ok 0;
while {$controlDisp < $maxU && $ok == 0} {
set ok [analyze 1]
set controlDisp [nodeDisp 2 1]
if { $ok != 0} {
puts "... trying an initial tangent iteration with Newton"
test NormDispIncr 1.0e-6 4000
algorithm ModifiedNewton -initial
set ok [analyze 1]
test NormDispIncr 1.0e-6 400
algorithm Newton
}
if {$ok != 0} {
puts "Trying Broyden .."
algorithm Broyden 8
set ok [analyze 1]
algorithm Newton
}
if {$ok != 0} {
puts "Trying NewtonWithLineSearch .."
algorithm NewtonLineSearch .8
set ok [analyze 1 ]
algorithm Newton
}
}
set accelSeries "Series -dt 0.01 -filePath -factor 1"; #
define acceleration vector from file (dt=0.01 is associated with the input
file gm)
pattern UniformExcitation 2 1 -accel $accelSeries; # define
where and how (pattern tag, dof) acceleration is applied
rayleigh 0. 0. 0. [expr 2*0.02/pow([eigen 1],0.5)]; # set damping
based on first eigen mode
时程分析换算法
set Nsteps [expr int($TmaxAnalysis/$DtAnalysis)];
set ok [analyze $Nsteps $DtAnalysis]; # actually perform
analysis; returns ok=0 if analysis was successful
if {$ok!= 0} { ; # if analysis was not successful.
# change some analysis parameters to achieve convergence
# performance is slower inside this loop
# Time-controlled analysis
set ok 0;
set controlTime [getTime];
while {$controlTime < $TmaxAnalysis && $ok == 0} {
set ok [analyze 1 $DtAnalysis]
set controlTime [getTime]
if {$ok != 0} {
puts "Trying Newton with Initial Tangent .."
test NormDispIncr $Tol 1000 0
algorithm Newton -initial
set ok [analyze 1 $DtAnalysis]
test $TestType $Tol $maxNumIter 0
algorithm $algorithmType
}
if {$ok != 0} {
puts "Trying Broyden .."
algorithm Broyden 8
set ok [analyze 1 $DtAnalysis]
algorithm $algorithmType
}
if {$ok != 0} {
puts "Trying NewtonWithLineSearch .."
algorithm NewtonLineSearch .8
set ok [analyze 1 $DtAnalysis]
algorithm $algorithmType
}
}
}; # end if ok !0
puts "Ground Motion Done. End Time: [getTime]"
(定义刚性板)
# define Rigid Floor Diaphragm
set RigidDiaphragm ON ; # options: ON, OFF. specify this before the
analysis parameters are set the constraints are handled differently.
set Xa [expr ($X2+$X1)/2]; # mid-span coordinate for rigid diaphragm
set Za [expr ($Z2+$Z1)/2];
# rigid-diaphragm nodes in center of each diaphram
node 1121 $Xa $Y2 $Za; # master nodes for rigid diaphragm -- story
2, bay 1, frame 1-2
node 1131 $Xa $Y3 $Za; # master nodes for rigid diaphragm -- story
3, bay 1, frame 1-2
node 1141 $Xa $Y4 $Za; # master nodes for rigid diaphragm -- story
4, bay 1, frame 1-2
# Constraints for rigid diaphragm master nodes
fix 1121 0 1 0 1 0 1
fix 1131 0 1 0 1 0 1
fix 1141 0 1 0 1 0 1
# ------------------------define Rigid Diaphram, dof 2 is normal to floor
set perpDirn 2;
rigidDiaphragm $perpDirn 1121 121 122 221 222; # level 2
rigidDiaphragm $perpDirn 1131 131 132 231 232; # level 3
rigidDiaphragm $perpDirn 1141 141 142 241 242; # level 4
三维的需要考虑扭转刚度,把扭转刚度定义为无穷大就不考虑扭转刚度啦
uniaxialMaterial Elastic $SecTagTorsion $Ubig
section Aggregator $ColSecTag $SecTagTorsion T -section
$ColSecTagFiber
section Aggregator $BeamSecTag $SecTagTorsion T -section
$BeamSecTagFiber
section Aggregator $GirdSecTag $SecTagTorsion T -section
$GirdSecTagFiber
三维坐标转化
set IDColTransf 1; # all columns
set IDBeamTransf 2; # all beams
set IDGirdTransf 3; # all girders
set ColTransfType Linear ; # options, Linear PDelta
Corotational
geomTransf $ColTransfType $IDColTransf 0 0 1 ; # only columns can
have PDelta effects (gravity effects)
geomTransf Linear $IDBeamTransf 0 0 1
geomTransf Linear $IDGirdTransf 1 0 0
零长度单元定义铰
同一个节点用两个不同的点来定义
再用equalDOF的命令确定这两个点上一样的自由度
再定义一个小的弹性模量
定义零长度单元来释放转动的约束
equalDOF 1 2 (1.2节点在X.Y方向一样)1 2
equalDOF 3 4 (3.4节点在Y方向一样) 2
uniaxialMaterial Elastic 3 1.0(定义材料3的弹性模量很小)
element zeroLength 2 1 2 -mat 3 -dir 3(定义零长度单元2,在1.2节点上是同一个点,采用弹性模量很小的材料来释放转动方向的约束)
constraints Plain
numberer Plain
system BandGeneral
test EnergyIncr 1.0e-6 200
algorithm Newton
analysis Static
##################################################################低周往复分析
integrator DisplacementControl 26 2 0.2
analyze 100
integrator DisplacementControl 26 2 -0.4
analyze 100
integrator DisplacementControl 26 2 0.6
analyze 100
integrator DisplacementControl 26 2 -0.8
analyze 100
integrator DisplacementControl 26 2 1.0
analyze 100
integrator DisplacementControl 26 2 -1.2
analyze 100
integrator DisplacementControl 26 2 1.4
analyze 100
integrator DisplacementControl 26 2 -1.4
analyze 100
integrator DisplacementControl 26 2 1.2
analyze 100
set Sa [lindex $8_85_6m [expr $t-1]]提取列表中的数
set SAlist {};建立一个列表
lappend SAlist $SA;把变量存储到列表中


发布评论