Dear all,
I have a double curvature RC column. It has rotation spring on the bottom, and shear spring on the top. However it can not run due to a error "Matrix singular U(i,i)=0. Do anyone can help me to solve this problem.
Thank you so much for your help.
#-------------------------
# create the model builder
#-------------------------
wipe; # clear memory of all past model definitions
file mkdir Data_T1-1_spring; # create data directory
model Basic -ndm 2 -ndf 3
#---------------
# create nodes
#---------------
node 1 0.0 0.0
node 2 0.0 31.5
node 3 0.0 63.0
node 4 0.0 94.5
node 5 0.0 126.0
node 11 0.0 0.0
node 55 0.0 126.0
#-----------------------------
# Constrain: fix node 1
#-----------------------------
fix 1 1 1 1
fix 5 0 0 0
equalDOF 1 11 1 2 3
#---------------------------------------
# Define materials for nonlinear columns
# --------------------------------------
set ConConfTag 1;
set ConUnconTag 2;
set LongBarTag 3;
set TransBarTag 4;
set RigidMatTag 5;
set ShearCurveTag 7;
set ShearTag 8;
set BcTag 9;
set ZeroLengthElementTag 10;
# CONCRETE tag f'c ec0 f'cu ecu
# Core concrete (confined)
uniaxialMaterial Concrete01 $ConConfTag -8.546 -0.004 -5.0 -0.014
# Cover concrete (unconfined)
uniaxialMaterial Concrete01 $ConUnconTag -6.309 -0.002 0.0 -0.006
# Reinforcing steel
set fy 68.603; # Yield stress of long bars
set fyt 67.733; # Yield stress of trans. bars
set E 29733.0; # Young's modulus
# tag fy E0 b
uniaxialMaterial Steel01 $LongBarTag $fy $E 0.01
# SPRING
uniaxialMaterial Elastic $RigidMatTag 999e9;
#---------------------------------
# Define cross-section for columns
# --------------------------------
# set some paramaters
set colWidth 31.5;
set colDepth 31.5;
set cover 2.2;
set As 1.23; # area of no. 10 bars
# 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 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 16 4 [expr -$y1] [expr $z1-$cover] $y1 $z1
patch rect 2 16 4 [expr -$y1] [expr -$z1] $y1 [expr $cover-$z1]
patch rect 2 16 4 [expr -$y1] [expr $cover-$z1] [expr $cover-$y1] [expr $z1-$cover]
patch rect 2 16 4 [expr $y1-$cover] [expr $cover-$z1] $y1 [expr $z1-$cover]
# Create the reinforcing fibers (left, middle, right)
layer straight 3 7 $As [expr $y1-$cover] [expr $z1-$cover] [expr $y1-$cover] [expr -$z1+$cover]
layer straight 3 2 $As [expr $y1-$cover-4.5] [expr $z1-$cover] [expr $y1-$cover-4.5] [expr -$z1+$cover]
layer straight 3 2 $As [expr $y1-$cover-4.5*2] [expr $z1-$cover] [expr $y1-$cover-4.5*2] [expr -$z1+$cover]
layer straight 3 2 $As 0.0 [expr $z1-$cover] 0.0 [expr -$z1+$cover]
layer straight 3 2 $As [expr -$y1+$cover+4.5*2] [expr $z1-$cover] [expr -$y1+$cover+4.5*2] [expr -$z1+$cover]
layer straight 3 2 $As [expr -$y1+$cover+4.5] [expr $z1-$cover] [expr -$y1+$cover+4.5] [expr -$z1+$cover]
layer straight 3 7 $As [expr -$y1+$cover] [expr $z1-$cover] [expr -$y1+$cover] [expr -$z1+$cover]
}
#-----------------------
# Define column elements
# ----------------------
# Geometry of column elements
# tag
geomTransf PDelta 1
# Number of integration points along length of element
set np 5
# Create the coulumns using Beam-column elements
# e tag ndI ndJ nsecs secID transfTag
set eleType dispBeamColumn
element $eleType 1 11 2 $np 1 1
element $eleType 2 2 3 $np 1 1
element $eleType 3 3 4 $np 1 1
element $eleType 4 4 55 $np 1 1
# zero-length sections for bar-slip
element zeroLengthSection 5 1 11 1 -orient 0 1 0 -1 0 0
# --------- Shear model material properties and zeroLength Element --------- #
#limitCurve RotationShearCurve $crvTag $eleTag $ndI $ndJ $rotAxis $Vn $Vr $Kdeg $rotLim
set ndI 1
set ndJ 2
set rotAxis 3
set Vn 71.94
set Vr 0.450
set Kdeg -6807
set rotLim 0.02606
limitCurve RotationShearCurve $ShearCurveTag 1 $ndI $ndJ $rotAxis $Vn $Vr $Kdeg $rotLim
#uniaxialMaterial PinchingLimitStateMaterial $matTag $nodeT $nodeB $driftAxis $Kelas $crvTyp $crvTag $eleTag $b $d $h $a $st $As $Acc $ld $db $rhot $f'c $fy $fyt
#set matTag $shearMatTag
set nodeT 11
set nodeB 55
set driftAxis 1
set Kelas -1
set crvTyp 2
set b 31.5
set d 31.5
set h 39.37
set a 39.37
set st 5.9
set As [expr 17.18]
set Acc [expr 992]
set ld 100.0 ;
set db 1.25
set rhot 0.0030
set fc 6.309
set fy 6.86
set fyt 6.77
# Calibrated Model
uniaxialMaterial PinchingLimitStateMaterial 20 $nodeT $nodeB $driftAxis $Kelas $crvTyp $ShearCurveTag 1 $b $d $h $a $st $As $Acc $ld $db $rhot $fc $fy $fyt
# zero-length length for shear response
element zeroLength 6 55 5 -mat 20 $RigidMatTag $RigidMatTag -dir 1 2 6
puts "shear spring ok"
#---------------------
# Recorder
# --------------------
set HSec [expr 31.5]
set BSec [expr 31.5]
set coverY [expr $HSec/2.0];
set coverZ [expr $BSec/2.0];
set coreY [expr $coverY-$cover];
set coreZ [expr $coverZ-$cover];
set numIntgrPts 5; # number of integration points for force-based element
recorder Node -file Data_T1-1_spring/DFree.out -time -node 5 -dof 1 2 disp; # displacements of free nodes
recorder Node -file Data_T1-1_spring/DBase.out -time -node 1 -dof 1 2 disp; # displacements of support nodes
recorder Node -file Data_T1-1_spring/RBase.out -time -node 1 -dof 1 2 reaction; # support reaction
recorder Drift -file Data_T1-1_spring/Drift.out -time -iNode 1 -jNode 2 -dof 1 -perpDirn 2 ; # lateral drift
recorder Element -file Data_T1-1_spring/FCol.out -time -ele 1 4 globalForce; # element forces -- column
recorder Element -file Data_T1-1_spring/ForceColSec1.out -time -ele 1 4 section 1 force; # Column section forces, axial and moment, node i
recorder Element -file Data_T1-1_spring/DefoColSec1.out -time -ele 1 4 section 1 deformation; # section deformations, axial and curvature, node i
recorder Element -file Data_T1-1_spring/ForceColSec$numIntgrPts.out -time -ele 1 4 section $numIntgrPts force; # section forces, axial and moment, node j
recorder Element -file Data_T1-1_spring/DefoColSec$numIntgrPts.out -time -ele 1 4 section 1 deformation; # section deformations, axial and curvature, node j
#recorder Element -xml stressstrain.out -time -ele $bcTag section $sectionFiberTag fiber 0 0 $IDconcCore stressStrain;
#recorder Element -xml stressstraincover.out -time -ele $bcTag section $sectionFiberTag fiber $coverY $coverZ $IDconcCore stressStrain;
#recorder Element -xml stressstrainreinf.out -time -ele $bcTag section $sectionFiberTag fiber $coreY 0 $IDreinf2 stressStrain;
#recorder Element -xml stressstrainreinf2.out -time -ele $bcTag section $sectionFiberTag fiber $coreY $coreZ $IDreinf stressStrain;
# record and plot node displacement
#recorder Node -file Data_T1-1_spring/topdisp.out -time -node 5 -dof 1 disp;
#recorder plot Data_T1-1_spring/topdisp.out Node5_Xdisp 10 10 300 300 -columns 2 1;
#---------------------
# Define Axial load
# --------------------
set P [expr -3332.0 ]
pattern Plain 1 Linear {
load 5 0.0 $P 0.0
}
# Gravity-analysis parameters -- load-controlled static analysis
set Tol 1.0e-8; # convergence tolerance for test
constraints Plain; # how it handles boundary conditions
numberer Plain; # renumber dof's to minimize band-width (optimization), if you want to
system BandGeneral; # how to store and solve the system of equations in the analysis
test NormDispIncr $Tol 6 ; # determine if convergence has been achieved at the end of an iteration step
algorithm Newton; # use Newton's solution algorithm: updates tangent stiffness at every iteration
set NstepGravity 1; # 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"
pattern Plain 2 Linear {
load 5 1.0 0.0 0.0
}
set IDctrlNode 5
set IDctrlDOF 1
set Tol 1.0e-3; # convergence tolerance for test
constraints Plain; # how it handles boundary conditions
numberer Plain; # renumber dof's to minimize band-width (optimization), if you want to
system BandGeneral; # how to store and solve the system of equations in the analysis
test NormDispIncr $Tol 20; # determine if convergence has been achieved at the end of an iteration step
algorithm Newton; # use Newton's solution algorithm: updates tangent stiffness at every iteration
#Cyclic displacement
#Drift ratio= 0.25%
set Dincr [expr 0.01]
set Dmax [expr 0.315]
set Nsteps [expr int(abs($Dmax/$Dincr))]; # number of pushover analysis steps
integrator DisplacementControl $IDctrlNode $IDctrlDOF $Dincr
analysis Static
analyze $Nsteps
set Dincr_2 [expr -0.01]
set Dmax_2 [expr -0.315]
set Nsteps_2 [expr 2*int(abs($Dmax_2/$Dincr_2))]; # number of pushover analysis steps
integrator DisplacementControl $IDctrlNode $IDctrlDOF $Dincr_2
analysis Static
analyze $Nsteps_2
puts "cyclic lateral displacement finished!"
# ------------------------------
# Finally perform the analysis
# ------------------------------
# display displacement shape of the column
recorder display "Displaced shape" 10 10 500 500 -file 1.bmp
prp 189. 189. 1;
vup 0 1 0;
vpn 0 0 1;
display 1 5 40
print -bmp -file model.bmp
recorder plot Data_T1-1/DFree.out 10 10 500 500 1 2
puts "Done"
The model did not work due to singular matrix
Moderators: silvia, selimgunay, Moderators