HELP/nonlinearbeamcolumnpush

[mpaksoy]

Dear Vesna,

I performed pushover analysis on a 50m column with nonlinearbeamcolumn element. To see what it differs I performed the analysis with 5 element and only with one element. When I push the 5 element model, at a top level displacement of 2.28m the analysis fails and opensees gives notices 'Bar fails'.
I am not sure is it a compression faiulere of the rebar?
In the 1 element model, the column can be easily pushed beyond the 2.28 level. Why is that?
When I checked the force displacement response and moment curvature response of the 2 models. I noticed they are different. The curvature level of the 5 element model at the end of analysis is 0.026 but the 1 element models' curvature is at a level of 0.011. Although it has been pushed to a higher displacement. I send you an e-mail with these files to vesna@berkeley with a topic name"VESNAFROMMILANO". Im confused. I will be glad if you have a look.

Thank you very much.
Regards,
Murathan

[vesna]

Murathan,

A message you get regarding bar failure is related to ReinforcingSteel material. The bar fails after buckling in tension.

When you have only one element you are not accounting for P-smal delta effect. Bur if you have 5 elements you will account for some portion of P-small delta. This is one reason for seeing more deformation if you model the column with 5 elements. Also, in case you have 5 elements and 5 integration points per element, the distance between the two bottom integration points becomes really small which will create curvature concentrations. The common way to avoid it is to define number of integration points such that the bottom two points encompass the plastic hinge region.

[wuhua.zeng]

Dear vesna :
>
>I am modeling a rectangle column confined with FRP,the result is "bar fail' and "failed to get compatible element forces" , I tried to change the algorithm\the element type \ the number of integration points\the uniaxialMaterial , but failed again. I spent much time on it .I don't know what I should do now. I will be appreciated if you could help me.
The following is my code.
# SET UP ----------------------------------------------------------------------------
wipe; # clear memory of all past model definitions
model BasicBuilder -ndm 2 -ndf 3; # 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 basic and system units
#source BuildRCrectSection.tcl; # procedure for definining RC fiber section
# define GEOMETRY -------------------------------------------------------------
set LCol [expr 1500*$mm]; # column length
set Weight [expr 571*$kN]; # 上部结构重量 weight
# define section geometry
set HCol [expr 400*$mm]; # Column Depth
set BCol $HCol;
# calculated parameters
set PCol [expr $Weight]; # nodal dead-load weight per column
set Mass [expr $PCol/$g]; # nodal mass
# calculated geometry parameters

# nodal coordinates:
node 1 0 0; # node#, X, Y
node 2 0 $LCol
node 3 0 50;
node 4 0 450;

# Single point constraints -- Boundary Conditions
fix 1 1 1 1; # node DX DY RZ
fix 2 0 0 0; # node DX DY RZ


# nodal masses:
mass 2 $Mass 0. 0.; # node#, Mx My Mz, Mass=Weight/g, neglect rotational inertia at nodes


# Define ELEMENTS & SECTIONS -------------------------------------------------------------
set ColSecTag 1; # assign a tag number to the column section
set ColSecTagbelow 2;


# MATERIAL parameters -------------------------------------------------------------------
set IDunconcU 1; # material ID tag -- unconfined cover concrete
set IDconc 2;
set IDreinf 3;
set IDconcinfrp 4; # material ID tag -- reinforcement
set IDfrp 5;
# nominal concrete compressive strength
set fc [expr -19.75*$MPa]; # CONCRETE Compressive Strength, ksi (+Tension, -Compression)
set Ec [expr 28000*$MPa]; # Concrete Elastic Modulus
set fc1 [expr -20.54*$MPa];
# unconfined concrete
set fc1U $fc; # UNCONFINED concrete (todeschini parabolic model), maximum stress
set eps1U -0.006; # strain at maximum strength of unconfined concrete
set fc2U [expr 0.2*$fc1U]; # crushing stress
set eps2U -0.033; # strain at crushing stress
set lambda 0.1; # ratio between unloading slope at $eps2 and initial slope $Ec
# tensile-strength properties
set ftU [expr -0.1*$fc1U]; # tensile strength +tension
set Ets [expr $ftU/0.002]; # tension softening stiffness
#------------
set fcc $fc1;
set ec0 -0.00624;
set fcu [expr 0.8*$fc1];
set ecu -0.032;
#set ft [expr -0.14*$fc1];
#set Et [expr $ft/0.002];
# infrp concrete
set fccinfrp [expr 1.2*$fc1];
set ec0infrp -0.0075;
set fcuinfrp [expr 0.8*$fccinfrp];
set ecuinfrp -0.032;
# -----------
set Fy [expr 300*$MPa]; # STEEL yield stress
set Es [expr 200000*$MPa]; # modulus of steel
# set Bs 0.000001; # strain-hardening ratio
# set R0 18.5; # control the transition from elastic to plastic branches
# set cR1 0.925; # control the transition from elastic to plastic branches
# set cR2 0.15; # control the transition from elastic to plastic branches
set Ef 259500;
set ef 0.0152;
set fu [expr 450*$MPa];
set esh 0.1;
set eult 0.5;
set Esh [expr ($fu-$Fy)/($eult-$esh)*$MPa];
#reinforcngSteel para~
#set GABuck
set lsr 93.75; #lsr= 钢筋长度/钢筋直径
set beta 1.0;
set r 0.4;
set gama 0.5;
#set DMbuck
#set lsr;
#set alpha;
#set CMFatigue
set Cf 0.26;
set alpha 0.506;
set Cd 0.389;

#set IsoHard
set a1 4.3; # Hardening constant (default = 4.3)
set limit 1.0; #Limit for the reduction of the yield plateau.
#% of original plateau length to remain (0.01 < limit < 1.0 )
#Limit =1.0, then no reduction takes place (default =0.01)


uniaxialMaterial Concrete02 $IDunconcU $fc1U $eps1U $fc2U $eps2U $lambda $ftU $Ets; # build cover concrete (unconfined)
# uniaxialMaterial ConfinedConcrete01 $IDconcinfrp S1 $fccinfrp $Ec -epscu $ecuinfrp -varub 350.0 6.0 250.0 210.0 $Es 0.0 1000.0 16.0 -wrap 35.0 200.4 400 3944 $Ef -stRatio 1.0;
#uniaxialMaterial ConfinedConcrete01 $IDconc S1 $fcc $Ec -epscu $ecu -varub 350.0 6.0 250.0 210.0 $Es 0.0 1000.0 16.0 -stRatio 1.0;
uniaxialMaterial Concrete01 $IDconc $fcc $ec0 $fcu $ecu;
uniaxialMaterial Concrete01 $IDconcinfrp $fccinfrp $ec0infrp $fcuinfrp $ecuinfrp;
uniaxialMaterial ElasticPP $IDfrp $Ef $ef 0;
uniaxialMaterial ReinforcingSteel $IDreinf $Fy $fu $Es $Esh $esh $eult -CMFatigue $Cf $alpha $Cd; #-GABuck $lsr $beta $r $gama -IsoHard $a1 $limit;
# FIBER SECTION properties -------------------------------------------------------------
# Formal arguments
# id - tag for the section that is generated by this procedure
# HSec - depth of section, along local-y axis
# BSec - width of section, along local-z axis
# cH - distance from section boundary to neutral axis of reinforcement
# cB - distance from section boundary to side of reinforcement
# coreID - material tag for the core patch
# coverID - material tag for the cover patches
# steelID - material tag for the reinforcing steel
# numBarsTop - number of reinforcing bars in the top layer
# numBarsBot - number of reinforcing bars in the bottom layer
# numBarsIntTot - TOTAL number of reinforcing bars on the intermediate layers, symmetric about z axis and 2 bars per layer-- needs to be an even integer
# barAreaTop - cross-sectional area of each reinforcing bar in top layer
# barAreaBot - cross-sectional area of each reinforcing bar in bottom layer
# barAreaInt - cross-sectional area of each reinforcing bar in intermediate layer
# nfCoreY - number of fibers in the core patch in the y direction
# nfCoreZ - number of fibers in the core patch in the z direction
# nfCoverY - number of fibers in the cover patches with long sides in the y direction
# nfCoverZ - number of fibers in the cover patches with long sides in the z direction
# The core concrete ends at the NA of the reinforcement
# The center of the section is at (0,0) in the local axis system
#
set cover [expr 35*$mm]; # rectangular-RC-Column cover
#set numBarsTopCol 2; # number of longitudinal-reinforcement bars on top layer
#set numBarsBotCol 2; # number of longitudinal-reinforcement bars on bottom layer
# set numBarsIntCol 0; # TOTAL number of reinforcing bars on the intermediate layers
set barAreaTop [expr 201.1*$mm2]; # longitudinal-reinforcement bar area
set barAreaBot [expr 201.1*$mm2]; # longitudinal-reinforcement bar area
#set barAreaIntCol [expr 113*$mm2]; # longitudinal-reinforcement bar area

set coverY [expr $HCol/2.0]; # The distance from the section z-axis to the edge of the cover concrete -- outer edge of cover concrete
set coverZ [expr $BCol/2.0]; # The distance from the section y-axis to the edge of the cover concrete -- outer edge of cover concrete
set coreY [expr $coverY-$cover]; # The distance from the section z-axis to the edge of the core concrete -- edge of the core concrete/inner edge of cover concrete
set coreZ [expr $coverZ-$cover]; # The distance from the section y-axis to the edge of the core concrete -- edge of the core concrete/inner edge of cover concrete
set coreYfrp [expr $coverY];
set coreZfrp [expr $coverZ];
set coverYfrp [expr $HCol/2.0+0.501];
set coverZfrp [expr $HCol/2.0+0.501]
#set numBarsInt [expr $numBarsIntTot/2]; # number of intermediate bars per side

# Define the fiber section
section Fiber $ColSecTag {; # Define the fiber section

# Define the core patch
patch quad $IDconc 30 30 -$coreY $coreZ -$coreY -$coreZ $coreY -$coreZ $coreY $coreZ

# Define the four cover patches
patch quad $IDunconcU 30 2 -$coverY $coverZ -$coreY $coreZ $coreY $coreZ $coverY $coverZ
patch quad $IDunconcU 30 2 -$coreY -$coreZ -$coverY -$coverZ $coverY -$coverZ $coreY -$coreZ
patch quad $IDunconcU 2 30 -$coverY $coverZ -$coverY -$coverZ -$coreY -$coreZ -$coreY $coreZ
patch quad $IDunconcU 2 30 $coreY $coreZ $coreY -$coreZ $coverY -$coverZ $coverY $coverZ

# define reinforcing layers
layer straight $IDreinf 3 $barAreaTop $coreY $coreZ $coreY -$coreZ; # top layer reinfocement
layer straight $IDreinf 3 $barAreaBot -$coreY $coreZ -$coreY -$coreZ; # bottom layer reinforcement
layer straight $IDreinf 3 $barAreaTop -$coreY $coreZ $coreY $coreZ;
layer straight $IDreinf 3 $barAreaTop -$coreY -$coreZ $coreY -$coreZ;
}; # end of fibersection definition

section Fiber $ColSecTagbelow {; # Define the fiber section

# Define the core patch
patch quad $IDconcinfrp 30 30 -$coreYfrp $coreZfrp -$coreYfrp -$coreZfrp $coreYfrp -$coreZfrp $coreYfrp $coreZfrp;

# Define the frp cover patches
patch quad $IDfrp 30 1 -$coverYfrp $coverZfrp -$coreYfrp $coreZfrp $coreYfrp $coreZfrp $coverYfrp $coverZfrp
patch quad $IDfrp 30 1 -$coreYfrp -$coreZfrp -$coverYfrp -$coverZfrp $coverYfrp -$coverZfrp $coreYfrp -$coreZfrp
patch quad $IDfrp 1 30 -$coverYfrp $coverZfrp -$coverYfrp -$coverZfrp -$coreYfrp -$coreZfrp -$coreYfrp $coreZfrp
patch quad $IDfrp 1 30 $coreYfrp $coreZfrp $coreYfrp -$coreZfrp $coverYfrp -$coverZfrp $coverYfrp $coverZfrp


# define reinforcing layers
#layer straight $IDreinf 2 $barAreaInt -$coreY $coreZ $coreY $coreZ; # intermediate skin reinf. +z
#layer straight $IDreinf $numBarsInt $barAreaInt -$coreY -$coreZ $coreY -$coreZ; # intermediate skin reinf. -z
layer straight $IDreinf 3 $barAreaTop $coreY $coreZ $coreY -$coreZ; # top layer reinfocement
layer straight $IDreinf 3 $barAreaBot -$coreY $coreZ -$coreY -$coreZ; # bottom layer reinforcement
layer straight $IDreinf 3 $barAreaTop -$coreY $coreZ $coreY $coreZ;
layer straight $IDreinf 3 $barAreaTop -$coreY -$coreZ $coreY -$coreZ;
}; # end of fibersection definition


# define geometric transformation: performs a linear geometric transformation of beam stiffness and resisting force from the basic system to the global-coordinate system
set ColTransfTag 1; # associate a tag to column transformation
set ColTransfType Linear ; # options, Linear PDelta Corotational
geomTransf $ColTransfType $ColTransfTag ; # only columns can have PDelta effects (gravity effects)


# element connectivity:
set numIntgrPts 4; # number of integration points for force-based element
element nonlinearBeamColumn 1 1 3 $numIntgrPts $ColSecTag $ColTransfTag; # self-explanatory when using variables
element nonlinearBeamColumn 2 4 2 $numIntgrPts $ColSecTag $ColTransfTag;
element nonlinearBeamColumn 3 3 4 $numIntgrPts $ColSecTagbelow $ColTransfTag;
# Define RECORDERS -------------------------------------------------------------
set xunconfc -200;
set yunconfc 0;
set xc 0;
set yc 0;
set xsteel -165;
set ysteel -165;
recorder Node -file $dataDir/DFree.out -time -node 2 -dof 1 disp; # displacements of free nodes

recorder Node -file $dataDir/Dreaction.out -time -node 1 -dof 1 reaction;

recorder Element -file $dataDir/FCol.out -time -ele 1 globalForce; # element forces -- column

recorder Element -file $dataDir/unconfinedconcrete.out -time -ele 1 section 1 fiber $xunconfc $yunconfc stressStrain;

recorder Element -file $dataDir/confinedconcrete.out -time -ele 1 section 1fiber $xc $yc stressStrain;

recorder Element -file $dataDir/steel.out -time -ele 1 section 1 fiber $xsteel $ysteel stressStrain;

pattern Plain 1 Linear {
load 2 0 -571000 0
}

# define GRAVITY -------------------------------------------------------------

# 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 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"








# The static cycle.tcl is:
set IDctrlNode 2; # node where displacement is read for displacement control
set IDctrlDOF 1; # degree of freedom of displacement read for displacement contro
# characteristics of cyclic analysis
set iDmax "2 4 6 8 10 14 18 22 26 30 34 38 42 46 50 54 58 62 66 70 74 78 82 "; # vector of displacement-cycle peaks, in terms of storey drift ratio
set Dincr 0.5; # displacement increment for pushover. you want this to be very small, but not too small to slow down the analysis
#set Fact LCol; # scale drift ratio by storey height for displacement cycles
set CycleType Full; # you can do Full / Push / Half cycles with the proc
set Ncycles 1; # specify the number of cycles at each peak

# create load pattern for lateral pushover load
source LibUnits.tcl;
#set Weight [expr 741.7*$kN];
#set Hload 741700; # define the lateral load as a proportion of the weight so that the pseudo time equals the lateral-load coefficient when using linear load pattern
#set iPushNode 2 ; # define nodes where lateral load is applied in static lateral analysis

pattern Plain 2 Linear {
load 2 571000 0.0 0.0
}


# ----------- set up analysis parameters
source LibAnalysisStaticParameters.tcl; # constraintsHandler,DOFnumberer,system-ofequations,convergenceTest,solutionAlgorithm,integrator

# --------------------------------- perform Static Cyclic Displacements Analysis
source LibGeneratePeaks.tcl
set fmt1 "%s Cyclic analysis: CtrlNode %.3i, dof %.1i, Disp=%.4f %s"; # format for screen/file output of DONE/PROBLEM analysis
foreach Dmax $iDmax {
set iDstep [GeneratePeaks $Dmax $Dincr $CycleType ]; # this proc is defined above
for {set i 1} {$i <= $Ncycles} {incr i 1} {
set zeroD 0
set D0 0.0
foreach Dstep $iDstep {
set D1 $Dstep
set Dincr [expr $D1 - $D0]
integrator DisplacementControl $IDctrlNode $IDctrlDOF $Dincr
analysis Static
# ----------------------------------------------first analyze command------------------------
set ok [analyze 1]
# ----------------------------------------------if convergence failure-------------------------
if {$ok != 0} {
# if analysis fails, we try some other stuff
# performance is slower inside this loop global maxNumIterStatic; # max no. of iterations performed before "failure to converge" is ret'd
if {$ok != 0} {
puts "Trying Newton with Initial Tangent .."
test NormDispIncr $Tol 2000 0
algorithm Newton -initial
set ok [analyze 1]
test $testTypeStatic $TolStatic $maxNumIterStatic 0
algorithm $algorithmTypeStatic
}
if {$ok != 0} {
puts "Trying Broyden .."
algorithm Broyden 8
set ok [analyze 1 ]
algorithm $algorithmTypeStatic
}
if {$ok != 0} {
puts "Trying NewtonWithLineSearch .."
algorithm NewtonLineSearch 0.8
set ok [analyze 1]
algorithm $algorithmTypeStatic
}
if {$ok != 0} {
set putout [format $fmt1 "PROBLEM" $IDctrlNode $IDctrlDOF [nodeDisp $IDctrlNode $IDctrlDOF] $LunitTXT]
puts $putout
return -1
}; # end if
}; # end if
# -----------------------------------------------------------------------------------------------------
set D0 $D1; # move to next step
}; # end Dstep
}; # end i
}; # end of iDmaxCycl
# -----------------------------------------------------------------------------------------------------
if {$ok != 0 } {
puts [format $fmt1 "PROBLEM" $IDctrlNode $IDctrlDOF [nodeDisp $IDctrlNode $IDctrlDOF] $LunitTXT]
} else {
puts [format $fmt1 "DONE" $IDctrlNode $IDctrlDOF [nodeDisp $IDctrlNode $IDctrlDOF] $LunitTXT]
}









#The static parameters Tcl is :
set systemTypeStatic BandGeneral; # try UmfPack for large model
system $systemTypeStatic

# TEST: # convergence test to
# Convergence TEST (http://opensees.berkeley.edu/OpenSees/m ... al/360.htm)
# -- Accept the current state of the domain as being on the converged solution path
# -- determine if convergence has been achieved at the end of an iteration step
# NormUnbalance -- Specifies a tolerance on the norm of the unbalanced load at the current iteration
# NormDispIncr -- Specifies a tolerance on the norm of the displacement increments at the current iteration
# EnergyIncr-- Specifies a tolerance on the inner product of the unbalanced load and displacement increments at the current iteration
# RelativeNormUnbalance --
# RelativeNormDispIncr --
# RelativeEnergyIncr --
variable TolStatic 1.e-10; # Convergence Test: tolerance
variable maxNumIterStatic 6; # Convergence Test: maximum number of iterations that will be performed before "failure to converge" is returned
variable printFlagStatic 0; # Convergence Test: flag used to print information on convergence (optional) # 1: print information on each step;
variable testTypeStatic EnergyIncr ; # Convergence-test type
test $testTypeStatic $TolStatic $maxNumIterStatic $printFlagStatic;
# for improved-convergence procedure:
variable maxNumIterConvergeStatic 2000;
variable printFlagConvergeStatic 0;


# Solution ALGORITHM: -- Iterate from the last time step to the current (http://opensees.berkeley.edu/OpenSees/m ... al/682.htm)
# Linear -- Uses the solution at the first iteration and continues
# Newton -- Uses the tangent at the current iteration to iterate to convergence
# ModifiedNewton -- Uses the tangent at the first iteration to iterate to convergence
# NewtonLineSearch --
# KrylovNewton --
# BFGS --
# Broyden --
variable algorithmTypeStatic Newton
algorithm $algorithmTypeStatic;

# Static INTEGRATOR: -- determine the next time step for an analysis (http://opensees.berkeley.edu/OpenSees/m ... al/689.htm)
# LoadControl -- Specifies the incremental load factor to be applied to the loads in the domain
# DisplacementControl -- Specifies the incremental displacement at a specified DOF in the domain
# Minimum Unbalanced Displacement Norm -- Specifies the incremental load factor such that the residual displacement norm in minimized
# Arc Length -- Specifies the incremental arc-length of the load-displacement path
# Transient INTEGRATOR: -- determine the next time step for an analysis including inertial effects
# Newmark -- The two parameter time-stepping method developed by Newmark
# HHT -- The three parameter Hilbert-Hughes-Taylor time-stepping method
# Central Difference -- Approximates velocity and acceleration by centered finite differences of displacement
integrator DisplacementControl $IDctrlNode $IDctrlDOF $Dincr

# ANALYSIS -- defines what type of analysis is to be performed (http://opensees.berkeley.edu/OpenSees/m ... al/324.htm)
# Static Analysis -- solves the KU=R problem, without the mass or damping matrices.
# Transient Analysis -- solves the time-dependent analysis. The time step in this type of analysis is constant. The time step in the output is also constant.
# variableTransient Analysis -- performs the same analysis type as the Transient Analysis object. The time step, however, is variable. This method is used when
# there are convergence problems with the Transient Analysis object at a peak or when the time step is too small. The time step in the output is also variable.
set analysisTypeStatic Static
analysis $analysisTypeStatic


Best wishes !

[vesna]

Is your model working if you do eigen analysis?

[wuhua.zeng]

vesna wrote:
> Is your model working if you do eigen analysis?

Now , my model can run successfully, but the result seems wrong. The strength din't deteriorate smoothly, the load-displacement curve looks zigzag. Would you please guild me?

Here is my code:

wipe; # clear memory of all past model definitions
model BasicBuilder -ndm 2 -ndf 3; # 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 basic and system units
#source BuildRCrectSection.tcl; # procedure for definining RC fiber section
# define GEOMETRY -------------------------------------------------------------
set LCol [expr 1500*$mm]; # column length
set Weight [expr 571*$kN]; # 上部结构重量 weight
# define section geometry
set HCol [expr 400*$mm]; # Column Depth
set BCol $HCol;
# calculated parameters
set PCol [expr $Weight]; # nodal dead-load weight per column
set Mass [expr $PCol/$g]; # nodal mass
# calculated geometry parameters
set ACol [expr $HCol*$BCol]; # cross-sectional area
set IzCol [expr 1./12*$BCol*pow($HCol,3)]; # Column moment of inertia
# nodal coordinates:
node 1 0 0; # node#, X, Y
node 2 0 450;
node 3 0 800;
node 4 0 1150;
node 5 0 $LCol;
node 6 0 50;
node 7 0 150;
node 8 0 250;
node 9 0 350;
# Single point constraints -- Boundary Conditions
fix 1 1 1 1; # node DX DY RZ
fix 5 0 0 0; # node DX DY RZ


# nodal masses:
mass 5 $Mass 0. 0.; # node#, Mx My Mz, Mass=Weight/g, neglect rotational inertia at nodes


# Define ELEMENTS & SECTIONS -------------------------------------------------------------
set ColSecTag 11; # assign a tag number to the column section
set ColSecTagbelow 22;
#set BeamSecTag 2; # assign a tag number to the beam section
# define section geometry
#set coverCol [expr 15.*$mm]; # Column cover to reinforcing steel NA.
#set numBarsCol 4; # number of longitudinal-reinforcement bars in each side of column section. (symmetric top & bot)
#set barAreaCol [expr 113*$mm2]; # area of longitudinal-reinforcement bars

# MATERIAL parameters -------------------------------------------------------------------
set IDunconcU 1; # material ID tag -- unconfined cover concrete
set IDconc 2;
set IDreinf 3; # material ID tag -- reinforcement
set IDconcinfrp 4;
set IDfrp 5;
# nominal concrete compressive strength
set fc [expr -19.43*$MPa]; # CONCRETE Compressive Strength, ksi (+Tension, -Compression)
#set Ec [expr 24289*$MPa]; # Concrete Elastic Modulus
# set fc1 [expr -30.52*$MPa];
# unconfined concrete
set fc1U $fc; # UNCONFINED concrete (todeschini parabolic model), maximum stress
set eps1U -0.002; # strain at maximum strength of unconfined concrete
set fc2U [expr 0.2*$fc1U]; # crushing stress
set eps2U -0.0035; # strain at crushing stress
set lambda 0.1; # ratio between unloading slope at $eps2 and initial slope $Ec
# tensile-strength properties
set ftU [expr -0.1*$fc1U]; # tensile strength +tension
set Ets [expr $ftU/0.002]; # tension softening stiffness
#confined concrete
set fcc [expr -19.82*$MPa];
set ec0 -0.002;
set fcu [expr 0.2*$fcc];
set ecu -0.011;
set fccinfrp [expr -23*$MPa];
set ec0infrp -0.0024;
set fcuinfrp [expr 0.24*$fccinfrp];
set ecuinfrp -0.036;
#set ft [expr -0.14*$fc1];
#set Et [expr $ft/0.002];
# -----------
set Fy [expr 375*$MPa]; # STEEL yield stress
set Es [expr 201000*$MPa]; # modulus of steel
set Bs 0.008; # strain-hardening ratio
set R0 18; # control the transition from elastic to plastic branches
set cR1 0.925; # control the transition from elastic to plastic branches
set cR2 0.15; # control the transition from elastic to plastic branches

# set fu [expr 452*$MPa];
# set esh 0.1;
# set eult 0.49;
# set Esh [expr ($fu-$Fy)/($eult-$esh)*$MPa];
#reinforcngSteel para~
#set GABuck
# set lsr 214; #lsr= 钢筋长度/钢筋直径
# set beta 0.35;
# set r 1;
# set gama 1;
#set DMbuck
#set lsr;
#set alpha;
#set CMFatigue
# set Cf 0.5;
# set alpha 0.206;
# set Cd 0.39;

#set IsoHard
# set a1 3.3; # Hardening constant (default = 4.3)
# set limit 1.0; #Limit for the reduction of the yield plateau.
#% of original plateau length to remain (0.01 < limit < 1.0 )
#Limit =1.0, then no reduction takes place (default =0.01)


uniaxialMaterial Concrete02 $IDunconcU $fc1U $eps1U $fc2U $eps2U $lambda $ftU $Ets; # build cover concrete (unconfined)
uniaxialMaterial Concrete01 $IDconc $fcc $ec0 $fcu $ecu;
uniaxialMaterial Concrete01 $IDconcinfrp $fccinfrp $ec0infrp $fcuinfrp $ecuinfrp;
# uniaxialMaterial Concrete02 $IDconc $fcc $ec0 $fcu $ecu $lambda $ft $Et;
#uniaxialMaterial Steel02 $IDreinf $Fy $Es $Bs $R0 $cR1 $cR2; # build reinforcement material
uniaxialMaterial Steel02 $IDreinf $Fy $Es $Bs $R0 $cR1 $cR2; #-IsoHard $a1 $limit;
uniaxialMaterial ElasticPP $IDfrp 259500 0.004 0 0;
uniaxialMaterial MinMax 6 $IDfrp -min 0 -max 0.003;

# Define the fiber section
section Fiber $ColSecTag {; # Define the fiber section

# Define the core patch
patch rect $IDconc 30 30 -165 -165 165 165 ;

# Define the four cover patches
patch rect $IDunconcU 30 6 -200 -200 200 -165 ;
patch rect $IDunconcU 30 6 -200 165 200 200 ;
patch rect $IDunconcU 6 30 -200 -165 -165 165 ;
patch rect $IDunconcU 6 30 165 -165 200 165 ;

# define reinforcing layers
#layer straight $IDreinf 2 $barAreaInt -$coreY $coreZ $coreY $coreZ; # intermediate skin reinf. +z
#layer straight $IDreinf $numBarsInt $barAreaInt -$coreY -$coreZ $coreY -$coreZ; # intermediate skin reinf. -z
layer straight $IDreinf 3 $barAreaInt -165 -165 165 -165; # top layer reinfocement
layer straight $IDreinf 3 $barAreaInt 165 -165 165 165; # bottom layer reinforcement
layer straight $IDreinf 3 $barAreaInt -165 165 165 165;
layer straight $IDreinf 3 $barAreaInt -165 -165 -165 165;
}; # end of fibersection definition
section Fiber $ColSecTagbelow {; # Define the fiber section

# Define the core patch
patch rect $IDconcinfrp 30 30 -200 -200 200 200 ;

# Define the four cover patches
patch rect $IDfrp 30 6 -200.334 -200.334 200.334 -200 ;
patch rect $IDfrp 30 6 -200.334 200 200.334 200.334;
patch rect $IDfrp 6 30 -200.334 -200 -200 200 ;
patch rect $IDfrp 6 30 200 -200 200.334 200 ;

# define reinforcing layers
#layer straight $IDreinf 2 $barAreaInt -$coreY $coreZ $coreY $coreZ; # intermediate skin reinf. +z
#layer straight $IDreinf $numBarsInt $barAreaInt -$coreY -$coreZ $coreY -$coreZ; # intermediate skin reinf. -z
layer straight $IDreinf 3 $barAreaInt -165 -165 165 -165; # top layer reinfocement
layer straight $IDreinf 3 $barAreaInt 165 -165 165 165; # bottom layer reinforcement
layer straight $IDreinf 3 $barAreaInt -165 165 165 165;
layer straight $IDreinf 3 $barAreaInt -165 -165 -165 165;
}; # end of fibersection definition


# define geometric transformation: performs a linear geometric transformation of beam stiffness and resisting force from the basic system to the global-coordinate system
set ColTransfTag 1; # associate a tag to column transformation
set ColTransfType Linear ; # options, Linear PDelta Corotational
geomTransf $ColTransfType $ColTransfTag ; # only columns can have PDelta effects (gravity effects)


# element connectivity:
set numIntgrPts 5; # number of integration points for force-based element
element nonlinearBeamColumn 1 1 6 2 $ColSecTag $ColTransfTag; # self-explanatory when using variables
element nonlinearBeamColumn 2 2 3 $numIntgrPts $ColSecTag $ColTransfTag;
element nonlinearBeamColumn 3 3 4 $numIntgrPts $ColSecTag $ColTransfTag;
element nonlinearBeamColumn 4 4 5 $numIntgrPts $ColSecTag $ColTransfTag;
element nonlinearBeamColumn 5 6 7 3 $ColSecTagbelow $ColTransfTag;
element nonlinearBeamColumn 6 7 8 3 $ColSecTagbelow $ColTransfTag;
element nonlinearBeamColumn 7 8 9 3 $ColSecTagbelow $ColTransfTag;
element nonlinearBeamColumn 8 9 2 3 $ColSecTagbelow $ColTransfTag;
# Define RECORDERS -------------------------------------------------------------
set xunconfc -200;
set yunconfc 0;
set xc 0;
set yc 0;
set xsteel -165;
set ysteel -165;
recorder Node -file $dataDir/DFree.out -time -node 5 -dof 1 disp; # displacements of free nodes

recorder Node -file $dataDir/Dreaction.out -time -node 1 -dof 1 reaction;

recorder Element -file $dataDir/FCol.out -time -ele 1 localForce; # element forces -- column

recorder Element -file $dataDir/sectionforce.out -time -ele 1 section 1 force;

recorder Element -file $dataDir/sectiondeformation.out -time -ele 1 section 1 deformation;

recorder Element -file $dataDir/unconfinedconcrete.out -time -ele 1 section 1 fiber $xunconfc $yunconfc stressStrain;

recorder Element -file $dataDir/confinedconcrete.out -time -ele 1 section 1 fiber $xc $yc stressStrain;

recorder Element -file $dataDir/steel.out -time -ele 1 section 1 fiber $xsteel $ysteel stressStrain;

pattern Plain 1 Linear {
load 5 0 -$Weight 0
}

# define GRAVITY -------------------------------------------------------------

# 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 EnergyIncr $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 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"



# we need to set up parameters that are particular to the model.
set IDctrlNode 5; # node where displacement is read for displacement control
set IDctrlDOF 1; # degree of freedom of displacement read for displacement contro
# characteristics of cyclic analysis
set iDmax "2 4 6 8 10 12 14 16 18 20 22 24 26 28 30 32 34 36 38 40 42 44 46 48 50 52 54 56 58 "; # vector of displacement-cycle peaks, in terms of storey drift ratio
set Dincr 0.5; # displacement increment for pushover. you want this to be very small, but not too small to slow down the analysis
set Fact LCol; # scale drift ratio by storey height for displacement cycles
set CycleType Full; # you can do Full / Push / Half cycles with the proc
set Ncycles 1; # specify the number of cycles at each peak

# create load pattern for lateral pushover load
source LibUnits.tcl;
#set Weight [expr 741.7*$kN];
#set Hload 741700; # define the lateral load as a proportion of the weight so that the pseudo time equals the lateral-load coefficient when using linear load pattern
#set iPushNode 2 ; # define nodes where lateral load is applied in static lateral analysis

pattern Plain 2 Linear {
load 5 571000 0.0 0.0
}


# ----------- set up analysis parameters
source LibAnalysisStaticParameters.tcl; # constraintsHandler,DOFnumberer,system-ofequations,convergenceTest,solutionAlgorithm,integrator

# --------------------------------- perform Static Cyclic Displacements Analysis
source LibGeneratePeaks.tcl
puts "LibGeneratePeaks has been done"
set fmt1 "%s Cyclic analysis: CtrlNode %.3i, dof %.1i, Disp=%.4f %s"; # format for screen/file output of DONE/PROBLEM analysis
foreach Dmax $iDmax {
set iDstep [GeneratePeaks $Dmax $Dincr $CycleType ]; # this proc is defined above
for {set i 1} {$i <= $Ncycles} {incr i 1} {
set zeroD 0
set D0 0.0
foreach Dstep $iDstep {
set D1 $Dstep
set Dincr [expr $D1 - $D0]
integrator DisplacementControl $IDctrlNode $IDctrlDOF $Dincr
analysis Static
# ----------------------------------------------first analyze command------------------------
set ok [analyze 1]
# ----------------------------------------------if convergence failure-------------------------
if {$ok != 0} {
# if analysis fails, we try some other stuff
# performance is slower inside this loop global maxNumIterStatic; # max no. of iterations performed before "failure to converge" is ret'd
if {$ok != 0} {
puts "Trying Newton with Initial Tangent .."
test NormDispIncr $Tol 2000 0
algorithm Newton -initial
set ok [analyze 1]
test $testTypeStatic $TolStatic $maxNumIterStatic 0
algorithm $algorithmTypeStatic
}
if {$ok != 0} {
puts "Trying Broyden .."
algorithm Broyden 8
set ok [analyze 1 ]
algorithm $algorithmTypeStatic
}
if {$ok != 0} {
puts "Trying NewtonWithLineSearch .."
algorithm NewtonLineSearch 0.8
set ok [analyze 1]
algorithm $algorithmTypeStatic
}
if {$ok != 0} {
set putout [format $fmt1 "PROBLEM" $IDctrlNode $IDctrlDOF [nodeDisp $IDctrlNode $IDctrlDOF] $LunitTXT]
puts $putout
return -1
}; # end if
}; # end if
# -----------------------------------------------------------------------------------------------------
set D0 $D1; # move to next step
}; # end Dstep
}; # end i
}; # end of iDmaxCycl
# -----------------------------------------------------------------------------------------------------
if {$ok != 0 } {
puts [format $fmt1 "PROBLEM" $IDctrlNode $IDctrlDOF [nodeDisp $IDctrlNode $IDctrlDOF] $LunitTXT]
} else {
puts [format $fmt1 "DONE" $IDctrlNode $IDctrlDOF [nodeDisp $IDctrlNode $IDctrlDOF] $LunitTXT]
}

[wuhua.zeng]

vesna wrote:
> Is your model working if you do eigen analysis?


I have tried to change the constituti verelation of concrete(change "set fcu [expr 0.2*$fcc]" to "set fcu [expr 0.8*$fcc]") ,and the load-displacement curve looks smooth, but the strength don't drop after a few cycls, it seems unrealistic.
When fcu = 0.8*$fcc and BS =0.0001(strain-hardening ratio of Steel02) ,the load-displacement curve looks smooth, and the strength flatten out,but it still not dropped.what's the reason?

[vesna]

If you do not see degradation of load-displcement curve try to push your model to a larger displacement.

The zig-zaged degradation you see on a pushover curve is something I would expect. With the materials you use it cannot be smooth.

[wuhua.zeng]

Thank you for your reply,vesna
The largest displacement I can get is 88mm, but I still can't see the degradation (set "fcu=0.8fcc").Is it because of the FRP material ?I am not sure if I have defined the right material .
Thank you!

[andy0071]

I am afraid that the model you build of FRP is not correct.

because of you write the material of the FRP as

uniaxialMaterial ElasticPP $IDfrp $Ef $ef 0;

but the material of FRP is not as same as steel which could keep the stress after the elastic transformation into plasticity. the FRP material will break and quit after it reach the ultimate load, so the curve will come down after the elastic transformation but not go horizontal. That means you could not use the "uniaxialMaterial ElasticPP" to model the FRP, is that the reason for your wrong result?

I don't know. Now I am doing the same analysis of FRP reinforce strucutre and could not use the basic material and Parallel material to build the FRP!!!
╮(╯▽╰)╭, who can tell me how to do it??? ~~~~(>_<)~~~~

[vesna]

There is material ConfinedDoncrete01 that accounts for FRP in the section. See if it is going to work for you:http://opensees.berkeley.edu/wiki/index ... 1_Material

[wuhua.zeng]

vesna wrote:
> There is material ConfinedDoncrete01 that accounts for FRP in the section.
> See if it is going to work for
> you:http://opensees.berkeley.edu/wiki/index ... 1_Material


Thank you! I tried ConfinedDoncrete01, and I got the degradation. Thank you!

[wuhua.zeng]

Than you for your reply,andy0071 . I'm afraid you are right