Can anyone help, i build this very simple model 2 story 1 bay with elasticBeamColumn and rotational spring to simulate the plastic hinge. I check and recheck numerous time. it still keep Pop up error of " Warning BandGenLinLapckSolver::Solve<> -LAPACK routine returned 4,7,7,7,7,7,7" and give me incorrect mode period.
Here is the code.
you can find the source of "rotSpring2DModIKModel.tcl" in "http://opensees.berkeley.edu/wiki/index ... KModel.tcl."
Thank you so much for checking!!
#---------------------------------------------------------------
# Dynamic Analysis
# Units: kips, inches, seconds
# Created by: Dong Wen Liang
###################################################################################################
# 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
source rotSpring2DModIKModel.tcl; # procedure for defining a rotational spring (zero-length element) for plastic hinges
# procedure for defining a rotational spring (zero-length element) to capture panel zone shear distortions
# procedure for defining 8 elements to create a rectangular panel zone
set xPixels 1024;
set yPixels 786;
set xLoc1 30;
set yLoc1 30;
set dAmp 2;
#DisplayModel2D NodeNumbers $dAmp $xLoc1 $yLoc1 $xPixels $yPixels or DeformedShape
###################################################################################################
# Define Analysis Type
###################################################################################################
# Define type of analysis: "pushover" = pushover; "dynamic" = dynamic
set analysisType "dynamic";
if {$analysisType == "pushover"} {
set dataDir Concentrated-PanelZone-Pushover-Output; # name of output folder
file mkdir $dataDir; # create output folder
}
if {$analysisType == "dynamic"} {
set dataDir Concentrated-PanelZone-Dynamic-Output; # name of output folder
file mkdir $dataDir; # create output folder
}
###################################################################################################
# Define Building Geometry, Nodes, Masses, and Constraints
###################################################################################################
# define structure-geometry parameters
set NStories 2; # number of stories
set NBays 1; # number of frame bays (excludes bay for P-delta column)
set WBay [expr 30.0*12.0]; # bay width in inches
set HStory1 [expr 13.0*12.0]; # 1st story height in inches
set HStory2 [expr 13.0*12.0 + $HStory1]; # 2nd story height in inches
set HStory3 [expr 13.0*12.0 + $HStory2]; # 3rd story height in inches
set HStoryTyp [expr 13.0*12.0]; # story height of other stories in inches
set HBuilding [expr ($NStories)* $HStory1]; # height of building
# calculate locations of beam-column joint centerlines:
set Pier1 0.0; # leftmost column line
set Pier2 [expr $Pier1 + $WBay];
set Pier3 [expr $Pier2 + $WBay];
set Floor1 0.0; # ground floor
set Floor2 [expr $Floor1 + $HStory1];
set Floor3 [expr $Floor2 + $Floor2];
# calculate panel zone dimensions
set Pzlat13 [expr 0]; # lateral dist from CL of beam-col joint to edge of panel zone (= half the column depth) of pier 1 floor 3
set Pzvert13 [expr 0]; # vert dist from CL of beam-col joint to edge of panel zone (= half the beam depth) of pier 1 floor 3
set Pzlat12 [expr 0]; # lateral dist from CL of beam-col joint to edge of panel zone (= half the column depth) of pier 1 floor 2
set Pzvert12 [expr 0]; # vert dist from CL of beam-col joint to edge of panel zone (= half the beam depth) of pier 1 floor 2
# calculate plastic hinge offsets from beam-column centerlines:
set Phlat12 [expr $Pzlat12]; # lateral dist from CL of beam-col joint to beam hinge
set Phvert12 [expr $Pzvert12]; # vert dist from CL of beam-col joint to column hinge (forms at edge of panel zone)
set Phlat13 [expr $Pzlat13];
set Phvert13 [expr $Pzvert13];
#------------------------------------------------------------------------------------
set Pzlat23 [expr 0];
set Pzvert23 [expr 0];
set Pzlat22 [expr 0];
set Pzvert22 [expr 0];
set Phlat22 [expr $Pzlat22]; # lateral dist from CL of beam-col joint to beam hinge
set Phvert22 [expr $Pzvert22]; # vert dist from CL of beam-col joint to column hinge (forms at edge of panel zone)
set Phlat23 [expr $Pzlat23];
set Phvert23 [expr $Pzvert23];
# define nodes and assign masses to beam-column intersections of frame
# command: node nodeID xcoord ycoord -mass mass_dof1 mass_dof2 mass_dof3
# nodeID convention: "xy" where x = Pier # and y = Floor #
node 11 $Pier1 $Floor1;
node 12 $Pier1 $Floor2;
node 13 $Pier1 $Floor3;
node 21 $Pier2 $Floor1;
node 22 $Pier2 $Floor2;
node 23 $Pier2 $Floor3;
node 31 $Pier3 $Floor1;
# define extra nodes for plastic hinge rotational springs
# nodeID convention: "xya" where x = Pier #, y = Floor #, a = location relative to beam-column joint
# "a" convention: 1,2 = left; 3,4 = right; (used for beams)
# "a" convention: 5,6 = below; 7,8 = above; (used for columns)
# column hinges at bottom of Story 1 (base)
node 111 $Pier1 $Floor1;
node 211 $Pier2 $Floor1;
node 311 $Pier3 $Floor1;
#------------------------------------------------------------------------------
# column & Beam hinges at top of (Pier 1,2nd floor)
node 121 $Pier1 [expr $Floor2 - $Phvert12];
node 123 $Pier1 [expr $Floor2 + $Phvert12];
node 125 [expr $Pier1 + $Pzlat12] $Floor2;
# column & Beam hinges at top of (Pier 1,3rd floor)
node 133 $Pier1 [expr $Floor3 + $Phvert13];
node 135 [expr $Pier1 + $Pzlat13] $Floor3;
# column & Beam hinges at top of (Pier 2,2nd floor)
node 221 $Pier2 [expr $Floor2 - $Phvert22];
node 225 [expr $Pier2 + $Pzlat22] $Floor2;
node 227 [expr $Pier2 - $Pzlat22] $Floor2;
# column & Beam hinges at top of (Pier 2,3rd floor)
node 235 [expr $Pier2 + $Pzlat23] $Floor3;
node 237 [expr $Pier2 - $Pzlat23] $Floor3;
# column & Beam hinges at top of (Pier 5,2th floor)
node 32 $Pier3 $Floor2;
node 321 $Pier3 $Floor2;# zero-stiffness spring will be used on p-delta column
node 322 $Pier3 $Floor2;# zero-stiffness spring will be used on p-delta column
# column & Beam hinges at top of (Pier 5,3th floor)
node 33 $Pier3 $Floor3;
node 331 $Pier3 $Floor3;
# calculate nodal masses -- lump floor masses at frame nodes
set g 386.2; # acceleration due to gravity
set WBuilding 6503.637; # total building weight kip
set NodalMass2 [expr 5.463 / 8.0]; # mass at each node on Floor 2 (2109.824 kip / 386.2 in/s^2)=5.463
set NodalMass3 [expr 5.463 / 8.0]; # mass at each node on Floor 3 (2109.824 kip / 386.2 in/s^2)=5.463
set NodalMass4 [expr 5.937 / 8.0]; # mass at each node on Floor 4
set Negligible 1e-9; # a very small number to avoid problems with zero
# define nodal masses: lump at beam-column joints in frame
# command: mass $nodeID $dof1mass $dof2mass $dof3mass
mass 12 $NodalMass2 $Negligible $Negligible; # Pier 1, Floor 2
mass 13 $NodalMass3 $Negligible $Negligible; # Pier 1, Floor 3
mass 22 $NodalMass2 $Negligible $Negligible; # Pier 2, Floor 2
mass 23 $NodalMass3 $Negligible $Negligible; # Pier 2, Floor 3
# 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 12 22 $dof1; # Floor 2: Pier 1 to Pier 2
equalDOF 12 32 $dof1; # Floor 3: Pier 1 to Pier 2
equalDOF 13 23 $dof1; # Floor 2: Pier 1 to Pier 3
equalDOF 13 33 $dof1; # Floor 2: Pier 1 to Pier 3
# 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 11 1 1 1;
fix 21 1 1 1;
fix 31 1 1 0;
###################################################################################################
# Define Section Properties and Elements #
###################################################################################################
# define material properties
set Es 29000.0; # steel Young's modulus
set Fy 50.0; # steel yield strength
set mod 1.1; # modified moment of inertia value
#Define columns properties
# W14x257
set Acol_11 75.6;
set Icol_11 [expr 3400.0 * $mod];
set Mycol_11 20750.0;
set dcol_11 16.4;
set bfcol_11 16.0;
set tfcol_11 1.89;
set twcol_11 1.18;
# W14x311
set Acol_12 91.4;
set Icol_12 [expr 4330.0 * $mod];
set Mycol_12 25300.0;
set dcol_12 17.1;
set bfcol_12 16.2;
set tfcol_12 2.26;
set twcol_12 1.41;
#Define Beams properties
# w33x118
set Abeam_22 34.7;
set Ibeam_22 [expr 5900.0 * $mod];
set Mybeam_22 17950.0;
set dbeam_22 32.9;
set bfbeam_22 11.5;
set tfbeam_22 0.74;
set twbeam_22 0.55;
# w30x116
set Abeam_23 34.2;
set Ibeam_23 [expr 4930.0 * $mod];
set Mybeam_23 16450.0;
set dbeam_23 30.0;
set bfbeam_23 10.5;
set tfbeam_23 0.85;
set twbeam_23 0.565;
# determine stiffness modifications to equate the stiffness of the spring-elastic element-spring subassembly to the stiffness of the actual frame member
# References: (1) Ibarra, L. F., and Krawinkler, H. (2005). "Global collapse of frame structures under seismic excitations," Technical Report 152,
# The John A. Blume Earthquake Engineering Research Center, Department of Civil Engineering, Stanford University, Stanford, CA.
# (2) Zareian, F. and Medina, R. A. (2010). “A practical method for proper modeling of structural damping in inelastic plane
# structural systems,” Computers & Structures, Vol. 88, 1-2, pp. 45-53.
# calculate modified section properties to account for spring stiffness being in series with the elastic element stiffness
set n 10.0; # stiffness multiplier for rotational spring
set Ks_col_F11 [expr $n * 6.0 * $Es * $Icol_11 /($HStory1-$Phvert12)]; #spring-elastic element for column
set Ks_col_F12 [expr $n * 6.0 * $Es * $Icol_11 /($HStory1-$Phvert12-$Phvert13)];
set Ks_col_F21 [expr $n * 6.0 * $Es * $Icol_12 /($HStory1-$Phvert12)];
set Ks_col_F22 [expr $n * 6.0 * $Es * $Icol_12 /($HStory1-$Phvert12-$Phvert13)];
set Ks_Beam_12 [expr $n * 6.0 * $Es * $Ibeam_22 /($WBay-$Phlat12-$Phlat22)]; #spring-elastic element for beam
set Ks_Beam_13 [expr $n * 6.0 * $Es * $Ibeam_23 /($WBay-$Phlat13-$Phlat23)];
# set up geometric transformation of elements
set IDColTransf 1; #all columns
#set IDBeamTransf 2; #all beams
#geomTransf Linear $IDBeamTransf; # Linear transformation 2
geomTransf PDelta $IDColTransf; # PDelta transformation 1
# 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 #
# pier 1 column
element elasticBeamColumn 111 111 125 75.6 $Es $Icol_11 $IDColTransf;
element elasticBeamColumn 112 121 135 75.6 $Es $Icol_11 $IDColTransf;
# pier 2 column
element elasticBeamColumn 211 211 225 91.4 $Es $Icol_12 $IDColTransf;
element elasticBeamColumn 212 221 235 91.4 $Es $Icol_12 $IDColTransf;
# define elastic beam elements
# element between plastic hinges: eleID convention = "2xy" where 2 = beam, x = Bay #, y = Floor #
# element between plastic hinge and panel zone: eleID convention = "2xya" where 2 = beam, x = Bay #, y = Floor #, a = loc in bay
# Floor 1 beam
element elasticBeamColumn 12 123 227 34.7 $Es $Ibeam_22 $IDColTransf;
# Floor 2 beam
element elasticBeamColumn 13 133 237 34.2 $Es $Ibeam_23 $IDColTransf;
# define p-delta columns and rigid links
set TrussMatID 600; # define a material ID
set Arigid 1000.0; # define area of truss section (make much larger than A of frame elements)
set Irigid 100000.0; # moment of inertia for p-delta columns (make much larger than I of frame elements)
uniaxialMaterial Elastic $TrussMatID $Es; # define truss material
# rigid links
# command: element truss $eleID $iNode $jNode $A $materialID
# eleID convention: 6xy, 6 = truss link, x = Bay #, y = Floor #
element truss 522 22 32 $Arigid $TrussMatID; # Floor 2
element truss 523 23 33 $Arigid $TrussMatID; # Floor 3
# p-delta columns
# eleID convention: 7xy, 7 = p-delta columns, x = Pier #, y = Story #
element elasticBeamColumn 511 31 321 $Arigid $Es $Irigid $IDColTransf; # Story 1
element elasticBeamColumn 512 322 331 $Arigid $Es $Irigid $IDColTransf; # Story 2
###################################################################################################
# Define Rotational Springs for Plastic Hinges, Panel Zones, and Leaning Columns
###################################################################################################
# define rotational spring properties and create spring elements using "rotSpring2DModIKModel" procedure
# rotSpring2DModIKModel creates a uniaxial material spring with a bilinear response based on Modified Ibarra Krawinkler Deterioration Model
# references provided in rotSpring2DModIKModel.tcl
# input values for Story 1 column springs
set McMy 1.05; # ratio of capping moment to yield moment, Mc / My
set LS 1000.0; # basic strength deterioration (a very large # = no cyclic deterioration)
set LK 1000.0; # unloading stiffness deterioration (a very large # = no cyclic deterioration)
set LA 1000.0; # accelerated reloading stiffness deterioration (a very large # = no cyclic deterioration)
set LD 1000.0; # post-capping strength deterioration (a very large # = no deterioration)
set cS 1.0; # exponent for basic strength deterioration (c = 1.0 for no deterioration)
set cK 1.0; # exponent for unloading stiffness deterioration (c = 1.0 for no deterioration)
set cA 1.0; # exponent for accelerated reloading stiffness deterioration (c = 1.0 for no deterioration)
set cD 1.0; # exponent for post-capping strength deterioration (c = 1.0 for no deterioration)
set th_pP 0.025; # plastic rot capacity for pos loading
set th_pN 0.025; # plastic rot capacity for neg loading
set th_pcP 0.3; # post-capping rot capacity for pos loading
set th_pcN 0.3; # post-capping rot capacity for neg loading
set ResP 0.4; # residual strength ratio for pos loading
set ResN 0.4; # residual strength ratio for neg loading
set th_uP 0.4; # ultimate rot capacity for pos loading
set th_uN 0.4; # ultimate rot capacity for neg loading
set DP 1.0; # rate of cyclic deterioration for pos loading
set DN 1.0; # rate of cyclic deterioration for neg loading
#Pier 1 floor 1 column spring
set a_mem [expr ($n+1.0)*($Mycol_11*($McMy-1.0)) / ($Ks_col_F11*$th_pP)]; # strain hardening ratio of spring
set b [expr ($a_mem)/(1.0+$n*(1.0-$a_mem))]; # modified strain hardening ratio of spring (Ibarra & Krawinkler 2005, note: Eqn B.5 is incorrect)
rotSpring2DModIKModel 91111 11 111 $Ks_col_F11 $b $b 20750.0 -20750.0 $LS $LK $LA $LD $cS $cK $cA $cD $th_pP $th_pN $th_pcP $th_pcN $ResP $ResN $th_uP $th_uN $DP $DN;
rotSpring2DModIKModel 91112 125 12 $Ks_col_F11 $b $b 20750.0 -20750.0 $LS $LK $LA $LD $cS $cK $cA $cD $th_pP $th_pN $th_pcP $th_pcN $ResP $ResN $th_uP $th_uN $DP $DN;
#Pier 1 floor 2 column spring
set a_mem [expr ($n+1.0)*($Mycol_11*($McMy-1.0)) / ($Ks_col_F12*$th_pP)];
set b [expr ($a_mem)/(1.0+$n*(1.0-$a_mem))];
rotSpring2DModIKModel 91121 12 121 $Ks_col_F12 $b $b 20750.0 -20750.0 $LS $LK $LA $LD $cS $cK $cA $cD $th_pP $th_pN $th_pcP $th_pcN $ResP $ResN $th_uP $th_uN $DP $DN;
rotSpring2DModIKModel 91122 135 13 $Ks_col_F12 $b $b 20750.0 -20750.0 $LS $LK $LA $LD $cS $cK $cA $cD $th_pP $th_pN $th_pcP $th_pcN $ResP $ResN $th_uP $th_uN $DP $DN;
#Pier 2 floor 2 column spring
set a_mem [expr ($n+1.0)*($Mycol_12*($McMy-1.0)) / ($Ks_col_F21*$th_pP)];
set b [expr ($a_mem)/(1.0+$n*(1.0-$a_mem))];
rotSpring2DModIKModel 92111 21 211 $Ks_col_F21 $b $b 25300.0 -25300.0 $LS $LK $LA $LD $cS $cK $cA $cD $th_pP $th_pN $th_pcP $th_pcN $ResP $ResN $th_uP $th_uN $DP $DN;
rotSpring2DModIKModel 92112 225 22 $Ks_col_F21 $b $b 25300.0 -25300.0 $LS $LK $LA $LD $cS $cK $cA $cD $th_pP $th_pN $th_pcP $th_pcN $ResP $ResN $th_uP $th_uN $DP $DN;
#Pier 2 floor 2 column spring
set a_mem [expr ($n+1.0)*($Mycol_12*($McMy-1.0)) / ($Ks_col_F22*$th_pP)];
set b [expr ($a_mem)/(1.0+$n*(1.0-$a_mem))];
rotSpring2DModIKModel 92131 22 221 $Ks_col_F22 $b $b 25300.0 -25300.0 $LS $LK $LA $LD $cS $cK $cA $cD $th_pP $th_pN $th_pcP $th_pcN $ResP $ResN $th_uP $th_uN $DP $DN;
rotSpring2DModIKModel 92132 235 23 $Ks_col_F22 $b $b 25300.0 -25300.0 $LS $LK $LA $LD $cS $cK $cA $cD $th_pP $th_pN $th_pcP $th_pcN $ResP $ResN $th_uP $th_uN $DP $DN;
region 1 -ele 91111 91112 91121 91122 92111 92112 92131 92132;
set th_pP 0.02;
set th_pN 0.02;
set th_pcP 0.16;
set th_pcN 0.16;
#Bay 1 Floor 1
set a_mem [expr ($n+1.0)*($Mybeam_22*($McMy-1.0)) / ($Ks_Beam_12*$th_pP)];
set b [expr ($a_mem)/(1.0+$n*(1.0-$a_mem))];
rotSpring2DModIKModel 9111 12 123 $Ks_Beam_12 $b $b 17950 -17950 $LS $LK $LA $LD $cS $cK $cA $cD $th_pP $th_pN $th_pcP $th_pcN $ResP $ResN $th_uP $th_uN $DP $DN;
rotSpring2DModIKModel 9112 227 22 $Ks_Beam_12 $b $b 17950 -17950 $LS $LK $LA $LD $cS $cK $cA $cD $th_pP $th_pN $th_pcP $th_pcN $ResP $ResN $th_uP $th_uN $DP $DN;
#Bay 1 Floor 2
set a_mem [expr ($n+1.0)*($Mybeam_23*($McMy-1.0)) / ($Ks_Beam_13*$th_pP)];
set b [expr ($a_mem)/(1.0+$n*(1.0-$a_mem))];
rotSpring2DModIKModel 9121 13 133 $Ks_Beam_13 $b $b 16450 -16450 $LS $LK $LA $LD $cS $cK $cA $cD $th_pP $th_pN $th_pcP $th_pcN $ResP $ResN $th_uP $th_uN $DP $DN;
rotSpring2DModIKModel 9122 237 23 $Ks_Beam_13 $b $b 16450 -16450 $LS $LK $LA $LD $cS $cK $cA $cD $th_pP $th_pN $th_pcP $th_pcN $ResP $ResN $th_uP $th_uN $DP $DN;
region 2 -ele 9111 9112 9121 9122;
# define p-delta column spring: zero-stiffness elastic spring
#Spring ID: "5xya" where 5 = leaning column spring, x = Pier #, y = Story #, a = location in story
# "a" convention: 1 = bottom of story, 2 = top of story
# rotLeaningCol ElemID ndR ndC
rotLeaningCol 5521 321 32; # top of Story 1
rotLeaningCol 5522 32 322; # bottom of Story 2
rotLeaningCol 5531 331 33; # top of Story 2
# create region for P-Delta column springsa
region 3 -ele 5521 5522 5531;
############################################################################
# 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 $nEigenI-1]]; # 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
DisplayModel2D DeformedShape $dAmp $xLoc1 $yLoc1 $xPixels $yPixels
Help!!! Where are the errors?
Moderators: silvia, selimgunay, Moderators
Re: Help!!! Where are the errors?
This error indicates the instability of your model. Check the boundary conditions and also make sure that you do not have "hanging" nodes.