Rigid Offsets

2 min read • 419 words

Geometric transformations are employed to model rigid joint offsets in frames.

In this example we investigate the offset argument for frame geometric transformations. The syntax is as follows:

m.geomTransf("Linear", 1, (0, 0, 1),
             offset={
                1: ( off/np.sqrt(2), off/np.sqrt(2), 0),
                2: (-off/np.sqrt(2),-off/np.sqrt(2), 0),
             }
)

geomTransf Linear 1  0 0 1  -offset {
  1  0.35355339059327373  0.35355339059327373 0
  2 -0.35355339059327373 -0.35355339059327373 0
}

Validation

Recall the stiffness of a linear Euler-Bernoulli Beam:

\boldsymbol{K}_e = \begin{bmatrix} EA/L_e & 0 & 0 \\ 0 & 4 EI/L_e & 2 EI/L_e \\ 0 & 2 EI/L_e & 4 EI/L_e \\ \end{bmatrix}

This is transformed into the angled configuration using the transformation:

\boldsymbol{A}_v = \begin{bmatrix} -dX/L_e & dY/L_e & 0 & dX/L_e \\ -dY/L_e^2 & -dX/L_e^2 & 0 & dY/L_e^2 \\ -dY/L_e^2 & -dX/L_e^2 & 1 & dY/L_e^2 \\ \end{bmatrix}

\boldsymbol{A}_o = \begin{bmatrix} 0 & 1 & 0 \\ 1 & 0 & -L_o/\sqrt{2} \\ 0 & 0 & 1 \\ 0 & 0 & L_o/\sqrt{2} \\ \end{bmatrix}

diamond.py (collapsed)

"""
               ^ P
               |
               .<|o
             .
           .
   y ^   .
     | .
    _._ ---> x
    o o

"""
import opensees.openseespy as ops
import numpy as np
import veux
import os

def diamond_solution(EI, EA, L, off=0):
    import numpy as np
    Le = L - 2*off
    dX = dY = Le/np.sqrt(2)
    Ag = np.array([[-dX/Le   ,  dY/Le   , 0,  dX/Le   ],
                   [-dY/Le**2, -dX/Le**2, 0,  dY/Le**2],
                   [-dY/Le**2, -dX/Le**2, 1,  dY/Le**2]])
    Ao = np.array([[0, 1,              0],
                   [1, 0,-off/np.sqrt(2)],
                   [0, 0,              1],
                   [0, 0, off/np.sqrt(2)]])

    A  = Ag@Ao
    Ke = np.array([[EA/Le,       0,       0],
                   [    0, 4*EI/Le, 2*EI/Le],
                   [    0, 2*EI/Le, 4*EI/Le]])

    K = A.T@Ke@A

    kuw = K[1:,0]
    return 1/np.linalg.solve(K, [1,0,0])
#   return K[0][0] - kuw.dot(np.linalg.solve(K[1:,1:], kuw))


def create_diamond():
    m = ops.Model(ndm=3, ndf=6)

    L   = 5
    off = 0.1*L
    X   = L / np.sqrt(2)
    E   = 30
    I   = 1
    A   = 500*I

    m.node(1, (0, 0, 0))
    m.node(2, (X, X, 0))
    m.fix(1, (0, 1, 1,  1, 1, 1))
    m.fix(2, (1, 0, 1,  1, 0, 0))


    m.geomTransf("Linear", 1, (0, 0, 1),
                  # jntOffset=( off/np.sqrt(2), off/np.sqrt(2), 0,
                  #            -off/np.sqrt(2),-off/np.sqrt(2), 0),
                offset={
                    1: ( off/np.sqrt(2), off/np.sqrt(2), 0),
                    2: (-off/np.sqrt(2),-off/np.sqrt(2), 0),
                }
    )

    if "CRD04" in os.environ:
        section, element = "ElasticFrame", "ForceFrame"
    else:
        section, element = "Elastic", "forceBeamColumn"

    m.section(section, 1, E=E, A=A, Iy=I, Iz=I, G=1, J=3/2*I)
    m.element(element, 1, (1, 2), transform=1, section=1, shear=0)

    m.pattern("Plain", 1, "Linear", loads={2: (0, 1, 0,   0, 0, 0)})
    P = -4*L**2/(E*I)
    m.integrator("LoadControl", P)
    m.analysis("Static")
#   m.test("NormDispIncr", 1e-8, 2)
    m.analyze(1)

#   m.print(json=True)

    print(P/diamond_solution(E*I, E*A, L, off=off)/L)
    print(m.nodeDisp(2, 2)/L)
    print(m.numIter())
    return m


if __name__ == "__main__":
    model = create_diamond()
    a = veux.create_artist(model)
    a.draw_outlines()
    a.draw_origin()
    a.draw_outlines(state=model.nodeDisp)
    a.draw_axes()
    veux.serve(a)

Pure flexure

Shell undergoing finite rotations