I performed a elastic pushover analysis by the following scrips.
But I could not get accurate result.
So, please check my scripts.
###################################################################################################
# Set Up & Source Definition
###################################################################################################
wipe all; # clear memory of past model definitions
model BasicBuilder -ndm 2 -ndf 3; # Define the model builder, ndm = #dimension, ndf = #dofs
source DisplayModel2D.tcl; # procedure for displaying a 2D perspective of model
source DisplayPlane.tcl; # procedure for displaying a plane in a model
###################################################################################################
# Define Analysis Type
###################################################################################################
# Define type of analysis: "pushover" = pushover; "dynamic" = dynamic
set analysisType "pushover";
if {$analysisType == "pushover"} {
set dataDir Concentrated-Pushover-Output; # name of output folder
file mkdir $dataDir; # create output folder
}
if {$analysisType == "dynamic"} {
set dataDir Concentrated-Dynamic-Output; # name of output folder
file mkdir $dataDir; # create output folder
}
###################################################################################################
# Define Building Geometry, Nodes, and Constraints
###################################################################################################
# define structure-geometry parameters
set NStories 3; # number of stories
set NBays 3; # number of frame bays
set WBay 6000.0; # bay width in (mm)
set HStory 4000.0; # 1st story height in (mm)
set HBuilding [expr $HStory*$NStories]; # height of building in(mm)
# calculate locations of beam/column joints:
set Pier1 0.0; # leftmost column line
set Pier2 [expr $Pier1 + $WBay];
set Pier3 [expr $Pier2 + $WBay];
set Pier4 [expr $Pier3 + $WBay];
set Floor1 0.0; # ground floor
set Floor2 [expr $Floor1 + $HStory];
set Floor3 [expr $Floor2 + $HStory];
set Floor4 [expr $Floor3 + $HStory];
# calculate joint offset distance for beam plastic hinges
#set phlat23 [expr 0.0]; # lateral dist from beam-col joint to loc of hinge on Floor 2
# calculate nodal masses -- lump floor masses at frame nodes
set g 9806.65; # acceleration due to gravity(mm/s^2)
set Floor2Mass 70000; # Mass of Floor 2 in (kg)
set Floor3Mass 70000; # Mass of Floor 3 in (kg)
set Floor4Mass 70000; # Mass of Floor 3 in (kg)
set MBuilding [expr $Floor2Mass + $Floor3Mass + $Floor4Mass];# total building Mass
set NodalMass2 [expr ($Floor2Mass) / (4.0)]; # mass at each node on Floor 2
set NodalMass3 [expr ($Floor3Mass) / (4.0)]; # mass at each node on Floor 3
set NodalMass4 [expr ($Floor4Mass) / (4.0)]; # mass at each node on Floor 4
set Negligible 1e-9; # a very smnumber to avoid problems with zero in (kg)
# define nodes and assign masses to beam-column intersections of frame
# command: node nodeID xcoord ycoord -mass mass_dof1 mass_dof2 mass_dof3
node 1 $Pier1 $Floor1;
node 2 $Pier2 $Floor1;
node 3 $Pier3 $Floor1;
node 4 $Pier4 $Floor1;
node 5 $Pier1 $Floor2 -mass $NodalMass2 $Negligible $Negligible;
node 6 $Pier2 $Floor2 -mass $NodalMass2 $Negligible $Negligible;
node 7 $Pier3 $Floor2 -mass $NodalMass2 $Negligible $Negligible;
node 8 $Pier4 $Floor2 -mass $NodalMass2 $Negligible $Negligible;
node 9 $Pier1 $Floor3 -mass $NodalMass3 $Negligible $Negligible;
node 10 $Pier2 $Floor3 -mass $NodalMass3 $Negligible $Negligible;
node 11 $Pier3 $Floor3 -mass $NodalMass3 $Negligible $Negligible;
node 12 $Pier4 $Floor3 -mass $NodalMass3 $Negligible $Negligible;
node 13 $Pier1 $Floor4 -mass $NodalMass4 $Negligible $Negligible;
node 14 $Pier2 $Floor4 -mass $NodalMass4 $Negligible $Negligible;
node 15 $Pier3 $Floor4 -mass $NodalMass4 $Negligible $Negligible;
node 16 $Pier4 $Floor4 -mass $NodalMass4 $Negligible $Negligible;
# constrain beam-column joints in a floor to have the same lateral displacement using the "equalDOF" command
# command: equalDOF $MasterNodeID $SlaveNodeID $dof1 $dof2...
set dof1 1; # constrain movement in dof 1 (x-direction)
equalDOF 5 6 $dof1; # Floor 2:
equalDOF 5 7 $dof1; # Floor 2:
equalDOF 5 8 $dof1; # Floor 2:
equalDOF 9 10 $dof1; # Floor 3:
equalDOF 9 11 $dof1; # Floor 3:
equalDOF 9 12 $dof1; # Floor 3:
equalDOF 13 14 $dof1; # Floor 4:
equalDOF 13 15 $dof1; # Floor 4:
equalDOF 13 16 $dof1; # Floor 4:
# assign boundary condidtions
# command: fix nodeID dxFixity dyFixity rzFixity
# fixity values: 1 = constrained; 0 = unconstrained
# fix the base of the building; pin P-delta column at base
fix 1 1 1 1;
fix 2 1 1 1;
fix 3 1 1 1;
fix 4 1 1 1;
###################################################################################################
# Define Section Properties and Elements
###################################################################################################
# define material properties
set Es 205000000.0; # steel Young's modulus((kg*mm/s^2)*1/mm^2)
# define column section
set Acol_c 17190; # cross-sectional area(mm^2)
set Icol_c 398000000; # moment of inertia(mm^4)
# define beam section
set Abeam_g 9543; # cross-sectional area(mm^2)
set Ibeam_g 329000000; # moment of inertia(mm^4)
# set up geometric transformations of element
set PDeltaTransf 1;
set LinearTransf 2;
geomTransf PDelta $PDeltaTransf; # columns
geomTransf Linear $LinearTransf; # beams
# define elastic column elements using "element" command
# command: element elasticBeamColumn $eleID $iNode $jNode $A $E $I $transfID
# eleID convention: "1xy" where 1 = col, x = Pier #, y = Story #
# Columns Story 1
element elasticBeamColumn 111 1 5 $Acol_c $Es $Icol_c $PDeltaTransf;
element elasticBeamColumn 121 2 6 $Acol_c $Es $Icol_c $PDeltaTransf;
element elasticBeamColumn 131 3 7 $Acol_c $Es $Icol_c $PDeltaTransf;
element elasticBeamColumn 141 4 8 $Acol_c $Es $Icol_c $PDeltaTransf;
# Columns Story 2
element elasticBeamColumn 112 5 9 $Acol_c $Es $Icol_c $PDeltaTransf;
element elasticBeamColumn 122 6 10 $Acol_c $Es $Icol_c $PDeltaTransf;
element elasticBeamColumn 132 7 11 $Acol_c $Es $Icol_c $PDeltaTransf;
element elasticBeamColumn 142 8 12 $Acol_c $Es $Icol_c $PDeltaTransf;
# Columns Story 3
element elasticBeamColumn 113 9 13 $Acol_c $Es $Icol_c $PDeltaTransf;
element elasticBeamColumn 123 10 14 $Acol_c $Es $Icol_c $PDeltaTransf;
element elasticBeamColumn 133 11 15 $Acol_c $Es $Icol_c $PDeltaTransf;
element elasticBeamColumn 143 12 16 $Acol_c $Es $Icol_c $PDeltaTransf;
# define elastic beam elements
# eleID convention: "2xy" where 2 = beam, x = Bay #, y = Floor #
# Beams Story 1
element elasticBeamColumn 212 5 6 $Abeam_g $Es $Ibeam_g $LinearTransf;
element elasticBeamColumn 222 6 7 $Abeam_g $Es $Ibeam_g $LinearTransf;
element elasticBeamColumn 232 7 8 $Abeam_g $Es $Ibeam_g $LinearTransf;
# Beams Story 2
element elasticBeamColumn 213 9 10 $Abeam_g $Es $Ibeam_g $LinearTransf;
element elasticBeamColumn 223 10 11 $Abeam_g $Es $Ibeam_g $LinearTransf;
element elasticBeamColumn 233 11 12 $Abeam_g $Es $Ibeam_g $LinearTransf;
# Beams Story 3
element elasticBeamColumn 214 13 14 $Abeam_g $Es $Ibeam_g $LinearTransf;
element elasticBeamColumn 224 14 15 $Abeam_g $Es $Ibeam_g $LinearTransf;
element elasticBeamColumn 234 15 16 $Abeam_g $Es $Ibeam_g $LinearTransf;
# display the model with the node numbers
DisplayModel2D NodeNumbers
############################################################################
# Recorders
############################################################################
# record drift histories
# drift recorder command: recorder Drift -file $filename -time -iNode $NodeI_ID -jNode $NodeJ_ID -dof $dof -perpDirn $Record.drift.perpendicular.to.this.direction
recorder Drift -file $dataDir/Drift-Story1.out -time -iNode 1 -jNode 5 -dof 1 -perpDirn 2;
recorder Drift -file $dataDir/Drift-Story2.out -time -iNode 5 -jNode 9 -dof 1 -perpDirn 2;
recorder Drift -file $dataDir/Drift-Story3.out -time -iNode 9 -jNode 13 -dof 1 -perpDirn 2;
recorder Drift -file $dataDir/Drift-Roof.out -time -iNode 1 -jNode 13 -dof 1 -perpDirn 2;
# record floor displacements
recorder Node -file $dataDir/Disp.out -time -node 1 5 9 13 -dof 1 disp;
# record base shear reactions
recorder Node -file $dataDir/Vbase.out -node 1 2 3 4 -dof 1 reaction;
# record story 1 column forces in global coordinates
recorder Element -file $dataDir/Fcol111.out -time -ele 111 force;
recorder Element -file $dataDir/Fcol112.out -time -ele 112 force;
recorder Element -file $dataDir/Fcol113.out -time -ele 113 force; #output date : Time - Shear force of i-node (Vi) - Axial force at i-node (Pi) - moment at i-node (Mi) - Vj - Pj - Mj
#######################################################################################
# #
# Analysis Section #
# #
#######################################################################################
############################################################################
# Eigenvalue Analysis
############################################################################
set pi [expr 2.0*asin(1.0)]; # Definition of pi
set nEigenI 1; # mode i = 1
set nEigenJ 2; # mode j = 2
set lambdaN [eigen [expr $nEigenJ]]; # eigenvalue analysis for nEigenJ modes
set lambdaI [lindex $lambdaN [expr 0]]; # eigenvalue mode i = 1
set lambdaJ [lindex $lambdaN [expr $nEigenJ-1]]; # eigenvalue mode j = 2
set w1 [expr pow($lambdaI,0.5)]; # w1 (1st mode circular frequency)
set w2 [expr pow($lambdaJ,0.5)]; # w2 (2nd mode circular frequency)
set T1 [expr 2.0*$pi/$w1]; # 1st mode period of the structure
set T2 [expr 2.0*$pi/$w2]; # 2nd mode period of the structure
puts "T1 = $T1 s"; # display the first mode period in the command window
puts "T2 = $T2 s"; # display the second mode period in the command window
#NodeEigenvector Command : nodeEigenvector $nodeTag $eigenvector $dof
set eigenvector111 [expr [nodeEigenvector 1 1 1] / [nodeEigenvector 13 1 1]]; #gives the value of eigenvector that corresponds to mode 1 at node 1 in dof 1
set eigenvector511 [expr [nodeEigenvector 5 1 1] / [nodeEigenvector 13 1 1]]; #gives the value of eigenvector that corresponds to mode 1 at node 5 in dof 1
set eigenvector911 [expr [nodeEigenvector 9 1 1] / [nodeEigenvector 13 1 1]]; #gives the value of eigenvector that corresponds to mode 1 at node 9 in dof 1
set eigenvector1311 [expr [nodeEigenvector 13 1 1] / [nodeEigenvector 13 1 1]]; #gives the value of eigenvector that corresponds to mode 1 at node 13 in dof 1
puts "EV111=$eigenvector111";
puts "EV511=$eigenvector511";
puts "EV911=$eigenvector911";
puts "EV1311=$eigenvector1311";
############################################################################
# Pushover Analysis #
############################################################################
if {$analysisType == "pushover"} {
puts "Running Pushover..."
# assign lateral loads and create load pattern:
set lat2 [expr 1*$eigenvector511]; # Load factor on each frame node in Floor 2
set lat3 [expr 1*$eigenvector911]; # Load factor on each frame node in Floor 3
set lat4 [expr 1*$eigenvector1311]; # Load factor on each frame node in Floor 4
pattern Plain 200 Linear {
load 5 $lat2 0.0 0.0;
load 6 $lat2 0.0 0.0;
load 7 $lat2 0.0 0.0;
load 8 $lat2 0.0 0.0;
load 9 $lat3 0.0 0.0;
load 10 $lat3 0.0 0.0;
load 11 $lat3 0.0 0.0;
load 12 $lat3 0.0 0.0;
load 13 $lat4 0.0 0.0;
load 14 $lat4 0.0 0.0;
load 15 $lat4 0.0 0.0;
load 16 $lat4 0.0 0.0;
}
# display deformed shape:
set ViewScale 5;
DisplayModel2D DeformedShape $ViewScale ; # display deformed shape, the scaling factor needs to be adjusted for each model
# displacement parameters
set IDctrlNode 13; # node where disp is read for disp control
set IDctrlDOF 1; # degree of freedom read for disp control (1 = x displacement)
set Dmax [expr 0.1*$HBuilding]; # maximum displacement of pushover: 10% roof drift
set Dincr [expr 0.1]; # displacement increment
# analysis commands
constraints Plain; # how it handles boundary conditions
numberer RCM; # renumber dof's to minimize band-width (optimization)
system BandSPD; # how to store and solve the system of equations in the analysis (large model: try UmfPack)
test NormDispIncr 1.0e-6 6; # type of convergence criteria with tolerance, max iterations
algorithm Newton; # use Newton's solution algorithm: updates tangent stiffness at every iteration
integrator LoadControl 1000000; # use Load-controlled analysis (=1kN)
analysis Static; # define type of analysis: static for pushover
set Nsteps 1500;# number of pushover analysis steps
set ok [analyze $Nsteps]; # this will return zero if no convergence problems were encountered
puts "Pushover complete"; # display this message in the command window
#integrator DisplacementControl $IDctrlNode $IDctrlDOF $Dincr; # use displacement-controlled analysis
#set Nsteps [expr int($Dmax/$Dincr)];# number of pushover analysis steps
}
############################################################################
# Time History/Dynamic Analysis #
############################################################################
if {$analysisType == "dynamic"} {
puts "Running dynamic analysis..."
# display deformed shape:
set ViewScale 500; # amplify(=to render larger) display of deformed shape
DisplayModel2D DeformedShape $ViewScale; # display deformed shape, the scaling factor needs to be adjusted for each model
# Rayleigh Damping
# calculate damping parameters
set zeta 0.02; # percentage of critical damping
set a1 [expr $zeta*2.0/($w1 + $w2)];
rayleigh 0.0 0.0 $a1 0.0; # assign stiffness proportional damping to frame beams & columns
# define ground motion parameters
set patternID 1; # load pattern ID
set GMdirection 1; # ground motion direction (1 = x)
set GMfile "GMwave.acc"; # ground motion filename
set dt 0.02; # timestep of input GM file
set Scalefact 0.001; # ground motion scaling factor
set TotalNumberOfSteps 2048; # number of steps in ground motion
set GMtime [expr $dt*$TotalNumberOfSteps + 10.0]; # total time of ground motion + 10 sec of free vibration
# define the acceleration series for the ground motion
# syntax: "Series -dt $timestep_of_record -filePath $filename_with_acc_history -factor $scale_record_by_this_amount
set accelSeries "Series -dt $dt -filePath $GMfile -factor [expr $Scalefact*$g]";
# create load pattern: apply acceleration to all fixed nodes with UniformExcitation
# command: pattern UniformExcitation $patternID $GMdir -accel $timeSeriesID
pattern UniformExcitation $patternID $GMdirection -accel $accelSeries;
# define dynamic analysis parameters
set dt_analysis 0.02; # timestep of analysis
wipeAnalysis; # destroy all components of the Analysis object, i.e.(=that is) any objects created with system, numberer, constraints, integrator, algorithm, and analysis commands
constraints Plain; # how it handles boundary conditions
numberer RCM; # renumber dof's to minimize band-width (optimization)
system UmfPack; # how to store and solve the system of equations in the analysis
test NormDispIncr 1.0e-8 50; # type of convergence criteria with tolerance, max iterations
algorithm NewtonLineSearch; # use NewtonLineSearch solution algorithm: updates tangent stiffness at every iteration and introduces line search to the Newton-Raphson algorithm to solve the nonlinear residual equation. Line search increases the effectiveness of the Newton method
integrator Newmark 0.5 0.25; # uses Newmark's average acceleration method to compute the time history
analysis Transient; # type of analysis: transient or static
set NumSteps [expr round(($GMtime + 0.0)/$dt_analysis)]; # number of steps in analysis
# perform the dynamic analysis and display whether analysis was successful
set ok [analyze $NumSteps $dt_analysis]; # ok = 0 if analysis was completed
if {$ok == 0} {
puts "Dynamic analysis complete";
} else {
puts "Dynamic analysis did not converge";
}
# output time at end of analysis
set currentTime [getTime]; # get current analysis time (after dynamic analysis)
puts "The current time is: $currentTime";
wipe all;
}
wipe all;
Please check my scripts.
Moderators: silvia, selimgunay, Moderators