Arch Instability

6 min read • 1,092 words

Several nonlinear static analysis methods are used to investigate instabilities in a shallow arch.

The files for this example are:

arch.py - This file contains the function arch_model which is used construct the model
IncrementalAnalysis.ipynb - This is the current Jupyter notebook file

arch.py (collapsed)

import numpy as np
import opensees.openseespy

# Create the model
def arch_model():

    # Define model parameters
    L      = 5000
    Rise   = 500
    Offset = 200

    # Compute radius
    R  = Rise/2 + (2*L)**2/(8*Rise)
    th = 2*np.arcsin(L/R)

    #
    # Build the model
    #
    model = opensees.openseespy.Model(ndm=2, ndf=3)

    # Create nodes
    ne  = 10
    nen =  2             # nodes per element
    nn  = ne*(nen-1)+1
    mid = (nn+1)//2      # midpoint node

    for i, angle in enumerate(np.linspace(-th/2, th/2, nn)):
        tag = i + 1

        # Compute x and add offset if midpoint
        x = R*np.sin(angle)
        if tag == mid:
            x -= Offset

        # Compute y
        y = R*np.cos(angle)

        # create the node
        model.node(tag, x, y)

    # Define elastic section properties
    E = 200
    A = 1e4
    I = 1e8
    model.section("ElasticFrame", 1, A=A, E=E, I=I)

    # Create elements
    transf = 1
    model.geomTransf("Corotational", transf)
    for i in range(ne):
        tag   = i+1
        nodes = (i+1, i+2)
        model.element("PrismFrame", tag, nodes, section=1, transform=transf)


    model.fix( 1, 1, 1, 0)
    model.fix(nn, 1, 1, 0)

    # Create a load pattern that scales linearly
    model.pattern("Plain", 1, "Linear")

    # Add a nodal load to the pattern
    model.load(mid, 0.0, -1.0, 0.0, pattern=1)


    # Select linear solver, and ensure determinants are
    # computed and stored.
    # model.system("ProfileSPD")
    # model.system("FullGeneral")
    # model.system("BandGeneral")
    model.system("Umfpack", det=True)

    model.test("NormUnbalance", 1e-6, 25, 9)
    model.algorithm("Newton")
    model.analysis("Static")


    return model, mid

In this problem, we investigate how this programmatic interface can be used to solve a highly nonlinear problem. Consider the shallow arch shown below (Clarke and Hancock, 1990):

This problem exhibits several critical points along the solution path, and requires the use of a sophisticated algorithm to properly traverse these. In this study, we will investigate how the OpenSees framework can be used to compose these algorithms.

We will perform the analysis by creating an OpenSees Model data structure, and using its methods to perform various tasks. A method is simply a function that is linked to a particular instance of a data structure. In this case, the Model data structure holds the geometry and state (i.e., the current values of solution variables) of our structural model. Some methods that you might see being used in this notebook include (you wont need to make any changes to the use of any of these):

Model.integrator(...) This method configures the iteration strategy to be performed in the next increment
Model.analyze(n) This method applies n increments, and between each increment, performs Newton-Raphson iterations.

In this notebook, the Model data structure is created by the arch_model helper function, which we load from the file arch.py:

from arch import arch_model

We’ll also find the following imports convenient:

import numpy as np
import matplotlib.pyplot as plt

try:
    import scienceplots
    plt.style.use(["steel"]) #(["ieee", "science", "notebook"])
except:
    pass

Solution Strategy

Our general strategy is implemented in the following function solve(). This function adopts an incremental approach where the load is applied in small steps of varying size. The arguments to the function are:

model: an OpenSees Model object
node : an integer indicating which node to collect results from.

Both of these arguments will be supplied by the arch_model helper function mentioned above.

def analyze(model, mid, increment, steps, dx, *args):
    # Initialize some variables
    xy = []      # Container to hold solution history (i.e., load factor and displacement at `node`)
    status = 0   # Convergence flag; Model.analyze() will return 0 if successful.

    # Configure the first load increment strategy; explained below
    increment(model, dx, *args)

    for step in range(steps):

        # 1. Perform Newton-Raphson iterations until convergence for 1 load
        #    increment
        status = model.analyze(1)

        # 2. Store the displacement and load factor
        xy.append([model.nodeDisp(mid, 2), model.getTime()])

        # 3. If the iterations failed, try cutting
        #    the increment arc-length in half
        if status != 0:
            dx /= 2
            increment(model, dx, *args)

    return np.array(xy).T

The strategies used by Clarke and Hancock are:

Solution 1: Iterative strategy: Constant load (Section 3.1)
Load incrementation strategy: Direct incrementation of the load parameter (Section 4.1.1)
Solution 2: Iterative strategy: Constant vertical displacement under the load, $u_{\mathrm{mid}}$ (Section 3.2)
Load incrementation strategy: Incrementation of the displacement component $u_{\mathrm{mid}}$ (Section 4.1.2)
Solution 3: Iterative strategy: Constant arc-length (Section 3.3)
Load incrementation strategy: Incrementation of the arc-length (Section 4.1.3)
Solution 4: Iterative strategy: Minimum unbalanced displacement norm (Section 3.5)
Load incrementation strategy: Incrementation of the arc-length (Section 4.1.3)
Solution 5: Iterative strategy: Constant weighted response (Section 3.7, equation (39))
Load incrementation strategy: Incrementation of the arc-length (Section 4.1.3)
Solution 6: Iterative strategy: Minimum unbalanced force norm (Section 3.6)
Load incrementation strategy: Using the current stiffness parameter (Section 4.2, equation (57))
Solution 7: Iterative strategy: Minimum unbalanced force norm (Section 3.6)
Load incrementation strategy: Incrementation of the arc-length (Section 4.1.3)
Solution 8: Iterative strategy: Constant arc-length (Section 3.3)
Load incrementation strategy: Using the current stiffness parameter (Section 4.2, equation (57))

def solution0(model, dx):
    model.integrator("LoadControl", 400.0)


def solution1(model, dx):
    Jd = 5
    model.integrator("LoadControl", dx, Jd, -800., 800.)


def solution2(model,  dx, *args):
    Jd  = 5
    mid, dof = args
    model.integrator("DisplacementControl", mid, dof, dx, Jd)


def norm_control(model, dx, *args):
    Jd  = 15
    model.integrator("MinUnbalDispNorm",  dx, Jd, -10, 10, det=True)


def arc_control(model, dx, *args,  a=0):
    model.integrator("ArcLength", dx, a, det=True, exp=0.0, reference="point")

fig, ax = plt.subplots()

# x, y = solution0(*arch_model(), 6, 400.0)
# ax.plot(-x, y, 'x', label="S0")

# x, y = analyze(*arch_model(), solution1, 6, 400.0)
# ax.plot(-x, y, 'x', label="S1")
# print(y)

model, mid = arch_model()
x, y = analyze(model, mid, solution2,    7, -150, *(mid, 2))
ax.plot(-x, y, 'o', label="S2")

x, y = analyze(*arch_model(), solution2,  536, -1.5, *(mid, 2))
ax.plot(-x, y, '-', label="S2")

# x, y = analyze(*arch_model(), arc_control, 9500, 0.5, 0)
# ax.plot(-x, y, "-", label="arc")


x, y = analyze(*arch_model(), arc_control, 110,  45)
ax.plot(-x, y, "x", label="arc")

# x, y = analyze(*arch_model(), arc_control,  80,  88, 0)
# ax.plot(-x, y, "+", label="arc")

# x, y = analyze(*arch_model(), arc_control,  80, 188, 0)
# ax.plot(-x, y, "*", label="arc")

# x, y = analyze(*arch_model(), arc_control, 8000, 0.8, 0)
# ax.plot(-x, y, "x", label="arc")

# x, y = analyze(*arch_model(), norm_control,  7000, 1.0)
# ax.plot(-x, y, "-", label="norm")


ax.set_title("Displacement vs Load Factor")
ax.set_ylabel(r"Load factor $\lambda$")
ax.set_xlabel("Displacement $u$")

ax.set_xlim([0, 1200])
ax.set_ylim([-800, 3000])
fig.legend();

Output:

<Figure size 2560x1920 with 1 Axes>

fix, ax = plt.subplots()
ax.plot(-x, '.');
ax.set_title("Displacement vs Step Number")
ax.set_ylabel("Displacement $u$")
ax.set_xlabel("Analysis step $n$");

Output:

<Figure size 2560x1920 with 1 Axes>

The following animation of the solution is created in Animating.ipynb

References

Clarke, M.J. and Hancock, G.J. (1990) ‘A study of incremental‐iterative strategies for non‐linear analyses’, International Journal for Numerical Methods in Engineering, 29(7), pp. 1365–1391. Available at: https://doi.org/10.1002/nme.1620290702 .

Moment-Curvature Analysis

Cable Stayed

On this page: