Convergence Problem with ElasticPP material

[AbRiahi]

Hi
I want to check the force-displacement relationship associated with a single truss element having ElasticPP material. I wrote a very simple code and I use the proposed method in the manual for improving convergence at the end of my code, but unforthunately it cannot compelet the analysis despite capturing the plastic region after yielding and I don't know why.
1) How I can tackle the convergence problem to obtain the pushover curve up to the target displacement? (The target displacement is 75mm, the yield point occurs at 8.067mm but it fails to converge after 9mm)
2) I don't know why the force-displacement curve does not include (0,0) point and starts at (3,74.39) instead?
3) Why the Opensees warns the message of no integrator?
The code and the results are presented as the following:

wipe; # clear opensees model
model basic -ndm 2 -ndf 3; # 2 dimensions, 3 dof per node
file mkdir data; # create data directory

# define GEOMETRY -------------------------------------------------------------

set L 4150.0; # Bay width, mm
set h1 2667.0; # 1st story col. height, mm
set h2 2700.0; # 2nd story col. height, mm
set h3 2700.0; # 3rd story col. height, mm
set Lev1 0.0;
set Lev2 [expr $Lev1+$h1];
set Lev3 [expr $Lev2+$h2];
set Lev4 [expr $Lev3+$h3];
set PinLev 200;

# nodal coordinates:

node 78 [expr 0.5*$L] $Lev4
set PT_Level 200;
node 120 [expr 0.5*$L] [expr $Lev1-$PT_Level]
puts "nodes are succesfully assigned"
fix 120 1 1 1; # node DX DY RZ
fix 78 1 0 1; # node DX DY RZ
puts "Supports are succesfully assigned at restrained nodes"
geomTransf Linear 1; # associate a tag to transformation
set E 2.05e5; # Steel modolus of elasticity, MPa
set A_PT 1.0;
uniaxialMaterial ElasticPP 1000 $E 0.001;
element truss 10 120 78 $A_PT 1000
test EnergyIncr 1.0e-9 10
numberer Plain
constraints Plain
algorithm Newton
analysis Static

recorder Node -file Data/Displ.out -time -node 78 -dof 1 2 3 disp;
recorder Node -file Data/React.out -time -node 120 -dof 1 2 3 reaction;
pattern Plain 1 "Linear" {
load 78 0.0 1.0 0.0
}
integrator DisplacementControl 78 2 3
analyze 25
set maxK 30.0;
puts "maxK: $maxK"
set numIncr 10;
puts "numIncr:$numIncr";
set dK [expr $maxK/$numIncr];
puts "dK: $dK";
integrator DisplacementControl 78 2 $dK
# Do the section analysis
set ok [analyze $numIncr]

# ----------------------------------------------if convergence failure-------------------------
set IDctrlNode 78
set IDctrlDOF 2
set Dmax $maxK
set Dincr $dK
set TolStatic 1.e-9;
set testTypeStatic EnergyIncr
set maxNumIterStatic 20
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
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]
}

puts "Done!"

The obtained results are the following files:
Displ.out
74.3922 0 3 0
148.784 0 6 0
167.382 0 6.75 0
185.98 0 7.5 0
204.578 0 8.25 0
205 0 9 0
React.out
74.3922 0 -74.3922 0
148.784 0 -148.784 0
167.382 0 -167.382 0
185.98 0 -185.98 0
204.578 0 -204.578 0
205 0 -205 0

[vesna]

1) I modified your scrip so it does what you want. You only need to adjust number of analysis steps to get to the target displacement.
2) The first point you record is the result after 1st analysis step
3) You need to specify integrator before analysis command


wipe; # clear opensees model
model basic -ndm 2 -ndf 3; # 2 dimensions, 3 dof per node
file mkdir data; # create data directory

# define GEOMETRY -------------------------------------------------------------

set L 4150.0; # Bay width, mm
set h1 2667.0; # 1st story col. height, mm
set h2 2700.0; # 2nd story col. height, mm
set h3 2700.0; # 3rd story col. height, mm
set Lev1 0.0;
set Lev2 [expr $Lev1+$h1];
set Lev3 [expr $Lev2+$h2];
set Lev4 [expr $Lev3+$h3];
set PinLev 200;

# nodal coordinates:

node 78 [expr 0.5*$L] $Lev4
set PT_Level 200;
node 120 [expr 0.5*$L] [expr $Lev1-$PT_Level]
puts "nodes are succesfully assigned"
fix 120 1 1 1; # node DX DY RZ
fix 78 1 0 1; # node DX DY RZ
puts "Supports are succesfully assigned at restrained nodes"
geomTransf Linear 1; # associate a tag to transformation
set E 2.05e5; # Steel modolus of elasticity, MPa
set A_PT 1.0;
uniaxialMaterial ElasticPP 1000 $E 0.001;
element truss 10 120 78 $A_PT 1000


recorder Node -file Data/Displ.out -time -node 78 -dof 1 2 3 disp;
recorder Node -file Data/React.out -time -node 120 -dof 1 2 3 reaction;

pattern Plain 1 "Linear" {
sp 78 2 1.0
}

test EnergyIncr 1.0e-9 10
numberer Plain
constraints Penalty 1.0e+18 1.0e+18
system BandSPD
integrator LoadControl 0.1
#integrator DisplacementControl 78 2 0.1
algorithm Newton
analysis Static

analyze 2000

print -node 78
print -node 120

[AbRiahi]

Dear vesna
Thanks a lot for you reply. Now I can run the file you sent without any convergence error, But I have some questions:

1) Are the commands you added like "constraints Penalty 1.0e+18 1.0e+18" and "system BandSPD" are responsible for removing the
convergence problem? If yes why they should be used for ElasticPP material and how they solve the convergence problem?(Because I did not have any convergence problem with elements such as truss,beamcolumn,...,even with nonlinear materials, by using commands like "constraints Plain" and "System BandGeneral,...") ? Are there any rules for selecting these parameters appropriately for any analysis type or element?
2) I want to use this ElasticPP as a material for an element in a more complicated structural model subjected to earthquake ground motion. Can the same options be selected for "constraint" and "system" commands in THA (Time History Analysis) in order to avoid convergence problems?( I think the manual suggested not to use "constraints Penalty " for such a case)

[vesna]

1) you did not have convergence problems but stability problems. I solved it by "sp" command instead of "load" command inside of load pattern and "load control" instead of "displacement control". To be able to apply displacement using "sp" I had to use penalty constraint. For a material that has post-yielding slope of zero I had to use system BandSPD.

2) for time history analysis you do not have to use penalty constraint.

[tazarv]

if you don't like the changes that Vesna has done for you, there is also a tricky way. Instead of ElasticPP you can use Steel01 with small strain-hardening ratio.
http://opensees.berkeley.edu/wiki/index ... 1_Material

[AbRiahi]

Thanks a lot for your kindness and reply
I should say that I like the way that vesna mentined, but I should say that I changed the integrator again to diplacement control and I did not notice any stability problem.
By the way , I think for the dynamic analysis I can use "Transformation constraint" and "BandSPD" ,as vesna mentioned, and I should not encounter any stability problem.
Best Regards

[AbRiahi]

Dear vesna
I used the same parameters that you proposed for the EPP case for a simple case of zerolength element with Elastic No Tension (ENT)material, but the opensees returns warning message. Should I make any modification in the file in order to complete the analysis? Here is the file:

wipe; # clear opensees model
model basic -ndm 2 -ndf 3; # 2 dimensions, 3 dof per node
file mkdir data; # create data directory

# define GEOMETRY -------------------------------------------------------------
set e1 150.0;
set e2 100.0;
set e3 100.0;
set L 4150.0; # Bay width, mm
set h1 2667.0; # 1st story col. height, mm
set h2 2700.0; # 2nd story col. height, mm
set h3 2700.0; # 3rd story col. height, mm
set Lev1 0.0;
set Lev2 [expr $Lev1+$h1];
set Lev3 [expr $Lev2+$h2];
set Lev4 [expr $Lev3+$h3];
set PinLev 200;

# nodal coordinates:
node 101 [expr 0.0-$e1] $Lev1
node 1001 [expr 0.0-$e1] $Lev1

# Single point constraints -- Boundary Conditions: 1=restrained & o=released
fix 1001 1 1 1; # node DX DY RZ

set kbh 7.8e6; #horizontal support base stiffness
set tagContHor 3;

uniaxialMaterial ENT $tagContHor $kbh; # Horizontal Direction:defining contact material elastic no tension material with high compressional stiffness
element zeroLength 1001 1001 101 -mat $tagContHor -dir 1

recorder Node -file Data/Displ.out -time -node 101 -dof 1 2 3 disp;
recorder Node -file Data/React.out -time -node 1001 -dof 1 2 3 reaction;
pattern Plain 1 "Linear" {
sp 101 1 1.0
}

test EnergyIncr 1.0e-9 10
numberer Plain
constraints Penalty 1.0e+18 1.0e+18
system BandSPD
integrator LoadControl -0.1
algorithm Newton
analysis Static
analyze 1000
integrator LoadControl 0.1
analyze 2000
print -node 101

[vesna]

You need to add constraints to node 101:

fix 101 0 1 1

[AbRiahi]

Thanks a lot for your attention