error 3X3 BRB steel frames

[xiaoshu]

Hello，everyone
I am pretty new to Opensees and am trying to model a 3X3 steel brace frames,but I have encountered some difficulties;It can't always converge;I have tried many methods to dubug,but I can't get the reasons;
Would you please help me with this? I do appreciate your help.I have also attached the code here:Thanks very much.
set NStories 3; # number of stories
set NBays 3; # number of bays
set WBay [expr 7.3*$m]; # width 0f bay
set Hstory1 [expr 5.5*$m]; # 1st story height
set HstoryTyp [expr 3.7*$m];#story height of other stories
set Hbuilding [expr ($NBays-1)*$HstoryTyp+$Hstory1]; #height of building;
# display data of the buidling
puts "Nstories=$NStories "
puts "NBays=$NBays";
puts "WBay=$WBay";
puts "Hstory1 = $Hstory1 m";
puts "HstoryTyp = $HstoryTyp m";
puts "Hbuilding = $Hbuilding m";
# calculate locations of beam/column joints:
set Pier1 0; # leftmost column line
set Pier2 [expr $Pier1+$WBay];
set Pier3 [expr $Pier2+$WBay];
set Pier4 [expr $Pier3+$WBay];
set Pier5 [expr $Pier4+$WBay]; #P-Delta column line
set Floor1 0; #ground floor
set Floor2 [expr $Floor1+$Hstory1];
set Floor3 [expr $Floor2+$HstoryTyp];
set Floor4 [expr $Floor3+$HstoryTyp];


# culculate nodal mass
set FlooroneWeight [expr 1538*$kN];
set FloortwoWeight [expr 1538*$kN];
set FloorthreeWeight [expr 1538*$kN];
set Wbuilding [expr $FlooroneWeight+$FloortwoWeight+$FloorthreeWeight]; #total building weight
set NodalMass2 [expr $FlooroneWeight/$g/4]; #weight of 2 floor
set NodalMass3 [expr $FloortwoWeight/$g/4]; #weight of 3 floor
set NodalMass4 [expr $FloorthreeWeight/$g/4]; #weight of 4 floor
set Negligible 1e-9; # a very small number to avoid problems with zero
puts "Wbuilding = $Wbuilding kg"

# difine some information for cross-section
set d1499 [expr 360*$mm];
set d1468 [expr 357*$mm];
set d1448 [expr 350*$mm];
set d1840 [expr 455*$mm];

# define nodes and assign masses to beam-column intersections of frame
# command: node nodeID xcoord ycoord -mass mass_dof1 mass_dof2 mass_dof3
# nodeID convention: "xy" where x = Pier # and y = Floor #
node 11 $Pier1 $Floor1;
node 21 $Pier2 $Floor1;
node 31 $Pier3 $Floor1;
node 41 $Pier4 $Floor1;
node 51 $Pier5 $Floor1;
node 12 $Pier1 $Floor2 -mass $NodalMass2 $Negligible $Negligible;
node 22 $Pier2 $Floor2 -mass $NodalMass2 $Negligible $Negligible;
node 32 $Pier3 $Floor2 -mass $NodalMass2 $Negligible $Negligible;
node 42 $Pier4 $Floor2 -mass $NodalMass2 $Negligible $Negligible;
node 52 $Pier5 $Floor2;

node 13 $Pier1 $Floor3 -mass $NodalMass3 $Negligible $Negligible;
node 23 $Pier2 $Floor3 -mass $NodalMass3 $Negligible $Negligible;
node 33 $Pier3 $Floor3 -mass $NodalMass3 $Negligible $Negligible;
node 43 $Pier4 $Floor3 -mass $NodalMass3 $Negligible $Negligible;
node 53 $Pier5 $Floor3;

node 14 $Pier1 $Floor4 -mass $NodalMass4 $Negligible $Negligible;
node 24 $Pier2 $Floor4 -mass $NodalMass4 $Negligible $Negligible;
node 34 $Pier3 $Floor4 -mass $NodalMass4 $Negligible $Negligible;
node 44 $Pier4 $Floor4 -mass $NodalMass4 $Negligible $Negligible;
node 54 $Pier5 $Floor4;


# define the nodes which are located at the half of the beam
node 22032 [expr 3*$WBay/2] $Floor2;
node 23033 [expr 3*$WBay/2] $Floor3;
node 24034 [expr 3*$WBay/2] $Floor4;

# define the nodes to model the concentrated plasticity.
# the nodes for the rotation spring of the column
node 126 $Pier1 $Floor2;
node 127 $Pier1 $Floor2;
node 136 $Pier1 $Floor3;
node 137 $Pier1 $Floor3;
node 146 $Pier1 $Floor4;

node 226 $Pier2 $Floor2;
node 227 $Pier2 $Floor2;
node 236 $Pier2 $Floor3;
node 237 $Pier2 $Floor3;
node 246 $Pier2 $Floor4;

node 326 $Pier3 $Floor2;
node 327 $Pier3 $Floor2;
node 336 $Pier3 $Floor3;
node 337 $Pier3 $Floor3;
node 346 $Pier3 $Floor4;

node 426 $Pier4 $Floor2;
node 427 $Pier4 $Floor2;
node 436 $Pier4 $Floor3;
node 437 $Pier4 $Floor3;
node 446 $Pier4 $Floor4;

# the nodes for the rotation spring of the beam
node 123 $Pier1 $Floor2;
node 222 $Pier2 $Floor2;
node 223 $Pier2 $Floor2;
node 322 $Pier3 $Floor2;
node 323 $Pier3 $Floor2;
node 422 $Pier4 $Floor2;

node 133 $Pier1 $Floor3;
node 232 $Pier2 $Floor3;
node 233 $Pier2 $Floor3;
node 332 $Pier3 $Floor3;
node 333 $Pier3 $Floor3;
node 432 $Pier4 $Floor3;

node 143 $Pier1 $Floor4;
node 242 $Pier2 $Floor4;
node 243 $Pier2 $Floor4;
node 342 $Pier3 $Floor4;
node 343 $Pier3 $Floor4;
node 442 $Pier4 $Floor4;

# extra nodes for springs in p-delta columns
node 526 $Pier5 $Floor2; # zero-stiffness spring will be used on p-delta column
node 527 $Pier5 $Floor2; # zero-stiffness spring will be used on p-delta column
node 536 $Pier5 $Floor3; # zero-stiffness spring will be used on p-delta column
node 537 $Pier5 $Floor3; # zero-stiffness spring will be used on p-delta column
node 546 $Pier5 $Floor4; # zero-stiffness spring will be used on p-delta column

# extra nodes for changing the material
node 1 $Pier1 [expr $Hstory1+1.5*$m];
node 2 $Pier2 [expr $Hstory1+1.5*$m];
node 3 $Pier3 [expr $Hstory1+1.5*$m];
node 4 $Pier4 [expr $Hstory1+1.5*$m];


# constrain beam-column joints in a floor to have the same lateral displacement using the "equalDOF" command
# Floor2
equalDOF 12 22 1;
equalDOF 12 32 1;
equalDOF 12 42 1;
equalDOF 12 52 1;
equalDOF 12 22032 1;


# Floor3
equalDOF 13 23 1;
equalDOF 13 33 1;
equalDOF 13 43 1;
equalDOF 13 53 1;
equalDOF 13 23033 1;


# Floor4
equalDOF 14 24 1;
equalDOF 14 34 1;
equalDOF 14 44 1;
equalDOF 14 54 1;
equalDOF 14 24034 1;

# assign boundary condidtions
# command: fix nodeID dxFixity dyFixity rzFixity
# fixity values: 1 = constrained; 0 = unconstrained
# pin all the column to the base
fix 11 1 1 0;
fix 21 1 1 0;
fix 31 1 1 0;
fix 41 1 1 0;
fix 51 1 1 0; # P-delta column is pinned

####################################################################################################
# Define Properties and Material
####################################################################################################

# Difine relative parameters for structure
set McMy 1.05;
set FyBeam [expr 248*$Mpa];
set FyCol [expr 345*$Mpa];
set Es [expr 2e5*$Mpa];
set epyBeam [expr $FyBeam/$Es];
set epyCol [expr $FyCol/$Es];
set mat1499 1;
set mat1468 2;
set mat1448 3;
set mat1840 4;
set n 10;

# Difine section for columns
set W1499 1;
set W1468 2;
set W1448 3;
set W1840 4;

# for W1499
set lp1499 [expr $d1499];
set A1499 [expr 18774*$mm2];
set Iz1499 [expr 46246.23*$cm4];
set IzColMod [expr $Iz1499/$d1499/2];#section modulus;
set MyCol [expr $FyCol*$IzColMod];
set Iz1499Modi [expr $Iz1499*($n+1.0)/$n];
set Ks_Col [expr $n*6.0*$Es*$Iz1499Modi/$Hstory1];
set a_mem [expr ($n+1.0)*($MyCol*($McMy-1.0)) / ($Ks_Col*0.0645)];
set b_col [expr ($a_mem)/(1.0+$n*(1.0-$a_mem))];
uniaxialMaterial Bilin $mat1499 $Ks_Col $b_col $b_col [expr 1.05*$MyCol] [expr -1.05*$MyCol] 1.63 1.63 1.63 1.63 1 1 1 1 0.0645 0.0645 0.1666 0.1666 0.4 0.4 0.4 0.4 1 1;
# Wsection $W1499 $mat1499 $d1499 [expr 370*$mm] [expr 20*$mm] [expr 12*$mm] 12 2 16 4;

# for W1468
set lp1468 [expr $d1468];
set A1468 [expr 12903.2*$mm2];
set Iz1468 [expr 30089.72*$cm4];
set IzColMod [expr $Iz1468/$d1468/2];#section modulus;
set MyCol [expr $FyCol*$IzColMod];
set Iz1468Modi [expr $Iz1468*($n+1.0)/$n];
set Ks_Col [expr $n*6.0*$Es*$Iz1468Modi/$HstoryTyp];
set a_mem [expr ($n+1.0)*($MyCol*($McMy-1.0)) / ($Ks_Col*0.0645)];
set b_col [expr ($a_mem)/(1.0+$n*(1.0-$a_mem))];
uniaxialMaterial Bilin $mat1468 $Ks_Col $b_col $b_col [expr 1.05*$MyCol] [expr -1.05*$MyCol] 1.4948 1.4948 1.4948 1.4948 1 1 1 1 0.0631 0.0631 0.1877 0.1877 0.4 0.4 0.4 0.4 1 1;
# Wsection $W1468 $mat1468 $d1468 [expr 255*$mm] [expr 18*$mm] [expr 10*$mm] 12 2 10 3;

# for W1448
set lp1448 [expr $d1448];
set A1448 [expr 9096.7*$mm2];
set Iz1448 [expr 20188.83*$cm4];
set IzColMod [expr $Iz1448/$d1448/2];#section modulus;
set MyCol [expr $FyCol*$IzColMod];
set Iz1448Modi [expr $Iz1448*($n+1.0)/$n];
set Ks_Col [expr $n*6.0*$Es*$Iz1448Modi/$HstoryTyp];
set a_mem [expr ($n+1.0)*($MyCol*($McMy-1.0)) / ($Ks_Col*0.0645)];
set b_col [expr ($a_mem)/(1.0+$n*(1.0-$a_mem))];
uniaxialMaterial Bilin $mat1448 $Ks_Col $b_col $b_col [expr 1.05*$MyCol] [expr -1.05*$MyCol] 1.3329 1.3329 1.3329 1.3329 1 1 1 1 0.0546 0.0546 0.1838 0.1838 0.4 0.4 0.4 0.4 1 1;
# Wsection $W1448 $mat1448 $d1448 [expr 204*$mm] [expr 15*$mm] [expr 9*$mm] 12 2 8 2;

# for W1840
set lp1840 [expr $d1840];
set A1840 [expr 7612.8*$mm2];
set Iz1840 [expr 25497.51*$cm4];
set IzBeamMod [expr $Iz1840/$d1840/2];#section modulus;
set MyBeam [expr $FyBeam*$IzBeamMod];
set Iz1840Modi [expr $Iz1840*($n+1.0)/$n];
set Ks_Beam [expr $n*6.0*$Es*$Iz1840Modi/$WBay];
set a_mem [expr ($n+1.0)*($MyBeam*($McMy-1.0)) / ($Ks_Beam*0.0498)];
set b_beam [expr ($a_mem)/(1.0+$n*(1.0-$a_mem))];
uniaxialMaterial Bilin $mat1840 $Ks_Beam $b_beam $b_beam [expr 1.05*$MyBeam] [expr -1.05*$MyBeam] 0.9580 0.9580 0.9580 0.9580 1 1 1 1 0.0498 0.0498 0.1781 0.1781 0.4 0.4 0.4 0.4 1 1;
# Wsection $W1840 $mat1840 $d1840 [expr 153*$mm] [expr 13*$mm] [expr 8*$mm] 16 2 7 2;

##################################################################
# Define Rotational Springs for Plastic Hinges and Elements
##################################################################
# set up geometric transformations of element
geomTransf PDelta 1; #for columns
geomTransf Linear 2; #for beam

# command: element beamWithHinges $eleID $iNode $jNode $secTagI $Lpi $secTagJ $Lpj $E $A $I $transfID
# eleID convention: "1xy" where 1 = col, x = Pier #, y = Story #
#for Column element
# Columns Story 1
element elasticBeamColumn 111 11 126 $A1499 $Es $Iz1499 1;
element elasticBeamColumn 121 21 226 $A1499 $Es $Iz1499 1;
element elasticBeamColumn 131 31 326 $A1499 $Es $Iz1499 1;
element elasticBeamColumn 141 41 426 $A1468 $Es $Iz1468 1;
# Columns Story 2
element elasticBeamColumn 1 127 1 $A1499 $Es $Iz1499 1;
element elasticBeamColumn 2 1 136 $A1468 $Es $Iz1468 1;
element elasticBeamColumn 3 227 2 $A1499 $Es $Iz1499 1;
element elasticBeamColumn 4 2 236 $A1468 $Es $Iz1468 1;
element elasticBeamColumn 5 327 3 $A1499 $Es $Iz1499 1;
element elasticBeamColumn 6 3 336 $A1468 $Es $Iz1468 1;
element elasticBeamColumn 7 427 4 $A1468 $Es $Iz1468 1;
element elasticBeamColumn 8 4 436 $A1448 $Es $Iz1448 1;
# Columns Story 3
element elasticBeamColumn 113 137 146 $A1468 $Es $Iz1468 1;
element elasticBeamColumn 123 237 246 $A1468 $Es $Iz1468 1;
element elasticBeamColumn 133 337 346 $A1468 $Es $Iz1468 1;
element elasticBeamColumn 143 437 446 $A1448 $Es $Iz1448 1;

#for beam
# Beams Story 1
element elasticBeamColumn 9 123 222 $A1840 $Es $Iz1840 2;
element elasticBeamColumn 10 223 22032 $A1840 $Es $Iz1840 2;
element elasticBeamColumn 11 22032 322 $A1840 $Es $Iz1840 2;
element elasticBeamColumn 12 323 422 $A1840 $Es $Iz1840 2;

# Beams Story 2
element elasticBeamColumn 13 133 232 $A1840 $Es $Iz1840 2;
element elasticBeamColumn 14 233 23033 $A1840 $Es $Iz1840 2;
element elasticBeamColumn 15 23033 332 $A1840 $Es $Iz1840 2;
element elasticBeamColumn 16 333 432 $A1840 $Es $Iz1840 2;

#Beams Story 3
element elasticBeamColumn 17 143 242 $A1840 $Es $Iz1840 2;
element elasticBeamColumn 18 243 24034 $A1840 $Es $Iz1840 2;
element elasticBeamColumn 19 24034 342 $A1840 $Es $Iz1840 2;
element elasticBeamColumn 20 343 442 $A1840 $Es $Iz1840 2;

# define p-delta columns and rigid links
set GravcolumnsmatTag 10;
set Egrav $Es;
uniaxialMaterial Elastic $GravcolumnsmatTag $Egrav;

# rigid links
# command: element truss $eleID $iNode $jNode $A $materialID
# eleID convention: 6xy, 6 = truss link, x = Bay #, y = Floor #
element truss 642 42 52 [expr 5*$m2] $GravcolumnsmatTag;
element truss 643 43 53 [expr 5*$m2] $GravcolumnsmatTag;
element truss 644 44 54 [expr 5*$m2] $GravcolumnsmatTag;

# p-delta columns
# eleID convention: 7xy, 7 = p-delta columns, x = Pier #, y = Story #
set AGCol [expr 5100*$in2];
set IzGCol [expr 100000*$in4];
geomTransf PDelta 3;

element elasticBeamColumn 752 51 526 $AGCol $Es $IzGCol 3;
element elasticBeamColumn 753 527 536 $AGCol $Es $IzGCol 3;
element elasticBeamColumn 754 537 546 $AGCol $Es $IzGCol 3;

# define p-delta column spring: zero-stiffness elastic spring
#Spring ID: "5xya" where 5 = leaning column spring, x = Pier #, y = Story #, a = location in story
# "a" convention: 1 = bottom of story, 2 = top of story
rotLeaningCol 5512 52 526;
rotLeaningCol 5521 52 527;
rotLeaningCol 5522 53 536;
rotLeaningCol 5531 53 537;
rotLeaningCol 5532 54 546;
# create region for P-Delta column springs
region 1 -ele 5512 5521 5522 5531 5532;

# Simulate BRB
# Difine BRB Geometry Parameters
set BRB1 20;
set BRB2 21;
set BRB3 22;

# simulate BRB with Mennegato-Pinto model
uniaxialMaterial Steel02 40 [expr 345*$Mpa] $Es 0.003 20 0.925 0.15 0.0005 0.01 0.0005 0.01;

# Define some nodes for constructing the BRB
#story 1
node 15 [expr 8.019*$m] [expr 1.09*$m];
node 16 [expr 8.1296*$m] [expr 1.256*$m];
node 17 [expr 9.125*$m] [expr 2.75*$m];
node 18 [expr 10.1204*$m] [expr 4.244*$m];
node 19 [expr 10.231*$m] [expr 4.41*$m];

node 105 [expr 13.881*$m] [expr 1.09*$m];
node 106 [expr 13.7704*$m] [expr 1.256*$m];
node 107 [expr 12.775*$m] [expr 2.75*$m];
node 108 [expr 11.7796*$m] [expr 4.244*$m];
node 109 [expr 11.669*$m] [expr 4.41*$m];

#story 2
node 25 [expr 8.072*$m] [expr 6.282*$m];
node 26 [expr 8.1422*$m] [expr 6.3532*$m];
node 27 [expr 9.125*$m] [expr 7.35*$m];
node 28 [expr 10.1078*$m] [expr 8.3468*$m];
node 29 [expr 10.178*$m] [expr 8.418*$m];

node 205 [expr 13.828*$m] [expr 6.282*$m];
node 206 [expr 13.7578*$m] [expr 6.3532*$m];
node 207 [expr 12.775*$m] [expr 7.35*$m];
node 208 [expr 11.7922*$m] [expr 8.3468*$m];
node 209 [expr 11.722*$m] [expr 8.418*$m];

#story 3
node 35 [expr 8.072*$m] [expr 9.982*$m];
node 36 [expr 8.1422*$m] [expr 10.0532*$m];
node 37 [expr 9.125*$m] [expr 11.05*$m];
node 38 [expr 10.1078*$m] [expr 12.0468*$m];
node 39 [expr 10.178*$m] [expr 12.118*$m];

node 305 [expr 13.828*$m] [expr 9.982*$m];
node 306 [expr 13.7578*$m] [expr 10.0532*$m];
node 307 [expr 12.775*$m] [expr 11.05*$m];
node 308 [expr 11.7922*$m] [expr 12.0468*$m];
node 309 [expr 11.722*$m] [expr 12.118*$m];

node 224 $Pier2 $Floor2;
node 324 $Pier3 $Floor2;
node 234 $Pier2 $Floor3;
node 334 $Pier3 $Floor3;

# Difine some sections for BRB
section Fiber $BRB1 {
patch circ 40 30 2 0 0 [expr 1.5*$in] [expr 5.135*$in] 0 360
}

section Fiber $BRB2 {
patch circ 40 30 2 0 0 [expr 1.5*$in] [expr 4.027*$in] 0 360
}

section Fiber $BRB3 {
patch circ 40 30 2 0 0 [expr 1.5*$in] [expr 3.794*$in] 0 360
}

section fiberSec 23 {
patch quadr 40 6 6 [expr -20*$in] [expr 1.5*$in] [expr -20*$in] [expr -1.5*$in] [expr 20*$in] [expr -1.5*$in] [expr -20*$in] [expr 1.5*$in]
}

#Transformation
geomTransf Corotational 8;

#Difine model
set ARGD [expr 9000*pow($in,2)];
set IzRGD [expr 82300000*pow($in,4)];
set Er [expr 100*$Es]
#story 1
element elasticBeamColumn 101 21 15 $ARGD $Er $IzRGD 8;# ‚Rigid’ offset
element nonlinearBeamColumn 102 15 16 3 23 8;# Gusset plate
element nonlinearBeamColumn 103 16 17 5 $BRB1 8;# Brace
element nonlinearBeamColumn 104 17 18 5 $BRB1 8;# Brace
element nonlinearBeamColumn 105 18 19 3 23 8;# Gusset plate
element elasticBeamColumn 106 19 22032 $ARGD $Er $IzRGD 8;# ‚Rigid’ offset

element elasticBeamColumn 107 31 105 $ARGD $Er $IzRGD 8;# ‚Rigid’ offset
element nonlinearBeamColumn 108 105 106 3 23 8;# Gusset plate
element nonlinearBeamColumn 109 106 107 5 $BRB1 8;# Brace
element nonlinearBeamColumn 110 107 108 5 $BRB1 8;# Brace
element nonlinearBeamColumn 1110 108 109 3 23 8;# Gusset plate
element elasticBeamColumn 112 109 22032 $ARGD $Er $IzRGD 8;# ‚Rigid’ offset

#story 2
element elasticBeamColumn 1130 224 25 $ARGD $Er $IzRGD 8;# ‚Rigid’ offset
element nonlinearBeamColumn 114 25 26 3 23 8;# Gusset plate
element nonlinearBeamColumn 115 26 27 5 $BRB2 8;# Brace
element nonlinearBeamColumn 116 27 28 5 $BRB2 8;# Brace
element nonlinearBeamColumn 117 28 29 3 23 8;# Gusset plate
element elasticBeamColumn 118 29 23033 $ARGD $Er $IzRGD 8;# ‚Rigid’ offset

element elasticBeamColumn 119 324 205 $ARGD $Er $IzRGD 8;# ‚Rigid’ offset
element nonlinearBeamColumn 120 205 206 3 23 8;# Gusset plate
element nonlinearBeamColumn 1210 206 207 5 $BRB2 8;# Brace
element nonlinearBeamColumn 122 207 208 5 $BRB2 8;# Brace
element nonlinearBeamColumn 1230 208 209 3 23 8;# Gusset plate
element elasticBeamColumn 124 209 23033 $ARGD $Er $IzRGD 8;# ‚Rigid’ offset

#story 3
element elasticBeamColumn 125 234 35 $ARGD $Er $IzRGD 8;# ‚Rigid’ offset
element nonlinearBeamColumn 126 35 36 3 23 8;# Gusset plate
element nonlinearBeamColumn 127 36 37 5 $BRB2 8;# Brace
element nonlinearBeamColumn 128 37 38 5 $BRB2 8;# Brace
element nonlinearBeamColumn 129 38 39 3 23 8;# Gusset plate
element elasticBeamColumn 130 39 24034 $ARGD $Er $IzRGD 8;# ‚Rigid’ offset

element elasticBeamColumn 1310 334 305 $ARGD $Er $IzRGD 8;# ‚Rigid’ offset
element nonlinearBeamColumn 132 305 306 3 23 8;# Gusset plate
element nonlinearBeamColumn 1330 306 307 5 $BRB2 8;# Brace
element nonlinearBeamColumn 134 307 308 5 $BRB2 8;# Brace
element nonlinearBeamColumn 135 308 309 3 23 8;# Gusset plate
element elasticBeamColumn 136 309 24034 $ARGD $Er $IzRGD 8;# ‚Rigid’ offset

# define the rotation springs
# the rotation springs of beam
proc rotSpringbeam {eleID nodeR nodeC } {
#Spring Stiffness
set K 1;
uniaxialMaterial Elastic $eleID $K
element zeroLength $eleID $nodeR $nodeC -mat $eleID -dir 6
equalDOF $nodeR $nodeC 1 2
}


proc rotSpringcolumn {eleID nodeR nodeC tag} {
element zeroLength $eleID $nodeR $nodeC -mat $tag -dir 6
equalDOF $nodeR $nodeC 1 2
}

# Beams story 1
rotSpringbeam 27 12 123;
rotSpringbeam 28 22 222;
rotSpringbeam 29 22 223;
rotSpringbeam 30 32 322;
rotSpringbeam 31 32 323;
rotSpringbeam 32 42 422;

# Beams story 2
rotSpringbeam 33 13 133;
rotSpringbeam 34 23 232;
rotSpringbeam 35 23 233;
rotSpringbeam 36 33 332;
rotSpringbeam 37 33 333;
rotSpringbeam 38 43 432;

# Beams story 3
rotSpringbeam 39 14 143;
rotSpringbeam 400 24 242;
rotSpringbeam 41 24 243;
rotSpringbeam 42 34 342;
rotSpringbeam 43 34 343;
rotSpringbeam 44 44 442;

# for rotational springs at the end of the brb
rotSpringbeam 45 22 224;
rotSpringbeam 46 32 324;
rotSpringbeam 47 23 234;
rotSpringbeam 48 33 334;

# create region for frame beam springs
region 2 -eleRange 27 44;

# Columns floor2
rotSpringcolumn 49 12 126 $mat1499;
rotSpringcolumn 50 12 127 $mat1499;
rotSpringcolumn 51 22 226 $mat1499;
rotSpringcolumn 52 22 227 $mat1499;
rotSpringcolumn 53 32 326 $mat1499;
rotSpringcolumn 54 32 327 $mat1499;
rotSpringcolumn 55 42 426 $mat1468;
rotSpringcolumn 56 42 427 $mat1468;

# Columns floor3
rotSpringcolumn 57 13 136 $mat1468;
rotSpringcolumn 58 13 137 $mat1468;
rotSpringcolumn 59 23 236 $mat1468;
rotSpringcolumn 60 23 237 $mat1468;
rotSpringcolumn 61 33 336 $mat1468;
rotSpringcolumn 62 33 337 $mat1468;
rotSpringcolumn 63 43 436 $mat1448;
rotSpringcolumn 64 43 437 $mat1448;

# Columns floor4
rotSpringcolumn 65 14 146 $mat1468;
rotSpringcolumn 66 24 246 $mat1468;
rotSpringcolumn 67 34 346 $mat1468;
rotSpringcolumn 68 44 446 $mat1448;

# create region for frame column springs
region 3 -eleRange 49 68;

DisplayModel2D NodeNumbers;
############################################################################
# Eigenvalue Analysis
############################################################################

set nEigenI 1; # mode 1
set nEigenJ 3; # mode 3
set lambdaN [eigen -fullGenLapack [expr $nEigenJ]];
set lambdaI [lindex $lambdaN [expr $nEigenI-1]]; # eigenvalue mode i
set lambdaJ [lindex $lambdaN [expr $nEigenJ-1]]; # eigenvalue mode j
set omegaI [expr pow($lambdaI,0.5)];
set omegaJ [expr pow($lambdaJ,0.5)];
set T1 [expr 2*$pi/$omegaI];
set T3 [expr 2*$pi/$omegaJ];
puts "T1 = $T1 s , T2 = $T3 s"

############################################################################
# Gravity Loads & Gravity Analysis
############################################################################

# define Gravity
# Gravity Loads
pattern Plain 1 Linear {
load 12 0 [expr -50*$kN] 0
load 22 0 [expr -50*$kN] 0
load 32 0 [expr -50*$kN] 0
load 42 0 [expr -50*$kN] 0
load 13 0 [expr -50*$kN] 0
load 23 0 [expr -50*$kN] 0
load 33 0 [expr -50*$kN] 0
load 43 0 [expr -50*$kN] 0
load 14 0 [expr -50*$kN] 0
load 24 0 [expr -50*$kN] 0
load 34 0 [expr -50*$kN] 0
load 44 0 [expr -50*$kN] 0

load 52 0 [expr -1538*$kN] 0
load 53 0 [expr -1538*$kN] 0
load 54 0 [expr -1538*$kN] 0
}
source LibAnalysisStaticParameters.tcl;
loadConst -time 0.0;
puts "Model Built"
puts "get Static Done!"

############################################################################
# Recorders
############################################################################
#General Recorders
#for RBase
recorder Node -file 3X3_Steel_brace_framesout/dynamicsout/11r.out -time -node 11 -dof 1 2 3 reaction;
recorder Node -file 3X3_Steel_brace_framesout/dynamicsout/21r.out -time -node 21 -dof 1 2 3 reaction;
recorder Node -file 3X3_Steel_brace_framesout/dynamicsout/31r.out -time -node 31 -dof 1 2 3 reaction;
recorder Node -file 3X3_Steel_brace_framesout/dynamicsout/41r.out -time -node 41 -dof 1 2 3 reaction;

# for rdrift
recorder Drift -file 3X3_Steel_brace_framesout/dynamicsout/1d.out -time -iNode 11 -jNode 14 -dof 1 -perpDirn 2;
recorder Drift -file 3X3_Steel_brace_framesout/dynamicsout/2d.out -time -iNode 21 -jNode 22 -dof 1 -perpDirn 2;
recorder Drift -file 3X3_Steel_brace_framesout/dynamicsout/3d.out -time -iNode 31 -jNode 34 -dof 1 -perpDirn 2;
recorder Drift -file 3X3_Steel_brace_framesout/dynamicsout/4d.out -time -iNode 41 -jNode 44 -dof 1 -perpDirn 2;

# for rdisplacement
recorder Node -file 3X3_Steel_brace_framesout/dynamicsout/14D.out -time -node 14 -dof 1 disp;
recorder Node -file 3X3_Steel_brace_framesout/dynamicsout/24D.out -time -node 24 -dof 1 disp;
recorder Node -file 3X3_Steel_brace_framesout/dynamicsout/34D.out -time -node 34 -dof 1 disp;
recorder Node -file 3X3_Steel_brace_framesout/dynamicsout/44D.out -time -node 44 -dof 1 disp;

# inner drift
recorder Drift -file 3X3_Steel_brace_framesout/dynamicsout/101d.out -time -iNode 11 -jNode 12 -dof 1 -perpDirn 2;
recorder Drift -file 3X3_Steel_brace_framesout/dynamicsout/102d.out -time -iNode 12 -jNode 13 -dof 1 -perpDirn 2;
recorder Drift -file 3X3_Steel_brace_framesout/dynamicsout/103d.out -time -iNode 13 -jNode 14 -dof 1 -perpDirn 2;
recorder Drift -file 3X3_Steel_brace_framesout/dynamicsout/201d.out -time -iNode 21 -jNode 22 -dof 1 -perpDirn 2;
recorder Drift -file 3X3_Steel_brace_framesout/dynamicsout/202d.out -time -iNode 22 -jNode 23 -dof 1 -perpDirn 2;
recorder Drift -file 3X3_Steel_brace_framesout/dynamicsout/203d.out -time -iNode 23 -jNode 24 -dof 1 -perpDirn 2;
recorder Drift -file 3X3_Steel_brace_framesout/dynamicsout/301d.out -time -iNode 31 -jNode 32 -dof 1 -perpDirn 2;
recorder Drift -file 3X3_Steel_brace_framesout/dynamicsout/302d.out -time -iNode 32 -jNode 33 -dof 1 -perpDirn 2;
recorder Drift -file 3X3_Steel_brace_framesout/dynamicsout/303d.out -time -iNode 33 -jNode 34 -dof 1 -perpDirn 2;
recorder Drift -file 3X3_Steel_brace_framesout/dynamicsout/401d.out -time -iNode 41 -jNode 42 -dof 1 -perpDirn 2;
recorder Drift -file 3X3_Steel_brace_framesout/dynamicsout/402d.out -time -iNode 42 -jNode 43 -dof 1 -perpDirn 2;
recorder Drift -file 3X3_Steel_brace_framesout/dynamicsout/403d.out -time -iNode 43 -jNode 44 -dof 1 -perpDirn 2;

# for damper
recorder Element -file 3X3_Steel_brace_framesout/dynamicsout/ele102.out -time -ele 102 localForce;
recorder Element -file 3X3_Steel_brace_framesout/dynamicsout/ele108.out -time -ele 108 localForce;
recorder Element -file 3X3_Steel_brace_framesout/dynamicsout/ele114.out -time -ele 114 localForce;
recorder Element -file 3X3_Steel_brace_framesout/dynamicsout/ele120.out -time -ele 120 localForce;
recorder Element -file 3X3_Steel_brace_framesout/dynamicsout/ele126.out -time -ele 126 localForce;
recorder Element -file 3X3_Steel_brace_framesout/dynamicsout/ele132.out -time -ele 132 localForce;
recorder Element -file 3X3_Steel_brace_framesout/dynamicsout/ele102d.out -time -ele 102 deformation;
recorder Element -file 3X3_Steel_brace_framesout/dynamicsout/ele108d.out -time -ele 108 deformation;
recorder Element -file 3X3_Steel_brace_framesout/dynamicsout/ele114d.out -time -ele 114 deformation;
recorder Element -file 3X3_Steel_brace_framesout/dynamicsout/ele120d.out -time -ele 120 deformation;
recorder Element -file 3X3_Steel_brace_framesout/dynamicsout/ele126d.out -time -ele 126 deformation;
recorder Element -file 3X3_Steel_brace_framesout/dynamicsout/ele132d.out -time -ele 132 deformation;

#######################################################################################
# #
# Analysis Section #
# #
#######################################################################################
############################################################################
# Time History/Dynamic Analysis
############################################################################
puts "Running dynamic analysis..."
# display deformed shape:
set ViewScale 5; # amplify display of deformed shape
DisplayModel2D DeformedShape $ViewScale; # display deformed shape, the scaling factor needs to be adjusted for each model
set GMdirection 1;
set Scalefact 1;
set GMfact [expr $Scalefact*$g];
set dt [expr 0.02*$sec];
set DtAnalysis [expr 0.02*$sec];
set TmaxAnalysis [expr 53.76*$sec];
set xDamp 0.02; # damping ratio
set n 10
set a0 [expr $xDamp*2.0*$omegaI*$omegaJ/($omegaI + $omegaJ)]; # mass damping coefficient based on first and second modes
set a1 [expr $xDamp*2.0/($omegaI + $omegaJ)]; # stiffness damping coefficient based on first and second modes
set a1_mod [expr $a1*(1.0+$n)/$n]; # modified stiffness damping coefficient used for n modified elements.
## use "region" command when defining mass proportional damping so that the stiffness proportional damping isn't canceled
region 4 -ele 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 111 121 131 141 113 123 133 143 -rayleigh 0.0 0.0 $a1_mod 0.0; # assign stiffness proportional damping to frame beams & columns w/ n modifications
region 6 -node 12 22 32 42 13 23 33 43 14 24 34 44 -rayleigh $a0 0.0 0.0 0.0; # assign mass proportional damping to structure (assign to nodes with mass)
set AccelSeries "Series -dt $dt -filePath El.tcl -factor $GMfact";
pattern UniformExcitation 5 $GMdirection -accel $AccelSeries ;
source LibAnalysisDynamicParameters.tcl;
constraints Transformation;
numberer RCM;
system BandGeneral;
set TolDynamic 1e-6;
set maxNumIterDynamic 50;
test $testTypeDynamic $TolDynamic $maxNumIterDynamic
set algorithmTypeDynamic ModifiedNewton
algorithm $algorithmTypeDynamic;
integrator Newmark 0.5 0.25;
analysis Transient;
set Nsteps [expr int($TmaxAnalysis/$DtAnalysis)];
set ok [analyze $Nsteps $DtAnalysis]; # actually perform analysis; returns ok=0 if analysis was successful

if {$ok != 0} { ; # if 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 ok [analyze 1 $DtAnalysis]
set controlTime [getTime]
set ok [analyze 1 $DtAnalysis]
if {$ok != 0} {
puts "Trying Newton with Initial Tangent .."
test NormDispIncr $TolDynamic 100 0
algorithm Newton -initial
set ok [analyze 1 $DtAnalysis]
test $testTypeDynamic $TolDynamic $maxNumIterDynamic 0
algorithm $algorithmTypeDynamic
}
if {$ok != 0} {
puts "Trying Broyden .."
algorithm Broyden 8
set ok [analyze 1 $DtAnalysis]
algorithm $algorithmTypeDynamic
}
if {$ok != 0} {
puts "Trying NewtonWithLineSearch .."
algorithm NewtonLineSearch .8
set ok [analyze 1 $DtAnalysis]
algorithm $algorithmTypeDynamic
}
}
}; # end if ok !0



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

[xiaoshu]

WARNING: CTestEnergyIncr::test() - failed to converge
after: 50 iterations
current EnergyIncr: 2.12244e+007 (max: 1e-006) Norm deltaX: 16.9499, Norm deltaR: 2.66389e+006
ModifiedNewton::solveCurrentStep() -the ConvergenceTest object failed in test()
DirectIntegrationAnalysis::analyze() - the Algorithm failed at time 0.02
OpenSees > analyze failed, returned: -3 error flag
WARNING: CTestEnergyIncr::test() - failed to converge
after: 50 iterations
current EnergyIncr: 2.12244e+007 (max: 1e-006) Norm deltaX: 16.9499, Norm deltaR: 2.66389e+006
ModifiedNewton::solveCurrentStep() -the ConvergenceTest object failed in test()
DirectIntegrationAnalysis::analyze() - the Algorithm failed at time 0.02
OpenSees > analyze failed, returned: -3 error flag
WARNING: CTestEnergyIncr::test() - failed to converge
after: 50 iterations
current EnergyIncr: 2.12244e+007 (max: 1e-006) Norm deltaX: 16.9499, Norm deltaR: 2.66389e+006
ModifiedNewton::solveCurrentStep() -the ConvergenceTest object failed in test()
DirectIntegrationAnalysis::analyze() - the Algorithm failed at time 0.02
OpenSees > analyze failed, returned: -3 error flag
Trying Newton with Initial Tangent ..
WARNING: CTestNormDispIncr::test() - failed to converge
after: 5 iterations current Norm: 1.#INF (max: 1e-006, Norm deltaR: 1.#QNAN)
NewtnRaphson::solveCurrentStep() -the ConvergenceTest object failed in test()
DirectIntegrationAnalysis::analyze() - the Algorithm failed at time 0.02
OpenSees > analyze failed, returned: -3 error flag
Trying Broyden ..
WARNING - ForceBeamColumn2d::update - failed to get compatible element forces & deformations for element: 128(dW: << 2.43406e-008)
Domain::update - domain failed in update
Newmark::update() - failed to update the domain
WARNING Broyden::solveCurrentStep() -the Integrator failed in update()
DirectIntegrationAnalysis::analyze() - the Algorithm failed at time 0.02
OpenSees > analyze failed, returned: -3 error flag
Trying NewtonWithLineSearch ..
WARNING: CTestEnergyIncr::test() - failed to converge
after: 50 iterations
current EnergyIncr: 1.02676e+008 (max: 1e-006) Norm deltaX: 155.792, Norm deltaR: 1.42081e+006
NewtonLineSearch::solveCurrentStep() -the ConvergenceTest object failed in test()
DirectIntegrationAnalysis::analyze() - the Algorithm failed at time 0.02

~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

I don't even know what the error means. Would you please help me with this? I do appreciate your help.I have also attached the code here:
Thanks very much.

[Jhno]

WARNING: CTestEnergyIncr::test() - failed to converge
after: 50 iterations

This mean that your model failed to meet the test criteria.
The rest of the message says the same, only you used an "adaptative" solution so once once criteria fails, it try another one and so forth.

This is your test command : test $testTypeDynamic $TolDynamic $maxNumIterDynamic
this is your adaptative solution
if {$ok != 0} { ; # if 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 ok [analyze 1 $DtAnalysis]
set controlTime [getTime]
set ok [analyze 1 $DtAnalysis]
if {$ok != 0} {
puts "Trying Newton with Initial Tangent .."
test NormDispIncr $TolDynamic 100 0
algorithm Newton -initial
set ok [analyze 1 $DtAnalysis]
test $testTypeDynamic $TolDynamic $maxNumIterDynamic 0
algorithm $algorithmTypeDynamic
}
if {$ok != 0} {
puts "Trying Broyden .."
algorithm Broyden 8
set ok [analyze 1 $DtAnalysis]
algorithm $algorithmTypeDynamic
}
if {$ok != 0} {
puts "Trying NewtonWithLineSearch .."
algorithm NewtonLineSearch .8
set ok [analyze 1 $DtAnalysis]
algorithm $algorithmTypeDynamic
}
}
}; # end if ok !0

I guess for someone who's pretty new to OpenSees as you say, you should start with a rather simpler model. Check the simple example using elastic model on the OS wiki and try understand before you use a more complex model. Good Luck !

[xiaoshu]

thanks Jhno!This is a very good suggestion! the parameters I used to construct the IMK model are calculated according to the corresponding paper(I think that these parameters are not likely to be wrong.)! I guess my model has some incompatibility issues,how can I debug it? I just feel very clueless.I have tried to use a 1 story X 1 Bay steel frame to debug,but it didn't work unless I change the way to construct the leaning column.Can you tell me how to solve the incompatibility issues?Most sincere thanks!

[selimgunay]

You can try another material like Steel01 first to check your model and make sure that there are not issues

[xiaoshu]

thanks selimgunay! I have just fighting for the final exams last two weeks.I will try your advice!