Moment-Curvature - NOT solving

[IsmBaha]

Hi all,

Here I have a fiber section and I tried to get the moment-curvature of the section. I used two examples provided in the manual and website. I followed same way as in the example. This is a shear wall section with two boundary and one web sections. Only boundaries are confined. I used rectangle patch. I really appreciate if you could help me what is wrong.

Here is the first code:

proc MomentCurvatureSW {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 section$secTag.out -time -node 2 -dof 3 disp

# Define constant axial load
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-9 10
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
}

# 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 Concrete01 1 -6.0 -0.004 -5.0 -0.014

# Cover concrete (unconfined)
uniaxialMaterial Concrete01 2 -5.0 -0.002 0.0 -0.006

# STEEL
# Reinforcing steel
set fy 70.0; # Yield stress
set E 29000.0; # Young's modulus
# tag fy E0 b
uniaxialMaterial Steel01 3 $fy $E 0.01

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

# set some paramaters
set colWidth 20
set colDepth [expr 44*12]

set cover 1.5
set As 0.60; # 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 $matTag $numSubdivY $numSubdivZ $yI $zI $yJ $zJ
patch rect 1 20 4 [expr 169.5] [expr 8.5] [expr 262.5] [expr -8.5]
patch rect 1 20 4 [expr -262.5] [expr 8.5] [expr -169.5] [expr -8.5]

# Create the concrete cover fibers (top, bottom, left, right)
patch rect 2 40 1 [expr -264.0] [expr 8.5] 264.0 10.0; # 1
patch rect 2 40 1 [expr -264.0] [expr -8.5] 264.0 [expr -8.5]; #2
patch rect 2 4 1 [expr -264.0] [expr -8.5] [expr -262.5] [expr 8.5]; # 3
patch rect 2 4 1 [expr 262.5] [expr -8.5] 264.0 [expr 8.5]; # 4
patch rect 2 40 4 [expr -169.5] [expr -8.5] 169.5 [expr 8.5]; # 5

# Create the reinforcing fibers (left, middle, right)
# # layer straight $matTag $numBars $areaBar $yStart $zStart $yEnd $zEnd
layer straight 3 16 $As [expr 262.5] [expr 8.5] [expr 169.5] [expr 8.5]; # 1
layer straight 3 16 $As [expr 262.5] [expr -8.5] [expr 169.5] [expr -8.5]; # 2
layer straight 3 2 $As [expr 262.5] [expr 0.0] [expr 169.5] [expr 0.0]; # 2-1
layer straight 3 16 $As [expr -169.5] [expr 8.5] [expr -262.5] [expr 8.5]; # 3
layer straight 3 16 $As [expr -169.5] [expr -8.5] [expr -262.5] [expr -8.5]; # 4
layer straight 3 2 $As [expr -169.5] [expr 0.0] [expr -262.5] [expr 0.0]; # 4-1
layer straight 3 56 $As [expr 169.5] [expr 8.5] [expr -169.5] [expr 8.5]; # 5
layer straight 3 56 $As [expr 169.5] [expr -8.5] [expr -169.5] [expr -8.5]; # 6
}

# Estimate yield curvature
# (Assuming no axial load and only top and bottom steel)
set d [expr $colDepth-$cover]; # d -- from top cover to bottom 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 -2903.0

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

# Call the section analysis procedure
MomentCurvatureSW 1 $P [expr $Ky*$mu] $numIncr

[IsmBaha]

Here is the second code:


proc MomentCurvatureSW2D { secTag axialLoad maxK {numIncr 100} } {
##################################################
# A procedure for performing section analysis (only does
# moment-curvature, but can be easily modified to do any mode
# of section reponse.)
#
# MHS
# October 2000
# modified to 2D and to improve convergence by Silvia Mazzoni, 2006
#
# Arguments
# secTag -- tag identifying section to be analyzed
# axialLoad -- axial load applied to section (negative is compression)
# maxK -- maximum curvature reached during analysis
# numIncr -- number of increments used to reach maxK (default 100)
#
# Sets up a recorder which writes moment-curvature results to file
# section$secTag.out ... the moment is in column 1, and curvature in column 2

# Define two nodes at (0,0)
node 1001 0.0 0.0
node 1002 0.0 0.0

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

# Define element
# tag ndI ndJ secTag
element zeroLengthSection 2001 1001 1002 $secTag

# Create recorder
recorder Node -file data/Mphi.out -time -node 1002 -dof 3 disp; # output moment (col 1) & curvature (col 2)

# Define constant axial load
pattern Plain 3001 "Constant" {
load 1002 $axialLoad 0.0 0.0
}

# Define analysis parameters
integrator LoadControl 0 1 0 0
system SparseGeneral -piv; # Overkill, but may need the pivoting!
test EnergyIncr 1.0e-9 10
numberer Plain
constraints Plain
algorithm Newton
analysis Static

# Do one analysis for constant axial load
analyze 1

loadConst -time 0.0

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

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

# Use displacement control at node 1002 for section analysis, dof 3
integrator DisplacementControl 1002 3 $dK 1 $dK $dK

# Do the section analysis
set ok [analyze $numIncr]

# ----------------------------------------------if convergence failure-------------------------
set IDctrlNode 1002
set IDctrlDOF 3
set Dmax $maxK
set Dincr $dK
set TolStatic 1.e-9;
set testTypeStatic EnergyIncr
set maxNumIterStatic 6
set algorithmTypeStatic Newton
if {$ok != 0} {
# if analysis fails, we try some other stuff, performance is slower inside this loop
set Dstep 0.0;
set ok 0
while {$Dstep <= 1.0 && $ok == 0} {
set controlDisp [nodeDisp $IDctrlNode $IDctrlDOF ]
set Dstep [expr $controlDisp/$Dmax]
set ok [analyze 1]; # this will return zero if no convergence problems were encountered
if {$ok != 0} {; # reduce step size if still fails to converge
set Nk 4; # reduce step size
set DincrReduced [expr $Dincr/$Nk];
integrator DisplacementControl $IDctrlNode $IDctrlDOF $DincrReduced
for {set ik 1} {$ik <=$Nk} {incr ik 1} {
set ok [analyze 1]; # 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 global maxNumIterStatic; # max no. of iterations performed before "failure to converge" is ret'd
puts "Trying Newton with Initial Tangent .."
test NormDispIncr $TolStatic 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} {; # stop if still fails to converge
puts [format $fmt1 "PROBLEM" $IDctrlNode $IDctrlDOF [nodeDisp $IDctrlNode $IDctrlDOF] $LunitTXT]
return -1
}; # end if
}; # end for
integrator DisplacementControl $IDctrlNode $IDctrlDOF $Dincr; # bring back to original increment
}; # end if
}; # end while loop
}; # end if ok !0
# -----------------------------------------------------------------------------------------------------
global LunitTXT; # load time-unit text
if { [info exists LunitTXT] != 1} {set LunitTXT "Length"}; # set blank if it has not been defined previously.

set fmt1 "%s Pushover analysis: CtrlNode %.3i, dof %.1i, Curv=%.4f /%s"; # format for screen/file output of DONE/PROBLEM analysis
if {$ok != 0 } {
puts [format $fmt1 "PROBLEM" $IDctrlNode $IDctrlDOF [nodeDisp $IDctrlNode $IDctrlDOF] $LunitTXT]
} else {
puts [format $fmt1 "DONE" $IDctrlNode $IDctrlDOF [nodeDisp $IDctrlNode $IDctrlDOF] $LunitTXT]
}

}

# 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 Concrete01 1 -6.0 -0.004 -5.0 -0.014

# Cover concrete (unconfined)
uniaxialMaterial Concrete01 2 -5.0 -0.002 0.0 -0.006

# STEEL
# Reinforcing steel
set fy 70.0; # Yield stress
set E 29000.0; # Young's modulus
# tag fy E0 b
uniaxialMaterial Steel01 3 $fy $E 0.01

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

# set some paramaters
set colWidth 20
set colDepth [expr 44*12]

set cover 1.5
set As 0.60; # 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 $matTag $numSubdivY $numSubdivZ $yI $zI $yJ $zJ
patch rect 1 20 4 [expr 169.5] [expr 8.5] [expr 262.5] [expr -8.5]
patch rect 1 20 4 [expr -262.5] [expr 8.5] [expr -169.5] [expr -8.5]

# Create the concrete cover fibers (top, bottom, left, right)
patch rect 2 40 1 [expr -264.0] [expr 8.5] 264.0 10.0; # 1
patch rect 2 40 1 [expr -264.0] [expr -8.5] 264.0 [expr -8.5]; #2
patch rect 2 4 1 [expr -264.0] [expr -8.5] [expr -262.5] [expr 8.5]; # 3
patch rect 2 4 1 [expr 262.5] [expr -8.5] 264.0 [expr 8.5]; # 4
patch rect 2 40 4 [expr -169.5] [expr -8.5] 169.5 [expr 8.5]; # 5

# Create the reinforcing fibers (left, middle, right)
# # layer straight $matTag $numBars $areaBar $yStart $zStart $yEnd $zEnd
layer straight 3 16 $As [expr 262.5] [expr 8.5] [expr 169.5] [expr 8.5]; # 1
layer straight 3 16 $As [expr 262.5] [expr -8.5] [expr 169.5] [expr -8.5]; # 2
layer straight 3 2 $As [expr 262.5] [expr 0.0] [expr 169.5] [expr 0.0]; # 2-1
layer straight 3 16 $As [expr -169.5] [expr 8.5] [expr -262.5] [expr 8.5]; # 3
layer straight 3 16 $As [expr -169.5] [expr -8.5] [expr -262.5] [expr -8.5]; # 4
layer straight 3 2 $As [expr -169.5] [expr 0.0] [expr -262.5] [expr 0.0]; # 4-1
layer straight 3 56 $As [expr 169.5] [expr 8.5] [expr -169.5] [expr 8.5]; # 5
layer straight 3 56 $As [expr 169.5] [expr -8.5] [expr -169.5] [expr -8.5]; # 6
}

# Estimate yield curvature
# (Assuming no axial load and only top and bottom steel)
set d [expr $colDepth-$cover]; # d -- from top cover to bottom 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 -2903.0

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

# Call the section analysis procedure
MomentCurvatureSW2D 1 $P [expr $Ky*$mu] $numIncr

[IsmBaha]

I found a coordinate mistake in both codes:

patch rect 2 40 1 [expr -264.0] [expr -8.5] 264.0 [expr -8.5]; #2
has to be
patch rect 2 40 1 [expr -264.0] [expr -10.0] 264.0 [expr -8.5]; #2

This was causing a node to hang on the air so discontinuity.

And I had to change the material properties. I used steel02 and concrete02 and they work.

I think I will have more similar questions on the process of mastering the opensees. Sorry for spamming the forum

[aahashib]

Hi. I am trying to get moment-curvature. But I am not understanding the following commands. What do the following lines mean? Thanks.

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

[IsmBaha]

- mu is ductility ratio = Displacement/Yield Displacement
- number of steps: in how many steps this target ductility (mu) will be achieved

Regards,

[aahashib]

Thanks