How to solve Convergence error of nonlinear static analysis

[omidgahremani]

Dear all
I am performing nonlinear static analysis on three steel structure with 3 ,6 and 12 stories that their plastic behavior modeled with concentrated plastic hinges just as done in the example represented for a 2 story steel frame in wiki opensees site.
for each structure at specified displacements convergence error happens and analysis cannot continue further and analysis fails.
for example for 3,6 and 12 structures the final displacements is 19 ,65 and 28 cm. I except the 12 story structure to have a higher displacement but converging error happening in 28 cm displacement and I don't know really what causes this error and How can I Overcome this problem.
I performed push over on the same structure using Sap2000 and the structure can have a displacement up to 90 cm without any error in analysis.
I've already checked the entire example represented for this process and I can't find any differences between creation of my models and the model given in the example and I don't know the cause of problem occured in my model.

I will be very grateful for your instructions and kind helps
best regards
gahremani

[fmk]

are you trying to do stuff in the script when the analysis fails. see analysis section of following basic example:

http://opensees.berkeley.edu/wiki/index ... r_Analysis

[omidgahremani]

Dear fmk
I've tried your method but I could not see any differences between old and new model's response ,
Here is the scripts I used before and the scripts that I used according your instruction ,
please contact me If you find anything that can help me with this:
____________________________
Old script :
____________________________
#all units are in meter-kg-Newton
constraints Plain
numberer RCM
system BandGeneral
test NormUnbalance 1.0e-6 400
algorithm Newton
integrator DisplacementControl 441 1 0.001
analysis Static
analyze 1000
__________________________
New script :
__________________________
#all units are in meter-kg-Newton
constraints Plain
numberer RCM
system BandGeneral
test NormDispIncr 1.0e-12 100
algorithm Newton
analysis Static
#new commands
set dU 0.001; # Displacement increment
# Change the integration scheme to be displacement control
# node dof init
integrator DisplacementControl 441 1 $dU
# Finally perform the analysis
# ------------------------------
# Set some parameters
set maxU 1; # Max displacement
set currentDisp 0.0;
set ok 0

while {$ok == 0 && $currentDisp < $maxU} {

set ok [analyze 1]

# if the analysis fails try initial tangent iteration
if {$ok != 0} {
puts "regular newton failed .. lets try an initail stiffness for this step"
test NormDispIncr 1.0e-12 1000
algorithm ModifiedNewton -initial
set ok [analyze 1]
if {$ok == 0} {puts "that worked .. back to regular newton"}
test NormDispIncr 1.0e-12 100
algorithm Newton
}
set currentDisp [nodeDisp 441 1]
}

if {$ok == 0} {
puts "Pushover analysis completed SUCCESSFULLY";
} else {
puts "Pushover analysis FAILED";
}
________________________________
P.s:
thank you very much for your time and attention
Please contact me if you find anything that can help me with this
I can send you the entire model scripts if it is needed.

[fmk]

is there enough iterations for initial to succeed .. it depends on your model .. put a 1 at end of test command:

if {$ok != 0} {
puts "regular newton failed .. lets try an initail stiffness for this step"
test NormDispIncr 1.0e-12 1000 1
algorithm ModifiedNewton -initial
set ok [analyze 1]
if {$ok == 0} {puts "that worked .. back to regular newton"}
test NormDispIncr 1.0e-12 100
algorithm Newton
}


and if that fails you can try other things like cutting step size, increasing tolerance for a step, different solution algorithms and integrators.

[omidgahremani]

Dear fmk
I don't know how to thank you for your help.
that worked .
what was that extra 1 at the end of test command for?
I used your command like this and it worked. however it took more time than further analysis:
#-----------------------------------------------------------------------------------------------------
#all units are in meter-kg-Newton
constraints Plain
numberer RCM
system BandGeneral
test NormDispIncr 1.0e-12 1000 1
algorithm Newton
analysis Static
#new commands
set dU 0.0001; # Displacement increment
# Change the integration scheme to be displacement control
# node dof init
integrator DisplacementControl 441 1 $dU
# Finally perform the analysis
# ------------------------------
# Set some parameters
set maxU 0.3; # Max displacement
set currentDisp 0.0;
set ok 0

while {$ok == 0 && $currentDisp < $maxU} {

set ok [analyze 1]

# if the analysis fails try initial tangent iteration
if {$ok != 0} {
puts "regular newton failed .. lets try an initail stiffness for this step"
test NormDispIncr 1.0e-12 10000 1
algorithm ModifiedNewton -initial
set ok [analyze 1]
if {$ok == 0} {puts "that worked .. back to regular newton"}
test NormDispIncr 1.0e-12 1000 1
algorithm Newton
}
set currentDisp [nodeDisp 441 1]
}

if {$ok == 0} {
puts "Pushover analysis completed SUCCESSFULLY";
} else {
puts "Pushover analysis FAILED";
}
#--------------------------------------------------------------------------------------------
Thank you a lot again.

[pkafando]

Dear Vesna,

I have tried to plot the base shear against the displacement but I am having some difficulties. I read in one of your post that I should plot the second column of Dfree.out and the Rbase.out but when I ran my file, my Rbase.out contains negative number, so I plotted without the negative number and compared it to a study done before on the same column but my results are slightly different. Please help me out. I need to understand why I have negative number? If my file contain errors that caused the results to be different? I will really appreciate any input.

Thank you

wipe
# This a texas standard single bent bridge for a 28' roadway.
# The Column is circular and the bent is square shape.
# In the Future we will set up a texas standard bridge with random node
# that can change at all times
#Units in kip/ft
# column.tcl
# No overlay load was considered in this model.
# This model is a 2D model

# SET UP ----------------------------------------------------------------------------
model basic -ndm 2 -ndf 3; # 2 spacial dimensions, 3 DOF's per nodey
set dataDir Data; # set up name for data directory
file mkdir $dataDir/; # create data directory
set GMdir "../GMfiles"; # ground-motion file directory
# -----------------------------------------------------
# --------------------------------------------------------------------------------------------------
# ------------- define system of units-----------------------------------------
#
# define UNITS ----------------------------------------------------------------------------
set in 1.; # define basic units -- output units
set kip 1.; # define basic units -- output units
set sec 1.; # define basic units -- output units
set LunitTXT "inch"; # define basic-unit text for output
set FunitTXT "kip"; # define basic-unit text for output
set TunitTXT "sec"; # define basic-unit text for output
set ft [expr 12.*$in]; # define engineering units
set ksi [expr $kip/pow($in,2)];
set psi [expr $ksi/1000.];
set lbf [expr $psi*$in*$in]; # pounds force
set pcf [expr $lbf/pow($ft,3)]; # pounds per cubic foot
set psf [expr $lbf/pow($ft,3)]; # pounds per square foot
set in2 [expr $in*$in]; # inch^2
set in4 [expr $in*$in*$in*$in]; # inch^4
set cm [expr $in/2.54]; # centimeter, needed for displacement input in MultipleSupport excitation
set PI [expr 2*asin(1.0)]; # define constants
#set g [expr 32.2*$ft/pow($sec,2)]; # gravitational acceleration
set g 386.4; # g.
set Ubig 1.e10; # a really large number

# Superstructure Weight+Rail+ Haunch; We neglected the bearing pad weight
# build important quantities from input paramters
set span [expr 120]; #forward and backward total span length
set triblength [expr $span/2]
set ngir 4; # the number of girder
set roadway 38; # roadway width
set girderw 0.851; # Standard Tx54 girder in ksi
set bentcap 3.5 ; #rectangular bent cap
set deckwidth 40; # width of the deck
set deck [expr 8.5/12] ; # Typical deck thickness in Texas
set trib [expr $deck*$triblength*$deckwidth]; # tributary area
set haunch 1.1; # factor for haunch dead load to be considered
set gamma 0.15; # the unit weight of concrete is in kcf)
set P1 [expr $haunch*$trib*$gamma]; # deck weight
set P2 [expr $ngir*$girderw*$triblength]; # girder weight kip/ft
set capweight [expr $bentcap*$bentcap*$roadway*$gamma]
set rwt [expr 2*0.382*$triblength]; # Standard T551 rail used in texas
#set P [expr $P1+$rwt+$P2+$capweight]; # superstructure weight

set P 800; # Load used by Choe and al.2008

# define section geometry
set Dcol [expr 5.5*$ft]; # Column Diameter
set LCol [expr 25*$ft]; # column length
set Weight [expr $P]
# calculated parameters
set Mass [expr $P/$g]; # nodal mass

# nodal coordinates:
node 1 0 0; # node#, X, Y
node 2 0 $LCol

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

# nodal masses:
mass 2 $Mass 1e-9 0.;

# we need to set up parameters that are particular to the model.
set IDctrlNode 2; # node where displacement is read for displacement control
set IDctrlDOF 1; # degree of freedom of displacement read for displacement control
set iSupportNode "1"; # define support node, if needed.

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

# Define Column geometry
set coverSec [expr 1.5*$in]; # Column cover to reinforcing steel NA.
set numBarsSec 40; # number of longitudinal-reinforcement bars in the column section
set barAreaSec [expr 1.56*$in*$in] ; # area of longitudinal-reinforcement bars #11
set fy [expr 65*$ksi]; # STEEL yield stress

# MATERIAL parameters -------------------------------------------------------------------
set IDconcC 1; # material ID tag -- confined core concrete
set IDconcU 2; # material ID tag -- unconfined cover concrete
set IDreinf 3; # material ID tag -- reinforcement

# nominal concrete compressive strength
set fc [expr 4.0*$ksi]; # CONCRETE Compressive Strength (+Tension, -Compression)
set E [expr (57*sqrt($fc*1000/$ksi))*$ksi]; # Concrete Elastic Modulus

# Define uniaxial materials

# Mander concrete model for confined and unconfined concrete
# requires rhos_trans = 4 Asp / (ds s)
# where Asp = cross sectional area of spiral reinforcement
# ds = center to center diameter inside spirals
# s = center to center spacing of spiral
set trans_bar 6.0
set As $barAreaSec
set rhos_trans 0.007
set pi [expr acos(-1.0)]
set epss 0.005
set eps0 0.002
set epsu 0.02
set f1 [expr 1.0/2.0*$rhos_trans*$fy]
set dcs [expr $Dcol-2.0*$coverSec-($trans_bar/8.0)]
set rho_cc [expr $As/($pi/4.0*pow($dcs,2))]
set spacing [expr $pi*pow($trans_bar/8.0,2)/$rhos_trans/$dcs]
set ke [expr (1.0-($spacing-$trans_bar/8.0)/2.0/$dcs)/(1.0-$rho_cc)]
set fcc [expr $fc*(-1.254+2.254*sqrt(1.0+7.94*$ke*$f1/$fc)-2.0*$ke*$f1/$fc)]
if { $fcc > 2.23*$fc } {
set fcc [expr 2.23*$fc]
}
set epsc [expr $eps0*(1.0+5.0*($fcc/$fc-1.0))]
set ecr [expr $E/($E-$fcc/$epsc)]
set fcu [expr $fcc*$epsu/$epsc*$ecr/($ecr-1.0+pow($epsu/$epsc,$ecr))]

# Cover concrete
# tag -f'c -epsco -f'cu -epscu
uniaxialMaterial Concrete01 1 -$fc -$eps0 0.0 -$epss
# Core concrete
uniaxialMaterial Concrete01 2 -$fcc -$epsc -$fcu -$epsu


# ----------- Reinforcing Steel --------
set Es [expr 29000.*$ksi]; # modulus of steel
set Esh [expr 2000.0*$ksi]; # Initial Hardening Modulus (estimated)
set Bs 0.01; # 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

uniaxialMaterial Steel02 $IDreinf $fy $Es $Bs $R0 $cR1 $cR2;

# Generate a circular reinforced concrete section
# with one layer of steel evenly distributed around the perimeter and a confined core.
# confined core.
# by: Michael H. Scott, 2003
#
#
# Notes
# The center of the reinforcing bars are placed at the inner radius
# The core concrete ends at the inner radius (same as reinforcing bars)
# The reinforcing bars are all the same size
# The center of the section is at (0,0) in the local axis system
# Zero degrees is along section y-axis
#
set ri 0.0; # inner radius of the section, only for hollow sections
set ro [expr $Dcol/2]; # overall (outer) radius of the section
set nfCoreR 96; # number of radial divisions in the core (number of "rings")
set nfCoreT 36; # number of theta divisions in the core (number of "wedges")
set nfCoverR 24; # number of radial divisions in the cover
set nfCoverT 36; # number of theta divisions in the cover


# Define the fiber section
section fiberSec $ColSecTag {
set rc [expr $ro-$coverSec]; # Core radius
patch circ $IDconcC $nfCoreT $nfCoreR 0 0 $ri $rc 0 360; # Define the core patch
patch circ $IDconcU $nfCoverT $nfCoverR 0 0 $rc $ro 0 360; # Define the cover patch
set theta [expr 360.0/$numBarsSec]; # Determine angle increment between bars
layer circ $IDreinf $numBarsSec $barAreaSec 0 0 $rc $theta 360; # Define the reinforcing layer
}


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

#Remember for Linear PDelta
# element connectivity:
set numIntgrPts 5; # number of integration points for force-based element
element nonlinearBeamColumn 1 1 2 $numIntgrPts $ColSecTag $ColTransfTag; # self-explanatory when using variables

# Column- Pile element connectivity:

# Create recorder
recorder Node -file $dataDir/DFree.out -time -node 2 -dof 1 2 3 disp; # displacements of free nodes
recorder Node -file $dataDir/DBase.out -time -node 1 -dof 1 2 3 disp; # displacements of support nodes
recorder Node -file $dataDir/RBase.out -time -node 1 -dof 1 2 3 reaction; # support reaction
recorder Drift -file $dataDir/Drift.out -time -iNode 1 -jNode 2 -dof 1 -perpDirn 2 ; # lateral drift
recorder Element -file $dataDir/FCol.out -time -ele 2 globalForce; # element forces -- column
recorder Element -file $dataDir/ForceColSec1.out -time -ele 1 section 1 force; # Column section forces, axial and moment, node i
recorder Element -file $dataDir/DefoColSec1.out -time -ele 1 section 1 deformation; # section deformations, axial and curvature, node i
recorder Element -file $dataDir/ForceColSec$numIntgrPts.out -time -ele 1 section $numIntgrPts force; # section forces, axial and moment, node j
recorder Element -file $dataDir/DefoColSec$numIntgrPts.out -time -ele 1 section 1 deformation; # section deformations, axial and curvature, node j
recorder Element -xml $dataDir/PlasticRotation.out -time -ele 1 plasticRotation; # section deformations, axial and curvature, node j

# define GRAVITY -------------------------------------------------------------
pattern Plain 1 Linear {
load 2 0 -$P 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

print -node 2

puts "Model Built"

# --------------------------------------------------------------------------------------------------
# Example 3. 2D Cantilever -- Static Pushover
# Silvia Mazzoni & Frank McKenna, 2006
# execute this file after you have built the model, and after you apply gravity
#

#source coltest.tcl;


# characteristics of pushover analysis
set Dmax [expr 0.05*$LCol]; # maximum displacement of pushover. push to 10% drift.
set Dincr [expr 0.001*$LCol]; # displacement increment for pushover. you want this to be very small, but not too small to slow down the analysis

# create load pattern for lateral pushover load
set Hload [expr $Weight]; # 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
pattern Plain 200 Linear {; # define load pattern -- generalized
load 2 $Hload 0.0 0.0
}

# STATIC-ANALYSIS parameters
# CONSTRAINTS handler -- Determines how the constraint equations are enforced in the analysis (http://opensees.berkeley.edu/OpenSees/m ... al/617.htm)
# Plain Constraints -- Removes constrained degrees of freedom from the system of equations (only for homogeneous equations)
# Lagrange Multipliers -- Uses the method of Lagrange multipliers to enforce constraints
# Penalty Method -- Uses penalty numbers to enforce constraints --good for static analysis with non-homogeneous eqns (rigidDiaphragm)
# Transformation Method -- Performs a condensation of constrained degrees of freedom
set constraintsType Plain; # default;
constraints $constraintsType

# DOF NUMBERER (number the degrees of freedom in the domain): (http://opensees.berkeley.edu/OpenSees/m ... al/366.htm)
# determines the mapping between equation numbers and degrees-of-freedom
# Plain -- Uses the numbering provided by the user
# RCM -- Renumbers the DOF to minimize the matrix band-width using the Reverse Cuthill-McKee algorithm
numberer Plain

# SYSTEM (http://opensees.berkeley.edu/OpenSees/m ... al/371.htm)
# Linear Equation Solvers (how to store and solve the system of equations in the analysis)
# -- provide the solution of the linear system of equations Ku = P. Each solver is tailored to a specific matrix topology.
# ProfileSPD -- Direct profile solver for symmetric positive definite matrices
# BandGeneral -- Direct solver for banded unsymmetric matrices
# BandSPD -- Direct solver for banded symmetric positive definite matrices
# SparseGeneral -- Direct solver for unsymmetric sparse matrices
# SparseSPD -- Direct solver for symmetric sparse matrices
# UmfPack -- Direct UmfPack solver for unsymmetric matrices
system BandGeneral

# 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
set Tol 1.e-8; # Convergence Test: tolerance
set maxNumIter 6; # Convergence Test: maximum number of iterations that will be performed before "failure to converge" is returned
set printFlag 0; # Convergence Test: flag used to print information on convergence (optional) # 1: print information on each step;
set TestType EnergyIncr; # Convergence-test type
test $TestType $Tol $maxNumIter $printFlag;

# 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
set algorithmType Newton
algorithm $algorithmType;

# 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.
analysis Static


# --------------------------------- perform Static Pushover Analysis
set Nsteps [expr int($Dmax/$Dincr)]; # number of pushover analysis steps
set ok [analyze $Nsteps]; # this will return zero if no convergence problems were encountered

if {$ok != 0} {
# if analysis fails, we try some other stuff, performance is slower inside this loop
set ok 0;
set controlDisp 0.0;
set D0 0.0; # analysis starts from zero
set Dstep [expr ($controlDisp-$D0)/($Dmax-$D0)]
while {$Dstep < 1.0 && $ok == 0} {
set controlDisp [nodeDisp $IDctrlNode $IDctrlDOF ]
set Dstep [expr ($controlDisp-$D0)/($Dmax-$D0)]
set ok [analyze 1 ]
if {$ok != 0} {
puts "Trying Newton with Initial Tangent .."
test NormDispIncr $Tol 2000 0
algorithm Newton -initial
set ok [analyze 1 ]
test $TestType $Tol $maxNumIter 0
algorithm $algorithmType
}
if {$ok != 0} {
puts "Trying Broyden .."
algorithm Broyden 8
set ok [analyze 1 ]
algorithm $algorithmType
}
if {$ok != 0} {
puts "Trying NewtonWithLineSearch .."
algorithm NewtonLineSearch .8
set ok [analyze 1 ]
algorithm $algorithmType
}
}; # end while loop
}; # end if ok !0

puts "Pushover Done. Control Disp=[nodeDisp $IDctrlNode $IDctrlDOF]"

recorder display "Column" 10 10 600 600 -wipe
vup 0 1 0
prp 0 0 200
display 1 0 50
print -node 2