Convergence problem

[mmasoudi]

I changed several parameters in order to achieve convergence. However, I found the earthquake records should be scaled down to converge. The earthquake loading is bidirectional (El-centro motion 1940 downloaded from peer database) and the unit of my model is kN, m.
Could you please help me on this issue?



wipe; # clear memory of all past model definitions
model BasicBuilder -ndm 3 -ndf 6; # Define the model builder, ndm=#dimension, ndf=#dofs
set dataDir Data; # set up name of data directory -- remove
file mkdir $dataDir; # create data directory
source DisplayPlane.tcl; # procedure for displaying a plane in model
source DisplayModel3D.tcl; # procedure for displaying a plane in model

# Define GEOMETRY -------------------------------------------------------------
# Define structure-geometry paramters
set LCol 4; # column height (parallel to Y axis)
set LBeam 6; # beam length (parallel to X axis)
set LGird 6; # girder length (parallel to Z axis)

# ------ frame configuration
set NStory 1; # number of stories above ground level
set NBay 1; # number of bays in X direction
set NBayZ 1; # number of bays in Z direction
puts "Number of Stories in Y: $NStory; Number of bays in X: $NBay; Number of bays in Z: $NBayZ"

# Define NODAL COORDINATES
# Ground level

node 111 0.0 0 0.0
node 112 6.0 0 0.0
node 211 0.0 0 6.0
node 212 6.0 0 6.0

# Floor 1

node 121 0.0 4 0.0 #numbering method, X-prallelGirdNo.,StoryNo.,Z-prallelGirdNo.
node 122 6.0 4 0.0
node 221 0.0 4 6.0
node 222 6.0 4 6.0

node 1101 0.9 4 0.0 #numbering method, FloorNo.,X-prallelGirdNo. if not 0,Z-prallelGirdNo. if not 0, cut place
node 1102 3.0 4 0.0
node 1103 5.1 4 0.0
node 1201 0.9 4 6.0
node 1202 3.0 4 6.0
node 1203 5.1 4 6.0

node 1011 0.0 4 0.9 #numbering method, FloorNo.,X-prallelGirdNo. if not 0,Z-prallelGirdNo. if not 0, cut place
node 1012 0.0 4 3.0
node 1013 0.0 4 5.1
node 1021 6.0 4 0.9
node 1022 6.0 4 3.0
node 1023 6.0 4 5.1

# Lumped nodal masses excluding self mass of elements
# End span masses
mass 1101 6.238 6.238 6.238 0.0 0.0 0.0;
mass 1103 6.238 6.238 6.238 0.0 0.0 0.0;
mass 1201 6.238 6.238 6.238 0.0 0.0 0.0;
mass 1203 6.238 6.238 6.238 0.0 0.0 0.0;
mass 1011 6.238 6.238 6.238 0.0 0.0 0.0;
mass 1013 6.238 6.238 6.238 0.0 0.0 0.0;
mass 1021 6.238 6.238 6.238 0.0 0.0 0.0;
mass 1023 6.238 6.238 6.238 0.0 0.0 0.0;

# Mid span masses
mass 1102 5.872 5.872 5.872 0.0 0.0 0.0;
mass 1202 5.872 5.872 5.872 0.0 0.0 0.0;
mass 1012 5.872 5.872 5.872 0.0 0.0 0.0;
mass 1022 5.872 5.872 5.872 0.0 0.0 0.0;

# Define Rigid Floor Diaphragm
rigidDiaphragm 2 121 122 221 222 1101 1102 1103 1201 1202 1203 1011 1012 1013 1021 1022 1023

# BOUNDARY CONDITIONS
fixY 0.0 1 1 1 1 1 1; # fix all Y=0.0 nodes

# define section tags:
set ColSecTagFiber 4
set BeamSecTagFiber 5
set SecTagTorsion 70

# Section Properties:
set HCol 0.3; # square-Column width
set BCol $HCol
set HBeam 0.4; # Beam depth -- perpendicular to bending axis
set BBeam 0.3; # Beam width -- parallel to bending axis


# MATERIAL parameters
uniaxialMaterial Concrete01 1 -27579.0 -0.0022 -22063.0 -0.005;
uniaxialMaterial Steel01 2 413685.5 2.000000E+8 0.01 0.0 1.0 0.0 1.0;

uniaxialMaterial Elastic 3 2300.0 0.0; # Reduced torsional Stiffness (0.2GJ) of the 30x30 RC section
uniaxialMaterial Elastic 4 4000.0 0.0; # Reduced torsional Stiffness (0.2GJ) of the 30x40 RC section


# FIBER SECTION properties
section Fiber 3030 {
patch quad 1 14 14 -0.150000 -0.150000 0.150000 -0.150000 0.150000 0.150000 -0.150000 0.150000;
layer straight 2 3 0.000314 -0.100000 0.100000 0.100000 0.100000;
layer straight 2 3 0.000314 -0.100000 -0.100000 0.100000 -0.100000;
layer straight 2 2 0.000314 0.100000 -0.000000 -0.100000 -0.000000;
}
section Fiber 3040 {
patch quad 1 14 14 -0.150000 -0.200000 0.150000 -0.200000 0.150000 0.200000 -0.150000 0.200000;
layer straight 2 4 0.000491 -0.100000 -0.150000 0.100000 -0.150000;
layer straight 2 4 0.000491 -0.100000 0.150000 0.100000 0.150000;
layer straight 2 2 0.000491 -0.100000 0.050000 -0.100000 -0.050000;
layer straight 2 2 0.000491 0.100000 0.050000 0.100000 -0.050000;
}

set ColSecTag 1;
set BeamSecTag 2;

section Aggregator $ColSecTag 3 T -section 3030;
section Aggregator $BeamSecTag 4 T -section 3040;

# Elements
set np 5; # number of Gauss integration points for nonlinear curvature distribution
set maxIters 1;
set tol 10e-16;

# set up geometric transformations of element
# separate columns and beams, in case of P-Delta analysis for columns
# in 3D model, assign vector vecxz
set IDColTransf 1; # all columns
set IDBeamXXTransf 2; # all beams
set IDBeamZZTransf 3; # all girders
set geomTransfType Linear;

geomTransf $geomTransfType $IDColTransf 0 0 1; # only columns can have PDelta effects (gravity effects)
geomTransf $geomTransfType $IDBeamXXTransf 0 -1 0;
geomTransf $geomTransfType $IDBeamZZTransf 0 -1 0;

set massDens30x30 0.216; # 0.30x0.30x1x2.4
set massDens30x40 0.288;

# Columns
element forceBeamColumn 111 111 121 $np $ColSecTag $IDColTransf -mass $massDens30x30 -iter $maxIters $tol;
element forceBeamColumn 121 112 122 $np $ColSecTag $IDColTransf -mass $massDens30x30 -iter $maxIters $tol;
element forceBeamColumn 211 211 221 $np $ColSecTag $IDColTransf -mass $massDens30x30 -iter $maxIters $tol;
element forceBeamColumn 221 212 222 $np $ColSecTag $IDColTransf -mass $massDens30x30 -iter $maxIters $tol;

# BeamXX
element forceBeamColumn 1101 121 1101 $np $BeamSecTag $IDBeamXXTransf -mass $massDens30x40 -iter $maxIters $tol;
element forceBeamColumn 1102 1101 1102 $np $BeamSecTag $IDBeamXXTransf -mass $massDens30x40 -iter $maxIters $tol;
element forceBeamColumn 1103 1102 1103 $np $BeamSecTag $IDBeamXXTransf -mass $massDens30x40 -iter $maxIters $tol;
element forceBeamColumn 1104 1103 122 $np $BeamSecTag $IDBeamXXTransf -mass $massDens30x40 -iter $maxIters $tol;

element forceBeamColumn 1201 221 1201 $np $BeamSecTag $IDBeamXXTransf -mass $massDens30x40 -iter $maxIters $tol;
element forceBeamColumn 1202 1201 1202 $np $BeamSecTag $IDBeamXXTransf -mass $massDens30x40 -iter $maxIters $tol;
element forceBeamColumn 1203 1202 1203 $np $BeamSecTag $IDBeamXXTransf -mass $massDens30x40 -iter $maxIters $tol;
element forceBeamColumn 1204 1203 222 $np $BeamSecTag $IDBeamXXTransf -mass $massDens30x40 -iter $maxIters $tol;

# BeamZZ
element forceBeamColumn 1011 121 1011 $np $BeamSecTag $IDBeamZZTransf -mass $massDens30x40 -iter $maxIters $tol;
element forceBeamColumn 1012 1011 1012 $np $BeamSecTag $IDBeamZZTransf -mass $massDens30x40 -iter $maxIters $tol;
element forceBeamColumn 1013 1012 1013 $np $BeamSecTag $IDBeamZZTransf -mass $massDens30x40 -iter $maxIters $tol;
element forceBeamColumn 1014 1013 221 $np $BeamSecTag $IDBeamZZTransf -mass $massDens30x40 -iter $maxIters $tol;

element forceBeamColumn 1021 122 1021 $np $BeamSecTag $IDBeamZZTransf -mass $massDens30x40 -iter $maxIters $tol;
element forceBeamColumn 1022 1021 1022 $np $BeamSecTag $IDBeamZZTransf -mass $massDens30x40 -iter $maxIters $tol;
element forceBeamColumn 1023 1022 1023 $np $BeamSecTag $IDBeamZZTransf -mass $massDens30x40 -iter $maxIters $tol;
element forceBeamColumn 1024 1023 222 $np $BeamSecTag $IDBeamZZTransf -mass $massDens30x40 -iter $maxIters $tol;

# Gravity loads
set uniLoad30x30 [expr -9.81*$massDens30x30];
set uniLoad30x40 [expr -9.81*$massDens30x40];

pattern Plain 1 Linear {

# Point loads
# End span loads
load 1101 0.0 -61.19478 0.0 0.0 0.0 0.0;
load 1103 0.0 -61.19478 0.0 0.0 0.0 0.0;
load 1201 0.0 -61.19478 0.0 0.0 0.0 0.0;
load 1203 0.0 -61.19478 0.0 0.0 0.0 0.0;
load 1011 0.0 -61.19478 0.0 0.0 0.0 0.0;
load 1013 0.0 -61.19478 0.0 0.0 0.0 0.0;
load 1021 0.0 -61.19478 0.0 0.0 0.0 0.0;
load 1023 0.0 -61.19478 0.0 0.0 0.0 0.0;

# Mid span loads
load 1102 0.0 -57.60432 0.0 0.0 0.0 0.0;
load 1202 0.0 -57.60432 0.0 0.0 0.0 0.0;
load 1012 0.0 -57.60432 0.0 0.0 0.0 0.0;
load 1022 0.0 -57.60432 0.0 0.0 0.0 0.0;

# Columns
eleLoad -ele 111 121 211 221 -type -beamUniform 0.0 0.0 $uniLoad30x30;

# BeamXX
eleLoad -ele 1101 1102 1103 1104 -type -beamUniform 0.0 $uniLoad30x40;
eleLoad -ele 1201 1202 1203 1204 -type -beamUniform 0.0 $uniLoad30x40;

# BeamZZ
eleLoad -ele 1011 1012 1013 1014 -type -beamUniform 0.0 $uniLoad30x40;
eleLoad -ele 1021 1022 1023 1024 -type -beamUniform 0.0 $uniLoad30x40;
}

# Define DISPLAY -------------------------------------------------------------
set xPixels 1200; # height of graphical window in pixels
set yPixels 800; # height of graphical window in pixels
set xLoc1 10; # horizontal location of graphical window (0=upper left-most corner)
set yLoc1 10; # vertical location of graphical window (0=upper left-most corner)
set dAmp 2; # scaling factor for viewing deformed shape, it depends on the dimensions of the model
DisplayModel3D NodeNumbers $dAmp $xLoc1 $yLoc1 $xPixels $yPixels

# Gravity load analysis
constraints Lagrange;
numberer Plain;
system BandGeneral;
test EnergyIncr 1.0e-08 10;
algorithm Newton;
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"

# Define RECORDERS -------------------------------------------------------------
recorder Node -file $dataDir/DFree.out -time -node 121 122 221 222 -dof 1 2 3 disp; # displacements of free node
recorder Node -file $dataDir/RBase.out -time -node 111 112 211 212 -dof 1 2 3 reaction; # support reaction

# Define DISPLAY -------------------------------------------------------------
set xPixels 1200; # height of graphical window in pixels
set yPixels 800; # height of graphical window in pixels
set xLoc1 10; # horizontal location of graphical window (0=upper left-most corner)
set yLoc1 10; # vertical location of graphical window (0=upper left-most corner)
set ViewScale 10; # scaling factor for viewing deformed shape, it depends on the dimensions of the model
DisplayModel3D DeformedShape $ViewScale $xLoc1 $yLoc1 $xPixels $yPixels

recorder plot $dataDir/DFree.out TopDispl 10 10 1200 400 -columns 1 2 -columns 1 4; # a window to plot the nodal displacements versus time

# ------------ define & apply damping
# RAYLEIGH damping parameters, Where to put M/K-prop damping, switches (http://opensees.berkeley.edu/OpenSees/m ... l/1099.htm)
# D=$alphaM*M + $betaKcurr*Kcurrent + $betaKcomm*KlastCommit + $beatKinit*$Kinitial
set xDamp 0.05; # damping ratio
set MpropSwitch 1.0;
set KcurrSwitch 0.0;
set KcommSwitch 1.0;
set KinitSwitch 0.0;
set omegaI [expr 2*3.14/0.6];
set omegaJ [expr 2*3.14/1.0];
set alphaM [expr $MpropSwitch*$xDamp*(2*$omegaI*$omegaJ)/($omegaI+$omegaJ)]; # M-prop. damping; D = alphaM*M
set betaKcurr [expr $KcurrSwitch*2.*$xDamp/($omegaI+$omegaJ)]; # current-K; +beatKcurr*KCurrent
set betaKcomm [expr $KcommSwitch*2.*$xDamp/($omegaI+$omegaJ)]; # last-committed K; +betaKcomm*KlastCommitt
set betaKinit [expr $KinitSwitch*2.*$xDamp/($omegaI+$omegaJ)]; # initial-K; +beatKinit*Kini
rayleigh $alphaM $betaKcurr $betaKinit $betaKcomm; # RAYLEIGH damping

timeSeries Path 101 -dt 0.01 -filePath ELC180.txt -factor 9.81; # define acceleration vector from file (dt=0.005 is associated with the input file gm)
timeSeries Path 102 -dt 0.01 -filePath ELC270.txt -factor 9.81;
pattern UniformExcitation 11 1 -accel 101; # define where and how (pattern tag, dof) acceleration is applied
pattern UniformExcitation 12 3 -accel 102;

set DynAlzTol 1.0e-6;
constraints Lagrange;
numberer RCM;
system BandGeneral;
test EnergyIncr $DynAlzTol 200;
algorithm ModifiedNewton;
integrator Newmark 0.5 0.25;
analysis Transient;

set DtAnalysis 0.01;
set TmaxAnalysis 40.00;
set Nsteps [expr int($TmaxAnalysis/$DtAnalysis)];
set ok [analyze $Nsteps $DtAnalysis]; # actually perform analysis; returns ok=0 if analysis was successful

if {$ok != 0} { ; # analysis was not successful.
# --------------------------------------------------------------------------------------------------
# change some analysis parameters to achieve convergence
# performance is slower inside this loop
# Time-controlled analysis
set ok 0;
set controlTime [getTime];
while {$controlTime < $TmaxAnalysis && $ok == 0} {
set controlTime [getTime]
set ok [analyze 1 $DtAnalysis]
if {$ok != 0} {
puts "Trying Newton with Initial Tangent .."
test NormDispIncr $DynAlzTol 1000 0
algorithm Newton -initial
set ok [analyze 1 $DtAnalysis]
test EnergyIncr $DynAlzTol 2000 0
algorithm ModifiedNewton
}
if {$ok != 0} {
puts "Trying Broyden .."
algorithm Broyden 8
set ok [analyze 1 $DtAnalysis]
algorithm ModifiedNewton
}
if {$ok != 0} {
puts "Trying NewtonWithLineSearch .."
algorithm NewtonLineSearch .8
set ok [analyze 1 $DtAnalysis]
algorithm ModifiedNewton
}
}
}; # end if ok !0

puts "Ground Motion Done. End Time: [getTime]"

[mmasoudi]

the model converged by increasing number of iterations to undertake to satisfy element compatibility.
Is it a proper way to achieve convergence?

[vesna]

Yes, that is absolutely fine.