I am facing a convergence issue with my 3D model. The following warnings and errors have been encountered:
WARNING BandGenLinLapackSolver::solve() -factorization failed, matrix singular U(i,i) = 0, i= 3
WARNING NewtonRaphson::solveCurrentStep() -the LinearSysOfEqn failed in solve()
StaticAnalysis::analyze() - the Algorithm failed at iteration: 0 with domain at load factor 0.1
OpenSees > analyze failed, returned: -3 error flag
I have made modifications to the model, which was originally based on the generic 3D frame provided by the OpenSees website. I am unsure about the primary cause of the issue in my model. If you could provide any insight into this matter, I would be very grateful. Thank you.
Code: Select all
# --------------------------------------------------------------------------------------------------;
# Example 8. 3D Steel W-section frame;
# Silvia Mazzoni & Frank McKenna, 2006;
#;
# SET UP ----------------------------------------------------------------------------;
wipe; # clear memory of all past model definitions;
model BasicBuilder -ndm 3 -ndf 6; # Define the model builder, ndm=#dimension, ndf=#dofs;
set dataDir Data; # set up name of data directory;
file mkdir $dataDir; # create data directory;
set GMdir "../GMfiles"; # ground-motion file directory;
source LibUnits.tcl; # define units;
source DisplayPlane.tcl; # procedure for displaying a plane in model;
source DisplayModel3D.tcl; # procedure for displaying 3D perspectives of model;
# ------ frame configuration;
set NStory 4; # number of stories above ground level;
set NBay 4; # number of bays in X direction;
set NBayZ 6; # number of bays in Z direction;
puts "Number of Stories in Y: $NStory, Number of bays in X: $NBay, Number of bays in Z: $NBayZ";
set NFrame [expr $NBayZ + 1]; # actually deal with frames in Z direction, as this is an easy extension of the 2d model;
# define GEOMETRY -------------------------------------------------------------;
# define structure-geometry paramters;
set LCol [expr 13*$ft]; # column height (parallel to Y axis);
set LBeam [expr 30*$ft]; # beam length (parallel to X axis);
set LGird [expr 30*$ft]; # girder length (parallel to Z axis);
# define NODAL COORDINATES;
set Dlevel 10000; # numbering increment for new-level nodes;
set Dframe 100; # numbering increment for new-frame nodes;
set Dzero3 990000000; # numbering increment for new zerolength nodes for 3;
set Dzero1 980000000; # numbering increment for new zerolength nodes for 1;
for {set frame 1} {$frame <=[expr $NFrame]} {incr frame 1} {;
set Z [expr ($frame-1)*$LGird];;
for {set level 1} {$level <=[expr $NStory+1]} {incr level 1} {;
set Y [expr ($level-1)*$LCol];;
for {set pier 1} {$pier <= [expr $NBay+1]} {incr pier 1} {;
set X [expr ($pier-1)*$LBeam];;
set nodeID [expr $frame*$Dframe+$level*$Dlevel+$pier];
node $nodeID $X $Y $Z; # actually define node;
if {$pier == 1} {;
set zeronodeID [expr $Dzero1 + $frame*$Dframe + $level*$Dlevel + $pier];
node $zeronodeID $X $Y $Z; # actually define node;
} elseif {$pier == [expr $NBay]} {;
if {$frame == 1 || $frame == [expr $NFrame]} {;
set zeronodeID [expr $Dzero3 + $frame*$Dframe + $level*$Dlevel + $pier];
node $zeronodeID $X $Y $Z; # actually define node;
} else {;
set zeronodeID1 [expr $Dzero1 + $frame*$Dframe + $level*$Dlevel + $pier];
node $zeronodeID1 $X $Y $Z; # actually define node;
set zeronodeID3 [expr $Dzero3 + $frame*$Dframe + $level*$Dlevel + $pier];
node $zeronodeID3 $X $Y $Z; # actually define node;
};
} elseif {$pier == [expr $NBay+1]} {;
if {$frame == 1 || $frame == [expr $NFrame]} {;
} else {;
set zeronodeID [expr $Dzero3 + $frame*$Dframe + $level*$Dlevel + $pier];
node $zeronodeID $X $Y $Z; # actually define node;
};
} else {;
set zeronodeID1 [expr $Dzero1 + $frame*$Dframe + $level*$Dlevel + $pier];
node $zeronodeID1 $X $Y $Z; # actually define node;
set zeronodeID3 [expr $Dzero3 + $frame*$Dframe + $level*$Dlevel + $pier];
node $zeronodeID3 $X $Y $Z; # actually define node;
};
};
};
};
# rigid diaphragm nodes;
set RigidDiaphragm ON ; # options: ON, OFF. specify this before the analysis parameters are set the constraints are handled differently.;
set Xa [expr ($NBay*$LBeam)/2]; # mid-span coordinate for rigid diaphragm;
set Za [expr ($NFrame-1)*$LGird/2];;
set iMasterNode "";
for {set level 2} {$level <=[expr $NStory+1]} {incr level 1} {;
set Y [expr ($level-1)*$LCol];;
# rigid-diaphragm nodes in center of each diaphram;
set MasterNodeID [expr 9900+$level];
node $MasterNodeID $Xa $Y $Za; # master nodes for rigid diaphragm;
fix $MasterNodeID 0 1 0 1 0 1; # constrain other dofs that dont belong to rigid diaphragm control;
lappend iMasterNode $MasterNodeID;
set perpDirn 2; # perpendicular to plane of rigid diaphragm;
for {set frame 1} {$frame <=[expr $NFrame]} {incr frame 1} {;
for {set pier 1} {$pier <= [expr $NBay+1]} {incr pier 1} {;
set nodeID [expr $level*$Dlevel+$frame*$Dframe+$pier];
rigidDiaphragm $perpDirn $MasterNodeID $nodeID; # define Rigid Diaphram,;
};
};
};
# determine support nodes where ground motions are input, for multiple-support excitation;
set iSupportNode "";
for {set frame 1} {$frame <=[expr $NFrame]} {incr frame 1} {;
set level 1;
for {set pier 1} {$pier <= [expr $NBay+1]} {incr pier 1} {;
set nodeID [expr $level*$Dlevel+$frame*$Dframe+$pier];
lappend iSupportNode $nodeID;
};
};
# BOUNDARY CONDITIONS;
fixY 0.0 1 1 1 0 1 0; # pin all Y=0.0 nodes;
# calculated MODEL PARAMETERS, particular to this model;
# Set up parameters that are particular to the model for displacement control;
set IDctrlNode [expr int((($NStory+1)*$Dlevel+$NFrame*$Dframe)+1)]; # node where displacement is read for displacement control;
set IDctrlDOF 1; # degree of freedom of displacement read for displacement control;
set LBuilding [expr $NStory*$LCol]; # total building height;
# Define SECTIONS -------------------------------------------------------------;
set n 10.0; # stiffness multiplier for rotational spring
set SectionType Elastic; # options: Elastic FiberSection ;
if {$SectionType == "Elastic"} {;
# material properties:;
set Es [expr 29000*$ksi]; # Steel Youngs Modulus;
set nu 0.3; # Poissons ratio;
set Gs [expr $Es/2./[expr 1+$nu]]; # Torsional stiffness Modulus;
set J $Ubig; # set large torsional stiffness;
# moment frame exterior; column sections: W14x257;
set Acol_M_Ex 75.6; # in^2 (Cross-sectional area)
set Icol_M_Ex 3400; # in^4 (Moment of inertia)
set Icol_z_M_Ex 1290; # in^4 (Moment of inertia)
set Jcol_M_Ex 79.1; # in^4 (Polar Moment of inertia)
set Ks_col_M_Ex [expr $n*6.0*$Es*$Icol_M_Ex/$LCol];# rotational stiffness of Story 1 column springs
# moment frame interior; column sections: W14x311;
set Acol_M_In 91.4; # in^2 (Cross-sectional area)
set Icol_M_In 4330; # in^4 (Moment of inertia)
set Icol_z_M_In 1610; # in^4 (Moment of inertia)
set Jcol_M_In 136; # in^4 (Polar Moment of inertia)
set Ks_col_M_In [expr $n*6.0*$Es*$Icol_M_In/$LCol];# rotational stiffness of Story 1 column springs
# gravity frame; column sections: W14x68;
set Acol_G 20; # in^2 (Cross-sectional area)
set Icol_G 722; # in^4 (Moment of inertia)
set Icol_z_G 121; # in^4 (Moment of inertia)
set Jcol_G 3.01; # in^4 (Polar Moment of inertia)
set Ks_G [expr $n*6.0*$Es*$Icol_G/$LCol]; # rotational stiffness of Story 1 column springs
# moment frame; beam sections on x-axis: W30x116;
set Abe_M 34.2; # in^2 (Cross-sectional area)
set Ibe_M 4930; # in^4 (Moment of inertia)
set Ibe_z_M 164; # in^4 (Moment of inertia)
set Jbe_M 6.43; # in^4 (Polar Moment of inertia)
set Ks_M [expr $n*6.0*$Es*$Ibe_M/$LGird]; # rotational stiffness of Floor 2 & 3 beam springs
# moment frame, top beam sections on x-axis: W24x62;
set Abe_M_T 18.2; # in^2 (Cross-sectional area)
set Ibe_M_T 1550; # in^4 (Moment of inertia)
set Ibe_z_M_T 34.5; # in^4 (Moment of inertia)
set Jbe_M_T 1.71; # in^4 (Polar Moment of inertia)
set Ks_M_T [expr $n*6.0*$Es*$Ibe_M/$LGird]; # rotational stiffness of Floor 2 & 3 beam springs
# gravity frame, beam sections on x-axis: W16x26;
set Abe_G 7.68; # in^2 (Cross-sectional area)
set Ibe_G 301; # in^4 (Moment of inertia)
set Ibe_z_G 44.2; # in^4 (Moment of inertia)
set Jbe_G 0.262; # in^4 (Polar Moment of inertia)
set Ks_G [expr $n*6.0*$Es*$Ibe_G/$LGird];# rotational stiffness of Floor 2 & 3 beam springs
set matIDhard 1; # material numbers for recorder (this stressstrain recorder will be blank, as this is an elastic section);
} else {;
puts "No section has been defined";
return -1;
};
set QdlCol_M_Ex [expr 257*$lbf/$ft]; # W-section weight per length;
set QdlCol_M_In [expr 311*$lbf/$ft]; # W-section weight per length;
set QdlCol_G [expr 68*$lbf/$ft]; # W-section weight per length;
set QBeam_M [expr 116*$lbf/$ft]; # W-section weight per length;
set QBeam_M_T [expr 62*$lbf/$ft]; # W-section weight per length;
set QBeam_G [expr 26*$lbf/$ft]; # W-section weight per length;
# --------------------------------------------------------------------------------------------------------------------------------;
# define ELEMENTS;
# set up geometric transformations of element;
# separate columns and beams, in case of P-Delta analysis for columns;
set IDColTransf 1; # all columns;
set IDBeamTransf 2; # all beams;
set IDGirdTransf 3; # all girds;
set ColTransfType PDelta ; # options for columns: Linear PDelta Corotational ;
geomTransf $ColTransfType $IDColTransf 0 0 1; # orientation of column stiffness affects bidirectional response.;
geomTransf Linear $IDBeamTransf 0 0 1;
geomTransf Linear $IDGirdTransf 1 0 0;
# columns;
set N0col [expr 10000-1]; # column element numbers;
set level 0;
for {set frame 1} {$frame <=[expr $NFrame]} {incr frame 1} {;
for {set level 1} {$level <=$NStory} {incr level 1} {;
for {set pier 1} {$pier <= [expr $NBay+1]} {incr pier 1} {;
set elemID [expr $N0col +$level*$Dlevel + $frame*$Dframe+$pier];
set nodeI [expr $level*$Dlevel + $frame*$Dframe+$pier];
set nodeJ [expr ($level+1)*$Dlevel + $frame*$Dframe+$pier];
if {$frame == 1 || $frame == $NFrame} {;
if {$pier == 1 || $pier == [expr $NBay+1]} {;
element elasticBeamColumn $elemID $nodeI $nodeJ $Acol_M_Ex $Es $Gs $Jcol_M_Ex $Icol_M_Ex $Icol_z_M_Ex $IDColTransf; # moment frame, exterior;
} else {;
element elasticBeamColumn $elemID $nodeI $nodeJ $Acol_M_In $Es $Gs $Jcol_M_In $Icol_M_In $Icol_z_M_In $IDColTransf; # moment frame, interior;
};
} else {;
element elasticBeamColumn $elemID $nodeI $nodeJ $Acol_G $Es $Gs $Jcol_G $Icol_G $Icol_z_G $IDColTransf; # gravity frame;
};
};
};
};
# girders -- parallel to X-axis;
source Pinch_grav.tcl;
source Bilin_moment.tcl;
uniaxialMaterial Elastic 999 $Es;
set N0gird 1000000; # beam element numbers;
for {set frame 1} {$frame <=[expr $NFrame]} {incr frame 1} {;
for {set level 2} {$level <=[expr $NStory+1]} {incr level 1} {;
for {set bay 1} {$bay <= $NBay} {incr bay 1} {;
set elemID [expr $N0gird +$level*$Dlevel + $frame*$Dframe+ $bay];
if {$level == [expr $NStory+1]} {;
if {$frame == 1 || $frame == $NFrame} {;
if {$bay == $NBay} {;
set nodeI [expr $level*$Dlevel + $frame*$Dframe+ $bay];
set nodeJ [expr $level*$Dlevel + $frame*$Dframe+ $bay+1];
element truss $elemID $nodeI $nodeJ 999 999; # truss;
} else {;
set zeronodeI [expr $Dzero1 + $level*$Dlevel + $frame*$Dframe+ $bay];
set zeronodeJ [expr $Dzero3 + $level*$Dlevel + $frame*$Dframe+ $bay+1];
element elasticBeamColumn $elemID $zeronodeI $zeronodeJ $Abe_M_T $Es $Gs $Jbe_M_T $Ibe_M_T $Ibe_z_M_T $IDBeamTransf; # girders for moment frame at top;
};
} else {;
set zeronodeI [expr $Dzero1 + $level*$Dlevel + $frame*$Dframe+ $bay];
set zeronodeJ [expr $Dzero3 + $level*$Dlevel + $frame*$Dframe+ $bay+1];
element elasticBeamColumn $elemID $zeronodeI $zeronodeJ $Abe_G $Es $Gs $Jbe_G $Ibe_G $Ibe_z_G $IDBeamTransf; # girders for gravity frame at top;
};
} else {;
if {$frame == 1 || $frame == $NFrame} {;
if {$bay == $NBay} {;
set nodeI [expr $level*$Dlevel + $frame*$Dframe+ $bay];
set nodeJ [expr $level*$Dlevel + $frame*$Dframe+ $bay+1];
element truss $elemID $nodeI $nodeJ 999 999; # truss;
} else {;
set zeronodeI [expr $Dzero1 + $level*$Dlevel + $frame*$Dframe+ $bay];
set zeronodeJ [expr $Dzero3 + $level*$Dlevel + $frame*$Dframe+ $bay+1];
element elasticBeamColumn $elemID $zeronodeI $zeronodeJ $Abe_M $Es $Gs $Jbe_M $Ibe_M $Ibe_z_M $IDBeamTransf; # girders for moment frame except top;
};
} else {;
set zeronodeI [expr $Dzero1 + $level*$Dlevel + $frame*$Dframe+ $bay];
set zeronodeJ [expr $Dzero3 + $level*$Dlevel + $frame*$Dframe+ $bay+1];
element elasticBeamColumn $elemID $zeronodeI $zeronodeJ $Abe_G $Es $Gs $Jbe_G $Ibe_G $Ibe_z_G $IDBeamTransf; # girders for gravity frame except top;
};
};
};
};
};
set N0zero 3000000; # zerolength element;
for {set frame 1} {$frame <=[expr $NFrame]} {incr frame 1} {;
for {set level 2} {$level <= [expr $NStory+1]} {incr level 1} {;
for {set pier 1} {$pier <= [expr $NBay+1]} {incr pier 1} {;
if {$level == [expr $NStory+1]} {;
if {$frame == 1 || $frame == $NFrame} {;
if {$pier == [expr $NBay+1]} {;
} elseif {$pier == 1} {;
set elemID [expr ($N0zero +$level*$Dlevel + $frame*$Dframe+ $pier)*10 + 1];
set nodeID [expr $frame*$Dframe+$level*$Dlevel+$pier];
set zeronodeID [expr $Dzero1 + $frame*$Dframe + $level*$Dlevel + $pier];
element zeroLength $elemID $nodeID $zeronodeID -mat 302 -dir 5;
equalDOF $nodeID $zeronodeID 1 2 3 4 6;
} elseif {$pier == $NBay} {;
set elemID [expr ($N0zero +$level*$Dlevel + $frame*$Dframe+ $pier)*10 + 3];
set nodeID [expr $frame*$Dframe+$level*$Dlevel+$pier];
set zeronodeID [expr $Dzero3 + $frame*$Dframe + $level*$Dlevel + $pier];
element zeroLength $elemID $zeronodeID $nodeID -mat 302 -dir 5;
equalDOF $zeronodeID $nodeID 1 2 3 4 6;
} else {;
set elemID1 [expr ($N0zero +$level*$Dlevel + $frame*$Dframe+ $pier)*10 + 1];
set elemID3 [expr ($N0zero +$level*$Dlevel + $frame*$Dframe+ $pier)*10 + 3];
set nodeID [expr $frame*$Dframe+$level*$Dlevel+$pier];
set zeronodeID1 [expr $Dzero1 + $frame*$Dframe + $level*$Dlevel + $pier];
set zeronodeID3 [expr $Dzero3 + $frame*$Dframe + $level*$Dlevel + $pier];
element zeroLength $elemID1 $nodeID $zeronodeID1 -mat 302 -dir 5;
equalDOF $nodeID $zeronodeID1 1 2 3 4 6;
element zeroLength $elemID3 $zeronodeID3 $nodeID -mat 302 -dir 5;
equalDOF $zeronodeID3 $nodeID 1 2 3 4 6;
};
} else {;
if {$pier == [expr $NBay+1]} {;
set elemID [expr ($N0zero +$level*$Dlevel + $frame*$Dframe+ $pier)*10 + 3];
set nodeID [expr $frame*$Dframe+$level*$Dlevel+$pier];
set zeronodeID [expr $Dzero3 + $frame*$Dframe + $level*$Dlevel + $pier];
element zeroLength $elemID $zeronodeID $nodeID -mat 401 -dir 5;
equalDOF $zeronodeID $nodeID 1 2 3 4 6;
} elseif {$pier == 1} {;
set elemID [expr ($N0zero +$level*$Dlevel + $frame*$Dframe+ $pier)*10 + 1];
set nodeID [expr $frame*$Dframe+$level*$Dlevel+$pier];
set zeronodeID [expr $Dzero1 + $frame*$Dframe + $level*$Dlevel + $pier];
element zeroLength $elemID $nodeID $zeronodeID -mat 401 -dir 5;
equalDOF $nodeID $zeronodeID 1 2 3 4 6;
} else {;
set elemID1 [expr ($N0zero +$level*$Dlevel + $frame*$Dframe+ $pier)*10 + 1];
set elemID3 [expr ($N0zero +$level*$Dlevel + $frame*$Dframe+ $pier)*10 + 3];
set nodeID [expr $frame*$Dframe+$level*$Dlevel+$pier];
set zeronodeID1 [expr $Dzero1 + $frame*$Dframe + $level*$Dlevel + $pier];
set zeronodeID3 [expr $Dzero3 + $frame*$Dframe + $level*$Dlevel + $pier];
element zeroLength $elemID1 $nodeID $zeronodeID1 -mat 401 -dir 5;
equalDOF $nodeID $zeronodeID1 1 2 3 4 6;
element zeroLength $elemID3 $zeronodeID3 $nodeID -mat 401 -dir 5;
equalDOF $zeronodeID3 $nodeID 1 2 3 4 6;
};
};
} else {;
if {$frame == 1 || $frame == $NFrame} {;
if {$pier == [expr $NBay+1]} {;
} elseif {$pier == 1} {;
set elemID [expr ($N0zero +$level*$Dlevel + $frame*$Dframe+ $pier)*10 + 1];
set nodeID [expr $frame*$Dframe+$level*$Dlevel+$pier];
set zeronodeID [expr $Dzero1 + $frame*$Dframe + $level*$Dlevel + $pier];
element zeroLength $elemID $nodeID $zeronodeID -mat 301 -dir 5;
equalDOF $nodeID $zeronodeID 1 2 3 4 6;
} elseif {$pier == $NBay} {;
set elemID [expr ($N0zero +$level*$Dlevel + $frame*$Dframe+ $pier)*10 + 3];
set nodeID [expr $frame*$Dframe+$level*$Dlevel+$pier];
set zeronodeID [expr $Dzero3 + $frame*$Dframe + $level*$Dlevel + $pier];
element zeroLength $elemID $zeronodeID $nodeID -mat 301 -dir 5;
equalDOF $zeronodeID $nodeID 1 2 3 4 6;
} else {;
set elemID1 [expr ($N0zero +$level*$Dlevel + $frame*$Dframe+ $pier)*10 + 1];
set elemID3 [expr ($N0zero +$level*$Dlevel + $frame*$Dframe+ $pier)*10 + 3];
set nodeID [expr $frame*$Dframe+$level*$Dlevel+$pier];
set zeronodeID1 [expr $Dzero1 + $frame*$Dframe + $level*$Dlevel + $pier];
set zeronodeID3 [expr $Dzero3 + $frame*$Dframe + $level*$Dlevel + $pier];
element zeroLength $elemID1 $nodeID $zeronodeID1 -mat 301 -dir 5;
equalDOF $nodeID $zeronodeID1 1 2 3 4 6;
element zeroLength $elemID3 $zeronodeID3 $nodeID -mat 301 -dir 5;
equalDOF $zeronodeID3 $nodeID 1 2 3 4 6;
};
} else {;
if {$pier == [expr $NBay+1]} {;
set elemID [expr ($N0zero +$level*$Dlevel + $frame*$Dframe+ $pier)*10 + 3];
set nodeID [expr $frame*$Dframe+$level*$Dlevel+$pier];
set zeronodeID [expr $Dzero3 + $frame*$Dframe + $level*$Dlevel + $pier];
element zeroLength $elemID $zeronodeID $nodeID -mat 401 -dir 5;
equalDOF $zeronodeID $nodeID 1 2 3 4 6;
} elseif {$pier == 1} {;
set elemID [expr ($N0zero +$level*$Dlevel + $frame*$Dframe+ $pier)*10 + 1];
set nodeID [expr $frame*$Dframe+$level*$Dlevel+$pier];
set zeronodeID [expr $Dzero1 + $frame*$Dframe + $level*$Dlevel + $pier];
element zeroLength $elemID $nodeID $zeronodeID -mat 401 -dir 5;
equalDOF $nodeID $zeronodeID 1 2 3 4 6;
} else {;
set elemID1 [expr ($N0zero +$level*$Dlevel + $frame*$Dframe+ $pier)*10 + 1];
set elemID3 [expr ($N0zero +$level*$Dlevel + $frame*$Dframe+ $pier)*10 + 3];
set nodeID [expr $frame*$Dframe+$level*$Dlevel+$pier];
set zeronodeID1 [expr $Dzero1 + $frame*$Dframe + $level*$Dlevel + $pier];
set zeronodeID3 [expr $Dzero3 + $frame*$Dframe + $level*$Dlevel + $pier];
element zeroLength $elemID1 $nodeID $zeronodeID1 -mat 401 -dir 5;
equalDOF $nodeID $zeronodeID1 1 2 3 4 6;
element zeroLength $elemID3 $zeronodeID3 $nodeID -mat 401 -dir 5;
equalDOF $zeronodeID3 $nodeID 1 2 3 4 6;
};
};
};
};
};
};
# beams -- parallel to Z-axis;
set N0beam 2000000; # gird element numbers;
for {set frame 1} {$frame <=[expr $NFrame-1]} {incr frame 1} {;
for {set level 2} {$level <=[expr $NStory+1]} {incr level 1} {;
for {set bay 1} {$bay <= $NBay+1} {incr bay 1} {;
set elemID [expr $N0beam + $level*$Dlevel +$frame*$Dframe+ $bay];
set nodeI [expr $level*$Dlevel + $frame*$Dframe+ $bay];
set nodeJ [expr $level*$Dlevel + ($frame+1)*$Dframe+ $bay];
element truss $elemID $nodeI $nodeJ 999 999; # truss;
};
};
};
# --------------------------------------------------------------------------------------------------------------------------------;
# Define GRAVITY LOADS, weight and masses;
# calculate dead load of frame, assume this to be an internal frame (do LL in a similar manner);
# calculate distributed weight along the beam length;
set GammaConcrete [expr 150*$pcf]; # Reinforced-Concrete floor slabs;
set Tslab [expr 6*$in]; # 6-inch slab;
set Lslab [expr 2*$LGird/2]; # slab extends a distance of $LGird/2 in/out of plane;
set DLfactor 2.0; # scale dead load up a little;
set Qslab [expr $GammaConcrete*$Tslab*$Lslab*$DLfactor]; ;
set QdlBeam_M [expr $Qslab + $QBeam_M]; # dead load distributed along beam (one-way slab);
set QdlBeam_T [expr $Qslab + $QBeam_M_T]; # dead load distributed along beam (one-way slab);
set QdlBeam_G [expr $Qslab + $QBeam_G]; # dead load distributed along beam (one-way slab);
set QdlGird $QBeam_M; # dead load distributed along girder;
set WeightCol_M_Ex [expr $QdlCol_M_Ex*$LCol]; # total Column weight;
set WeightCol_M_In [expr $QdlCol_M_In*$LCol]; # total Column weight;
set WeightCol_G [expr $QdlCol_G*$LCol]; # total Column weight;
set WeightBeam_M [expr $QdlBeam_M*$LBeam]; # total Beam weight;
set WeightBeam_T [expr $QdlBeam_T*$LBeam]; # total Beam weight;
set WeightBeam_G [expr $QdlBeam_G*$LBeam]; # total Beam weight;
# assign masses to the nodes that the columns are connected to ;
# each connection takes the mass of 1/2 of each element framing into it (mass=weight/$g);
set iFloorWeight "";
set WeightTotal 0.0;
set sumWiHi 0.0; # sum of storey weight times height, for lateral-load distribution;
for {set frame 1} {$frame <=[expr $NFrame]} {incr frame 1} {;
if {$frame == 1 || $frame == $NFrame} {;
set GirdWeightFact 1; # 1x1/2girder on exterior frames;
} else {;
set GirdWeightFact 2; # 2x1/2girder on interior frames;
};
for {set level 2} {$level <=[expr $NStory+1]} {incr level 1} { ; ;
set FloorWeight 0.0;
if {$level == [expr $NStory+1]} {;
set ColWeightFact 1; # one column in top story;
} else {;
set ColWeightFact 2; # two columns elsewhere;
};
for {set pier 1} {$pier <= [expr $NBay+1]} {incr pier 1} {;;
if {$pier == 1 || $pier == [expr $NBay+1]} {;
set BeamWeightFact 1; # one beam at exterior nodes;
if {$frame == 1 || $frame == $NFrame} {;
if {$level == [expr $NStory+1]} {;
set WeightNode [expr $ColWeightFact*$WeightCol_M_Ex/2 + $BeamWeightFact*$WeightBeam_T/2];
} else {;
set WeightNode [expr $ColWeightFact*$WeightCol_M_Ex/2 + $BeamWeightFact*$WeightBeam_M/2];
};
} else {;
set WeightNode [expr $ColWeightFact*$WeightCol_G/2 + $BeamWeightFact*$WeightBeam_G/2];
};
} else {;
set BeamWeightFact 2; # two beams elewhere;
if {$frame == 1 || $frame == $NFrame} {;
if {$level == [expr $NStory+1]} {;
set WeightNode [expr $ColWeightFact*$WeightCol_M_In/2 + $BeamWeightFact*$WeightBeam_T/2];
} else {;
set WeightNode [expr $ColWeightFact*$WeightCol_M_In/2 + $BeamWeightFact*$WeightBeam_M/2];
};
} else {;
set WeightNode [expr $ColWeightFact*$WeightCol_G/2 + $BeamWeightFact*$WeightBeam_G/2];
};
};
set MassNode [expr $WeightNode/$g];;
set nodeID [expr $level*$Dlevel+$frame*$Dframe+$pier];
mass $nodeID $MassNode 0. $MassNode 0. 0. 0.; # define mass;
set FloorWeight [expr $FloorWeight+$WeightNode];;
};
lappend iFloorWeight $FloorWeight;
set WeightTotal [expr $WeightTotal+ $FloorWeight];
set sumWiHi [expr $sumWiHi+$FloorWeight*($level-1)*$LCol]; # sum of storey weight times height, for lateral-load distribution;
};
};
set MassTotal [expr $WeightTotal/$g]; # total mass;
# --------------------------------------------------------------------------------------------------------------------------------;
# LATERAL-LOAD distribution for static pushover analysis;
# calculate distribution of lateral load based on mass/weight distributions along building height;
# Fj = WjHj/sum(WiHi) * Weight at each floor j;
set iFj "";
for {set level 2} {$level <=[expr $NStory+1]} {incr level 1} { ; ;
set FloorWeight [lindex $iFloorWeight [expr $level-1-1]];;
set FloorHeight [expr ($level-1)*$LCol];;
lappend iFj [expr $FloorWeight*$FloorHeight/$sumWiHi*$WeightTotal]; # per floor;
};
set iNodePush $iMasterNode; # nodes for pushover/cyclic, vectorized;
set iFPush $iFj; # lateral load for pushover, vectorized;
# Define RECORDERS -------------------------------------------------------------;
set FreeNodeID [expr $NFrame*$Dframe+($NStory+1)*$Dlevel+($NBay+1)]; # ID: free node;
set SupportNodeFirst [lindex $iSupportNode 0]; # ID: first support node;
set SupportNodeLast [lindex $iSupportNode [expr [llength $iSupportNode]-1]]; # ID: last support node;
set FirstColumn [expr $N0col + 1*$Dframe+1*$Dlevel +1]; # ID: first column;
recorder Node -file $dataDir/DFree.out -time -node $FreeNodeID -dof 1 2 3 disp; # displacements of free node;
recorder Node -file $dataDir/DBase.out -time -nodeRange $SupportNodeFirst $SupportNodeLast -dof 1 2 3 disp; # displacements of support nodes;
recorder Node -file $dataDir/RBase.out -time -nodeRange $SupportNodeFirst $SupportNodeLast -dof 1 2 3 reaction; # support reaction;
recorder Drift -file $dataDir/DrNode.out -time -iNode $SupportNodeFirst -jNode $FreeNodeID -dof 1 -perpDirn 2; # lateral drift;
recorder Element -file $dataDir/Fel1.out -time -ele $FirstColumn localForce; # element forces in local coordinates;
recorder Element -file $dataDir/ForceEle1sec1.out -time -ele $FirstColumn section 1 force; # section forces, axial and moment, node i;
recorder Element -file $dataDir/DefoEle1sec1.out -time -ele $FirstColumn section 1 deformation; # section deformations, axial and curvature, node i;
# Define DISPLAY -------------------------------------------------------------;
DisplayModel3D DeformedShape ; # options: DeformedShape NodeNumbers ModeShape;
# GRAVITY -------------------------------------------------------------;
# define GRAVITY load applied to beams and columns -- eleLoad applies loads in local coordinate axis;
pattern Plain 101 Linear {;
for {set frame 1} {$frame <=[expr $NFrame]} {incr frame 1} {;
for {set level 1} {$level <=$NStory} {incr level 1} {;
for {set pier 1} {$pier <= [expr $NBay+1]} {incr pier 1} {;
set elemID [expr $N0col + $level*$Dlevel +$frame*$Dframe+$pier];
if {$frame == 1 || $frame == $NFrame} {;
if {$pier == [expr $NBay+1]} {;
} else {;
eleLoad -ele $elemID -type -beamUniform 0. 0. -$QdlCol_M_Ex; # COLUMNS
};
} else {;
eleLoad -ele $elemID -type -beamUniform 0. 0. -$QdlCol_M_Ex; # COLUMNS
};
};
};
};
for {set frame 1} {$frame <=[expr $NFrame]} {incr frame 1} {;
for {set level 2} {$level <=[expr $NStory+1]} {incr level 1} {;
for {set bay 1} {$bay <= $NBay-1} {incr bay 1} {;
set elemID [expr $N0gird + $level*$Dlevel +$frame*$Dframe+ $bay];
eleLoad -ele $elemID -type -beamUniform -$QdlBeam_M 0.; # Girders;
};
};
};
for {set frame 1} {$frame <=[expr $NFrame-1]} {incr frame 1} {;
for {set level 2} {$level <=[expr $NStory+1]} {incr level 1} {;
for {set bay 1} {$bay <= $NBay+1} {incr bay 1} {;
set elemID [expr $N0beam + $level*$Dlevel +$frame*$Dframe+ $bay];
};
};
};
};
# apply GRAVITY-- # apply gravity load, set it constant and reset time to zero, load pattern has already been defined;
puts goGravity;
# Gravity-analysis parameters -- load-controlled static analysis;
set Tol 1.0e-8; # convergence tolerance for test;
variable constraintsTypeGravity Plain; # default;;
if { [info exists RigidDiaphragm] == 1} {;
if {$RigidDiaphragm=="ON"} {;
variable constraintsTypeGravity Lagrange; # large model: try Transformation;
}; # if rigid diaphragm is on;
}; # if rigid diaphragm exists;
constraints Penalty 1.0e14 1.0e14
numberer RCM; # renumber dofs to minimize band-width (optimization), if you want to;
system BandGeneral ; # how to store and solve the system of equations in the analysis (large model: try UmfPack);
test EnergyIncr $Tol 6 ; # determine if convergence has been achieved at the end of an iteration step;
algorithm Newton; # use Newtons solution algorithm: updates tangent stiffness at every iteration;
set NstepGravity 10; # apply gravity in 10 steps;
set DGravity [expr 1./$NstepGravity]; # first load increment;;
integrator LoadControl $DGravity; # determine the next time step for an analysis;
analysis Static; # define type of analysis static or transient;
analyze $NstepGravity; # apply gravity;
# ------------------------------------------------- maintain constant gravity loads and reset time to zero;
loadConst -time 0.0;
# -------------------------------------------------------------;
puts "Model Built";
All the best,