interesting error in factorization only for certain loads

Forum for OpenSees users to post questions, comments, etc. on the use of the OpenSees interpreter, OpenSees.exe

Moderators: silvia, selimgunay, Moderators

Post Reply
acv
Posts: 6
Joined: Wed Mar 05, 2014 9:32 am
Location: NTUA

interesting error in factorization only for certain loads

Post by acv »

hello everyone,
i am using the moment curvature example http://opensees.berkeley.edu/wiki/index ... re_Example.
i also changed the section and used a while loop to apply different axial loads.

for axial load=809.37 or for axial load= 1618.74 (=2*809.37) the analysis returns this message
NewtonRaphson:..... the convergence test object failed in test<>
....Error 1 returned in factorization dgstrf

for axial load=810 and for axial load= 1619 the analysis doesn't return those errors.

is there something that i should change in my model??

set v 3.0
puts "v equals $v"
set v 3.0
set sum [expr $v + 2.0]
puts "sum equals $sum"; # print the sum
# units: kip, in

# Remove existing model
wipe

# Create ModelBuilder (with two-dimensions and 3 DOF/node)
model BasicBuilder -ndm 2 -ndf 3

# Define materials for nonlinear columns
# ------------------------------------------
# CONCRETE tag f'c ec0 f'cu ecu
# Core concrete (confined)
uniaxialMaterial ElasticMultiLinear 1 -strain -0.0152 -0.0151 -0.015 -0.00679 -0.004535 -0.003391 -0.0024 -0.001326 -0.00063 0.00001 0.0002 -stress 0 0 -4.307 -5.108 -5.19 -5.033 -4.613 -3.46 -2.021 0 0
# Cover concrete (unconfined)
uniaxialMaterial ElasticMultiLinear 2 -strain -0.006182 -0.006181 -0.0037 -0.00212 -0.001554 -0.001083 -0.000637 0.0001 0.0022 -stress 0 0 -3.182 -3.991 -3.849 -3.247 -2.164 0 0

# STEEL
# Reinforcing steel
set fy 75.0; # Yield stress
set E 29000.0; # Young's modulus
# tag fy E0 b
uniaxialMaterial ElasticMultiLinear 3 -strain -0.062 -0.061 -0.060 -0.015 -0.0053 -0.00258 0.00258 0.060 0.061 0.062 -stress 0 0 -1 -10 -77 -75 75 100 0 0

# ------------------------------------------

# set some paramaters
set colWidth 19.685
set colDepth 19.685

set cover 1.0
set As 0.488; # area of no. 7 bars

# some variables derived from the parameters
set y1 [expr $colDepth/2.0]
set z1 [expr $colWidth/2.0]
# Define cross-section for nonlinear columns

section Fiber 1 {

# Create the concrete core fibers
patch rect 1 50 10 [expr $cover-$y1] [expr $cover-$z1] [expr $y1-$cover] [expr $z1-$cover]

# Create the concrete cover fibers (top, bottom, left, right)
patch rect 2 40 10 [expr -$y1] [expr $z1-$cover] $y1 $z1
patch rect 2 40 10 [expr -$y1] [expr -$z1] $y1 [expr $cover-$z1]
patch rect 2 40 10 [expr -$y1] [expr $cover-$z1] [expr $cover-$y1] [expr $z1-$cover]
patch rect 2 40 10 [expr $y1-$cover] [expr $cover-$z1] $y1 [expr $z1-$cover]

# Create the reinforcing fibers (left, middle, right)
layer straight 3 5 $As [expr $y1-$cover] [expr $z1-$cover] [expr $y1-$cover] [expr $cover-$z1]
layer straight 3 5 $As [expr $cover-$y1] [expr $z1-$cover] [expr $cover-$y1] [expr $cover-$z1]

}

# Estimate yield curvature
# (Assuming no axial load and only top and bottom steel)
set d [expr $colDepth-$cover] ;# d -- from cover to rebar
set epsy [expr $fy/$E] ;# steel yield strain
set Ky [expr $epsy/(0.7*$d)]

# Print estimate to standard output
puts "Estimated yield curvature: $Ky"

# Set axial load
set P -1618.74

set mu 15; # Target ductility for analysis
set numIncr 100; # Number of analysis increments


proc MomentCurvature {secTag axialLoad maxK {numIncr 100} } {
# Define two nodes at (0,0)
node 1 0.0 0.0
node 2 0.0 0.0

# Fix all degrees of freedom except axial and bending
fix 1 1 1 1
fix 2 0 1 0

# Define element
# tag ndI ndJ secTag
element zeroLengthSection 1 1 2 $secTag

# Create recorder
recorder Node -file section2.out -time -node 2 -dof 3 disp
recorder Element -file Element1.out -time -ele 1 section fiber 8.65 0.1 1 strain
recorder Element -file Element2.out -time -ele 1 section fiber 8.55 0.1 1 strain
pattern Plain 1 "Constant" {
load 2 $axialLoad 0.0 0.0
}

# Define analysis parameters
integrator LoadControl 0.0
system SparseGeneral -piv; # Overkill, but may need the pivoting!
test NormUnbalance 1.0e-11 10000
numberer Plain
constraints Plain
algorithm Newton
analysis Static

# Do one analysis for constant axial load
analyze 1

# Define reference moment
pattern Plain 2 "Linear" {
load 2 0.0 0.0 1.0
}

# Compute curvature increment
set dK [expr $maxK/$numIncr]

# Use displacement control at node 2 for section analysis
integrator DisplacementControl 2 3 $dK

# Do the section analysis
analyze $numIncr
}
# Call the section analysis procedure
MomentCurvature 1 $P [expr $Ky*$mu] $numIncr
Post Reply