Pull Request Proposal for OpenSeesPyDoc

pyExamples/nonlinear_mdof.py

@@ -0,0 +1,130 @@
+'''
+=========================================================================================================================
+================================================== NonLinear MDOF =======================================================
+=========================================================================================================================
+
+By M. Eng. Joseph Jaramillo, National University of Engineering
+e-mail: [email protected]
+Date - 22/10/2024
+
+This example models a multi-degree-of-freedom (MDOF) damped system commonly used in earthquake engineering of a three-story 
+building. It conducts a nonlinear dynamic (time history)  analysis using the El Centro 1940 earthquake as the input ground
+motion.
+'''
+
+import openseespy.opensees as ops
+import numpy as np
+import plot_nonlinear_mdof_graphs
+# Base units
+m = 1                        # Meters
+s = 1                        # Seconds
+kN = 1                       # Kilo Newtons
+
+# Derivated units
+g = 9.81*m/s**2              # Gravity
+cm = 1e-2*m                  # Centimeter
+mm = 1e-3*m                  # Milimeter
+Ton = kN*s**2/m              # Ton
+
+# Parameters - data
+N = 3                        # N° DOF
+h = 0.05                     # Damping ratio
+dt = 0.02                    # Time step
+dt_out = 0.001               # Output time step
+tFinal = 35                  # Analysis stop time
+m1 = 0.1*Ton                 # Mass / floor 
+m2 = 0.1*Ton
+m3 = 0.1*Ton 
+Py1 = 0.55*kN                 # Yielding strength / floor
+Py2 = 0.45*kN                 
+Py3 = 0.30*kN                 
+K1 = 60*kN/m                 # Stiffness / floor
+K2 = 50*kN/m
+K3 = 30*kN/m
+b = 0.01                     # Strain-hardening ratio
+
+######## Model ###############
+ops.wipe()					                     # clear memory of all past model definitions
+ops.model('basic', '-ndm', 1, '-ndf', 1) 		 # Define the model builder, ndm=#dimension, ndf=#dofs
+
+# Create nodes
+ops.node(0, 0)
+ops.node(1, 0, '-mass', m1)
+ops.node(2, 0, '-mass', m2)
+ops.node(3, 0, '-mass', m3)
+
+# Define boundary condition
+ops.fix(0, 1)
+
+# Material definition
+ops.uniaxialMaterial('Steel01', 1, Py1, K1, b)
+ops.uniaxialMaterial('Steel01', 2, Py2, K2, b)
+ops.uniaxialMaterial('Steel01', 3, Py3, K3, b)
+
+# Element definition
+ops.element('zeroLength', 1, 0, 1, '-mat',    1   , '-dir', 1, '-doRayleigh', 1)
+ops.element('zeroLength', 2, 1, 2, '-mat',    2   , '-dir', 1, '-doRayleigh', 1)
+ops.element('zeroLength', 3, 2, 3, '-mat',    3   , '-dir', 1, '-doRayleigh', 1)
+
+# Set Rayleigh damping
+w1, w2, w3 = np.array(ops.eigen('-fullGenLapack', 3))**0.5
+a0 = 2*h*w1*w2/(w1+w2)
+a1 = 2*h/(w1+w2)
+ops.rayleigh(a0, .0, .0, a1)                     # RAYLEIGH damping
+
+# Natural periods
+print('\nNatural periods:       Natural frequencies:')
+print(f'T1 = {2*np.pi/w1:.3f} [s]    |    f1 = {w1/(2*np.pi):.3f} [Hz]')
+print(f'T2 = {2*np.pi/w2:.3f} [s]    |    f2 = {w2/(2*np.pi):.3f} [Hz]')
+print(f'T3 = {2*np.pi/w3:.3f} [s]    |    f3 = {w3/(2*np.pi):.3f} [Hz]')
+
+# Define the dynamic analysis 
+load_tag = 1
+patter_tag = 1
+direc = 1 
+ops.timeSeries('Path', load_tag, '-dt', dt, '-filePath', r'./el_centro.th', '-factor', g)                                           # Reading ground motion
+ops.pattern('UniformExcitation', patter_tag, direc, '-accel', load_tag)
+
+# Define output data files
+rD_path =  r'./Relative_disp.out'
+ops.recorder('Node', '-file', r'./Relative_disp.out', '-time', '-dT', dt_out, '-node', 1, 2, 3, '-dof', 1, 'disp')                  # Relative displacements with respect to the ground
+
+rA_path =  r'./Relative_accel.out'
+ops.recorder('Node', '-file', rA_path, '-time', '-dT', dt_out, '-node', 1, 2, 3, '-dof', 1, 'accel')                                # Relative accelerations with respect to the ground
+
+aA_path =  r'./Absolute_accel.out'
+ops.recorder('Node', '-file', aA_path, '-timeSeries', load_tag, '-time', '-dT', dt_out, '-node', 0, 1, 2, 3, '-dof', 1, 'accel')    # Absolute accelerations
+
+eF_path =  r'./Element_force.out'
+ops.recorder('Element', '-file', eF_path, '-time', '-dT', dt_out, '-ele', 1, 2, 3, 'localForce')                                    # Local spring force
+
+# Run the dynamic analysis
+Gamma = 0.5
+Beta = 0.25
+tol = 1.0e-12
+itrs = 100
+ops.wipeAnalysis()
+ops.algorithm('Newton')
+ops.system('BandGen')
+ops.numberer('Plain')
+ops.constraints('Plain')
+ops.integrator('Newmark', Gamma, Beta)
+ops.analysis('Transient')
+ops.test('NormUnbalance', tol, itrs)
+num_steps = int(tFinal/dt_out+1)
+ops.analyze(num_steps, dt_out)
+ops.wipe()
+
+rD = np.genfromtxt(rD_path, usecols=[1, 2, 3]).T
+
+rA = np.genfromtxt(rA_path, usecols=[1, 2, 3]).T
+aA = np.genfromtxt(aA_path, usecols=[1, 2, 3, 4]).T
+
+eF = np.genfromtxt(eF_path, usecols=[1, 2, 3]).T
+
+
+# Matplotlib plots 
+rD /= mm
+K = [K1, K2, K3]
+Py = [Py1, Py2, Py3]
+plot_nonlinear_mdof_graphs.plot(N, rD, rA, aA, eF, dt_out, num_steps, tFinal, K, Py)

pyExamples/plot_nonlinear_mdof_graphs.py

@@ -0,0 +1,124 @@
+import numpy as np
+import matplotlib.pyplot as plt
+from matplotlib.ticker import MultipleLocator
+from matplotlib import rcParams
+
+rcParams['font.family'] = 'Times New Roman'
+rcParams['font.size'] = 10
+rcParams['lines.linewidth'] = 0.8
+rcParams['lines.color'] = 'C0'
+rcParams['legend.frameon'] = False
+rcParams['axes.linewidth'] = 0.5
+rcParams['xtick.major.width'] = 0.9
+rcParams['xtick.minor.width'] = 0.9
+rcParams['ytick.major.width'] = 0.9
+rcParams['ytick.minor.width'] = 0.9
+rcParams["savefig.format"] = 'jpg'
+rcParams["savefig.transparent"] = False
+rcParams["grid.alpha"] = 0.35
+rcParams["grid.color"] = 'k'
+rcParams["grid.linewidth"] = 0.5
+rcParams['figure.titleweight'] = 'bold'
+rcParams['savefig.dpi'] = 600
+
+
+def plot(N, rD, rA, aA, eF, dt_out, num_steps, tFinal, K, Py):
+
+    z = lambda x: np.min(x) if abs(np.min(x)) > abs(np.max(x)) else np.max(x)
+    t = np.arange(0, num_steps) * dt_out
+
+    # Plotting absolute accelerations
+    fig, axs = plt.subplots(N+1, figsize=(8, 5), sharex=True, sharey=True)
+    fig.suptitle('Absolute Accelerations')
+
+    gm = aA[0]
+    peak = z(gm)
+    axs[0].plot(t, gm, f'k', label=f'Ground   peak: {peak:.2f} [m/s2]')
+
+    for i, ax in enumerate(axs[1:]):
+        x = aA[i+1]
+        peak = z(x)
+        ax.plot(t, x, f'C{i+1}', label=f'Floor {i+1}   peak: {peak:.2f} [m/s2]')
+
+    for ax in axs:
+        ax.set_xlim(.0, tFinal)
+        ax.set_ylim(-8, 8)
+        ax.grid()
+        ax.legend(loc=(0.72, 1.0))
+        ax.label_outer()
+        ax.xaxis.set_major_locator(MultipleLocator(5))
+        ax.yaxis.set_major_locator(MultipleLocator(4))
+
+    axs[1].set_ylabel('Acceleration [m/s2]', loc='top')
+    axs[1].yaxis.set_label_coords(-0.05, .5)
+    axs[3].set_xlabel('Time [s]')
+
+    plt.subplots_adjust(0.08, 0.1, 0.97, 0.9, 0.1, 0.4)
+    plt.savefig(r'./nonlinear_mdof_abs_accel.jpg')
+
+
+    # Plotting relative accelerations
+    fig, axs = plt.subplots(N, figsize=(8, 3.8))
+    fig.suptitle('Relative Accelerations')
+
+    for i, ax in enumerate(axs):
+        peak = z(rA[i])
+        ax.plot(t, rA[i], f'C{i+1}', label=f'Floor {i+1}   peak: {peak:.2f} [m/s2]')
+        ax.set_xlim(.0, tFinal)
+        ax.set_ylim(-8, 8)
+        ax.xaxis.set_major_locator(MultipleLocator(5))
+        ax.yaxis.set_major_locator(MultipleLocator(4))
+        ax.grid()
+        ax.legend(loc=(0.72, 1.0))
+        ax.label_outer()
+    axs[1].set_ylabel('Acceleration [m/s2]')
+    axs[2].set_xlabel('Time [s]')
+    plt.subplots_adjust(0.08, 0.12, 0.97, 0.9, 0.1, 0.4)
+    plt.savefig(r'./nonlinear_mdof_rel_accel.jpg')
+
+    # Plotting relative displacements
+    fig, axs = plt.subplots(N, figsize=(8, 3.8))
+    fig.suptitle('Relative Displacements')
+
+    for i, ax in enumerate(axs):
+        peak = z(rD[i])
+        ax.plot(t, rD[i], f'C{i+1}', label=f'Floor {i+1}   peak: {peak:.2f} [mm]')
+        ax.set_xlim(.0, tFinal)
+        ax.set_ylim(-80, 80)
+        ax.xaxis.set_major_locator(MultipleLocator(5))
+        ax.yaxis.set_major_locator(MultipleLocator(40))
+        ax.grid()
+        ax.legend(loc=(0.72, 1.0))
+        ax.label_outer()
+    axs[1].set_ylabel('Displacement [mm]')
+    axs[2].set_xlabel('Time [s]')
+    plt.subplots_adjust(0.08, 0.12, 0.97, 0.9, 0.1, 0.4)
+    plt.savefig(r'./nonlinear_mdof_rel_disp.jpg')
+
+
+    # Plotting Hysteresis by floor
+    fig, axs = plt.subplots(1, N, figsize=(10, 3.5))
+    fig.suptitle('Floor Hysteresis')
+
+    axs[0].plot(rD[0,:], eF[0,:], 'C1', label=f'Floor 1 (K={K[0]} [kN/m]; Py={Py[0]:.2f} [kN])')
+    axs[1].plot(rD[1,:] - rD[0,:], eF[1,:], 'C2', label=f'Floor 2 (K={K[1]} [kN/m]; Py={Py[1]:.2f} [kN])')
+    axs[2].plot(rD[2,:] - rD[1,:], eF[2,:], 'C3', label=f'Floor 3 (K={K[2]} [kN/m]; Py={Py[2]:.2f} [kN])')
+
+    for ax in axs:
+        # ax.set_xlim(.0, tFinal)
+        ax.set_ylim(-0.6, 0.6)
+        ax.set_xlim(-30, 30)
+        ax.grid()
+        ax.legend(loc=(0.06, 1.0))
+        ax.label_outer()
+        ax.set_xlabel('Displacement [mm]')
+        ax.xaxis.set_major_locator(MultipleLocator(10))
+        ax.yaxis.set_major_locator(MultipleLocator(0.2))
+
+    axs[0].set_ylabel('Force [kN]')
+    plt.subplots_adjust(0.06, 0.13, 0.97, 0.86, 0.1)
+    plt.savefig(r'./nonlinear_mdof_hysteresis.jpg')
+
+    plt.show()
+
+

src/NonlinearMDOF.rst

@@ -0,0 +1,31 @@
+.. include:: sub.txt
+
+==============
+Nonlinear MDOF
+==============
+
+#. The source code is developed by `Joseph Jaramillo <https://github.com/jjaramillod93>`_ from National University of Engineering.
+#. The source code is shown below, which can be downloaded :download:`here </pyExamples/nonlinear_mdof.py>`.
+#. The ground motion file :download:`here </pyExamples/el_centro.th>`.
+#. And also the python plotting script :download:`here </pyExamples/plot_nonlinear_mdof_graphs.py>`.
+#. Make sure the `numpy`_ and `matplotlib`_ packages are installed in your Python distribution.
+#. Run the source code in your favorite Python program and should see:
+
+::
+
+   Natural periods:       Natural frequencies:
+   T1 = 0.628 [s]    |    f1 = 1.592 [Hz]
+   T2 = 0.257 [s]    |    f2 = 3.898 [Hz]
+   T3 = 0.162 [s]    |    f3 = 6.164 [Hz]
+
+.. image:: /_static/nonlinear_mdof_abs_accel.jpg
+
+.. image:: /_static/nonlinear_mdof_rel_accel.jpg
+
+.. image:: /_static/nonlinear_mdof_rel_disp.jpg
+
+.. image:: /_static/nonlinear_mdof_hysteresis.jpg
+
+.. literalinclude:: /pyExamples/nonlinear_mdof.py
+   :linenos:
+

src/earthquake.rst

@@ -1,21 +1,23 @@
-=====================
- Earthquake Examples
-=====================
-
-
-#. :doc:`Canti2DEQ`
-#. :doc:`RCFrameEarthquake`
-#. :doc:`exampleNamedSpacedNonlinearSDOF`
-#. :doc:`exampleRotDSpectra`
-
-
-
-.. toctree::
-   :maxdepth: 1
-   :hidden:
-
-   Canti2DEQ
-   RCFrameEarthquake
-   exampleNamedSpacedNonlinearSDOF
-   exampleRotDSpectra
-
+=====================
+ Earthquake Examples
+=====================
+
+
+#. :doc:`Canti2DEQ`
+#. :doc:`RCFrameEarthquake`
+#. :doc:`exampleNamedSpacedNonlinearSDOF`
+#. :doc:`exampleRotDSpectra`
+#. :doc:`NonlinearMDOF`
+
+
+
+.. toctree::
+   :maxdepth: 1
+   :hidden:
+
+   Canti2DEQ
+   RCFrameEarthquake
+   exampleNamedSpacedNonlinearSDOF
+   exampleRotDSpectra
+   NonlinearMDOF
+