G2 Example 2

G2Example2.py

from G2AnalysisSDOFNonlinearDynamic import *
from G2MaterialBilinear import *
import numpy as np
import matplotlib.pyplot as plt

# SI units: N, m, kg, sec.

# Natural periods considered
minPeriod = 0.1
maxPeriod = 5.0
numPeriodPoints = 50
naturalPeriods = np.linspace(minPeriod, maxPeriod, numPeriodPoints)

# Select natural period index for which to plot ductility demands
plotIndex = int(numPeriodPoints/5)

# Strengths considered
minRelativeStrength = 0.5
maxRelativeStrength = 1.0
numStrengthPoints = 6
sValues = np.linspace(minRelativeStrength, maxRelativeStrength, numStrengthPoints)

# Fix stiffness and let mass be calculated (elasto-plastic material)
E = 1e4
alpha = 0.0

# Damping ratio
dampingRatio = 0.05

# Load ground motion
groundMotion = "GroundMotionElCentro.txt"
dt = 0.02
gmScaling = 1
gm = []
f = open(groundMotion, "r")
lines = f.readlines()
for oneline in lines:
    splitline = oneline.split()
    for j in range(len(splitline)):
        value = 9.81 * float(splitline[j]) * gmScaling
        gm.append(value)
gm = np.array(gm)
duration = dt * (len(gm)-1)

# First do linear analysis to pick up ue for all natural periods
uy = 1e6
fy = E * uy
ueStorage = np.zeros(numPeriodPoints)
for i in range(numPeriodPoints):

    # Mass from natural period
    M = (naturalPeriods[i]/2/np.pi)**2 * E

    # Run analysis
    material = bilinearMaterial(['Bilinear', E, fy, alpha])
    t, uTrack, vTrack, aTrack, dudx, dvdx, dadx, dnl = nonlinearDynamicSDOFAnalysis(duration, dt, material, M, dampingRatio, gm, gmScaling, [])

    # Pick up maximum response
    ueStorage[i] = np.max(np.abs(uTrack))

# Loop over natural periods, and for each, loop over relative strengths and calculate ductility demand
ratioStorage = np.zeros((numPeriodPoints, numStrengthPoints))
for i in range(numPeriodPoints):

    # Mass from natural period
    M = (naturalPeriods[i]/2/np.pi)**2 * E

    # Store results if plots are requested
    if i == plotIndex:
        Tn = naturalPeriods[i]
        equalDisplacementDuctilityStorage = []
        ductilityDemandStorage = []

    # Loop over relative strengths
    for j in range(numStrengthPoints):

        # Yield displacement and yield strength
        uy = sValues[j]*ueStorage[i]
        fy = E * uy

        # Run analysis
        material = bilinearMaterial(['Bilinear', E, fy, alpha])
        t, uTrack, vTrack, aTrack, dudx, dvdx, dadx, dnl = nonlinearDynamicSDOFAnalysis(duration, dt, material, M, dampingRatio, gm, gmScaling, [])

        # Pick up maximum response
        um = np.max(np.abs(uTrack))

        # Ductility demands
        equalDisplacementDuctility = ueStorage[i]/uy
        ductilityDemand = um/uy

        # Append to storage
        if i == plotIndex:
            equalDisplacementDuctilityStorage.append(equalDisplacementDuctility)
            ductilityDemandStorage.append(ductilityDemand)

        ratioStorage[i, j] = ductilityDemand/equalDisplacementDuctility

    if i == plotIndex:

        # Plot "spectrum" as function of s
        plt.ion()
        plt.figure()
        plt.autoscale(True)
        plt.grid(True)
        plt.xlabel("Relative Strength")
        plt.ylabel("Ductility Demand")
        plt.title("Ductility Demand ($T_n=%.2fsec.$)" % Tn)
        plt.plot(sValues, ductilityDemandStorage, 'k-', linewidth=1.0, label="Actual ductility demand")
        plt.plot(sValues, equalDisplacementDuctilityStorage, 'k--', linewidth=1.0, label="Equal displacement rule")
        plt.legend(loc='upper right')
        print('\n'"Click somewhere in the plot to continue...")
        plt.waitforbuttonpress()

        # Another plot
        plt.ion()
        plt.figure()
        plt.autoscale(True)
        plt.grid(True)
        plt.xlabel("Equal Displacement Ductility")
        plt.ylabel("Actual Ductility Demand")
        plt.title("Ductility Demand ($T_n=%.2fsec.$)" % Tn)
        plt.plot(equalDisplacementDuctilityStorage, ductilityDemandStorage, 'k-', linewidth=1.0, label="Actual ductility demand")
        plt.plot(equalDisplacementDuctilityStorage, equalDisplacementDuctilityStorage, 'k--', linewidth=1.0, label="Equal displacement rule")
        plt.legend(loc='upper left')
        print('\n'"Click somewhere in the plot to continue...")
        plt.waitforbuttonpress()

# Plot spectrum
colorList = ["black", "blue", "orange", "green", "red", "purple", "brown", "pink", "gray", "olive", "cyan", "magenta", "yellow"]
plt.ion()
plt.figure()
plt.autoscale(True)
plt.grid(True)
plt.xlabel("Natural Period, $T_n$")
plt.ylabel("Ductility Ratio")
plt.title("Ductility Demand Spectrum")
for j in range(numStrengthPoints):
    plt.plot(naturalPeriods, ratioStorage[:, j], '-', color=colorList[j], linewidth=1.0, label=("s=%.1f" % (sValues[j])))
plt.legend(loc='upper right')
print('\n'"Click somewhere in the plot to continue...")
plt.waitforbuttonpress()