Code: Select all
# -*- coding: utf-8 -*-
"""
Created on Thu Jul 27 17:22:45 2023
@author: praf_
"""
# height is reguired to calculate drift from displacement. This was for experimental specimen
DriftAmp = [0.25,-0.25,0.25,-0.25
,0.5,-0.5,0.5,-0.5,
0.75,-0.75,0.75,-0.75,
1,-1,1,-1,
1.5,-1.5,1.5,-1.5,
2,-2,2,-2,
2.5,-2.5,2.5,-2.5,
3,-3, 3,-3,
4,-4,5,-0.1]
RunCyclic2Converge(controlNode, DriftAmp, fsheight, controlDOF=1)
# ---------RunCyclic2Converge.tcl----------------------------------------------------------------
# Changed into a procedure from Opensees Web-site. Slightly modified
# to stop analyzing if the loadfactor (pseudo-time) falls below zero.
# Also, it accepts an optional input, "dof" which defaults to 1.
# this is the dof in the POSITIVE OR NEGATIVE direction of which we push.
# Modified from Vam
import openseespy.opensees as ops
import numpy as np
def RunCyclic2Converge(controlNode, Drift, height, controlDOF=1):
# Set some parameters
# Arrays for plotting
currentDisp = ops.nodeDisp(controlNode, controlDOF)
arr=np.array(Drift) # converting to array to multipy to calculate Displacement from Drift
arr =height/100*arr # Displacement-Drift*height/100 max for each cycle
Disp = arr.tolist()
# ------------------------------------------------------------------------------
# Main analysis - performed as a number of analyses equal to the number of steps
# ------------------------------------------------------------------------------
ops.system("BandGeneral") # Overkill, but may need the pivoting!
ops.test("NormDispIncr", 1.0e-7, 10, 0)
ops.numberer("Plain")
ops.constraints("Penalty",1.0e14, 1.0e14)
ops.algorithm("Newton")
ops.analysis("Static")
for Dmax in Disp:
if Dmax > 0:
Dincr = 0.001
elif Dmax < 0:
Dincr =-0.001
ops.integrator("DisplacementControl", controlNode, controlDOF, Dincr)
ops.analysis("Static")
if Dincr > 0:
ok = 0
while (ok == 0 and currentDisp < Dmax):
ok = ops.analyze(1)
if ok != 0:
print( " ")
print ("Trying Newton with Initial Tangent ..")
ops.test("NormDispIncr", 1.0e-6, 2000, 0)
ops.algorithm('Newton', '-initial')
ok = ops.analyze(1)
ops.test("NormDispIncr", 1.0e-6, 6, 2)
ops.algorithm('Newton')
if ok != 0:
print( " ")
print ("Trying Broyden ..")
ops.test("NormDispIncr", 1.0e-6, 2000, 0)
ops.algorithm('Broyden', 8)
ok = ops.analyze(1)
ops.algorithm('Newton')
if ok != 0:
print( " ")
print ("Trying NewtonWithLineSearch ..")
ops.algorithm('NewtonLineSearch', 0.8)
ok = ops.analyze(1)
ops.algorithm('Newton')
currentDisp = ops.nodeDisp(controlNode, controlDOF)
factorLoad =ops.getTime()
#outputs["Load"].append(factorLoad)
#outputs["Displacement"].append(currentDisp)
#outputs["Drift"].append(currentDisp/height*100)
elif Dincr < 0 :
ok = 0
while (ok == 0 and currentDisp > Dmax):
ok = ops.analyze(1)
if ok != 0:
print( " ")
print ("Trying Newton with Initial Tangent ..")
ops.test("NormDispIncr", 1.0e-6, 2000, 0)
ops.algorithm('Newton', '-initial')
ok = ops.analyze(1)
ops.test("NormDispIncr", 1.0e-6, 6, 2)
ops.algorithm('Newton')
if ok != 0:
print( " ")
print ("Trying Broyden ..")
ops.test("NormDispIncr", 1.0e-6, 2000, 0)
ops.algorithm('Broyden', 8)
ok = ops.analyze(1)
ops.algorithm('Newton')
if ok != 0:
print( " ")
print ("Trying NewtonWithLineSearch ..")
ops.algorithm('NewtonLineSearch', 0.8)
ok = ops.analyze(1)
ops.algorithm('Newton')
currentDisp = ops.nodeDisp(controlNode, controlDOF)
factorLoad =ops.getTime()
# outputs["Load"].append(factorLoad)
# outputs["Displacement"].append(currentDisp)
# outputs["Drift"].append(currentDisp/height*100)
errorFileName = "ncerror.out"
# output success/failure results to a file
if errorFileName != "":
with open(errorFileName, "w") as errorFileID:
errorFileID.write(str(ok))
if ok != 0:
print("DispControl Analysis FAILED")
else:
print("DispControl Analysis SUCCESSFUL")