Soil–Structure Interaction

This example originates from the GitHub repository maintained by Professor Quan Gu of Xiamen University. Refer to: OpenSeesXMU

[1]:
import openseespy.opensees as ops
import opstool as opst
import opstool.vis.pyvista as opsvis
import matplotlib.pyplot as plt

The TCL script has been translated into a Python script using opstool.pre.tcl2py. For details, see: model.py

[2]:
from model import Model

Model()

Visualize the model

Here, notebook=True and jupyter_backend="jupyterlab" are used solely for the purpose of documentation generation. You can disable them during actual use.

[3]:
opsvis.set_plot_props(notebook=True)  # notebook=False should be used
fig = opst.vis.pyvista.plot_model()
fig.show(jupyter_backend="jupyterlab")  # fig.show() should be used for notebook=False
OPSTOOL ::  Model data has been saved to _OPSTOOL_ODB/ModelData-None.nc!
../../../_images/examples_post_SSI-GUQUAN_ssi_7_1.png

Gravity analysis

[4]:
ops.constraints("Transformation")
ops.numberer("RCM")
ops.test("NormDispIncr", 1e-06, 25, 2)
ops.integrator("LoadControl", 1, 1, 1, 1)
ops.algorithm("Newton")
ops.system("BandGeneral")
ops.analysis("Static")
ops.analyze(3)
print("soil gravity nonlinear analysis completed ...")
soil gravity nonlinear analysis completed ...
CTestNormDispIncr::test() - iteration: 6 current Norm: 2.68477e-11 (max: 1e-06, Norm deltaR: 6.07089e-11)
CTestNormDispIncr::test() - iteration: 1 current Norm: 1.49834e-16 (max: 1e-06, Norm deltaR: 6.91778e-11)
CTestNormDispIncr::test() - iteration: 1 current Norm: 1.47098e-16 (max: 1e-06, Norm deltaR: 8.07018e-11)

Earthquake analysis

[5]:
ops.timeSeries("Path", 1, "-dt", 0.01, "-filePath", "elcentro.txt", "-factor", 3)
ops.pattern("UniformExcitation", 1, 1, "-accel", 1)
[6]:
ops.wipeAnalysis()
ops.constraints("Transformation")
ops.test("NormDispIncr", 1e-06, 25)
ops.algorithm("Newton")
ops.numberer("RCM")
ops.system("BandGeneral")
ops.integrator("Newmark", 0.55, 0.275625)
ops.analysis("Transient")
[7]:
ODB = opst.post.CreateODB(odb_tag=1)  # Create ODB object

for _ in range(2400):
    ops.analyze(1, 0.005)
    ODB.fetch_response_step()  # Fetch response for the current step
ODB.save_response(zlib=True)  # Save response
OPSTOOL ::  All responses data with _odb_tag = 1 saved in _OPSTOOL_ODB/RespStepData-1.nc!

Post-processing

Frame Element Response

[8]:
FrameResp = opst.post.get_element_responses(odb_tag=1, ele_type="Frame")
OPSTOOL ::  Loading Frame response data from _OPSTOOL_ODB/RespStepData-1.nc ...
[9]:
f = FrameResp["sectionForces"].sel(eleTags=7, secDofs="MZ", secPoints=1)
d = FrameResp["sectionDeformations"].sel(eleTags=7, secDofs="MZ", secPoints=1)

plt.plot(d, f)
plt.show()
../../../_images/examples_post_SSI-GUQUAN_ssi_17_0.png

Nodal response

[10]:
NodalResp = opst.post.get_nodal_responses(odb_tag=1)
OPSTOOL ::  Loading all response data from _OPSTOOL_ODB/RespStepData-1.nc ...
[11]:
time = NodalResp.time
disp = NodalResp["disp"].sel(nodeTags=1, DOFs="UX")

plt.plot(time, disp)
plt.show()
../../../_images/examples_post_SSI-GUQUAN_ssi_20_0.png

Plane Soil Element

[12]:
PlaneResp = opst.post.get_element_responses(odb_tag=1, ele_type="Plane")
OPSTOOL ::  Loading Plane response data from _OPSTOOL_ODB/RespStepData-1.nc ...
[13]:
s = PlaneResp["Stresses"].sel(eleTags=37, GaussPoints=1, stressDOFs="sigma12")
e = PlaneResp["Strains"].sel(eleTags=37, GaussPoints=1, strainDOFs="eps12")

plt.plot(e, s, c="blue")
plt.xlabel("eps12")
plt.ylabel("sigma12")
plt.title("Stress-Strain Curve for Plane Element 37 at Gauss Point 1")
plt.show()
../../../_images/examples_post_SSI-GUQUAN_ssi_23_0.png
[ ]: