Rigid Offsets

This example tests the effect of rigid joint offsets in a diamond-shaped frame. The analytical solution is known, and the test is used to verify the correct treatment of offset geometry in a 3D frame element.

This test case consists of a simple two-node frame element arranged in a diamond pattern. Rigid offsets at both ends rotate the element relative to the applied force direction. This configuration is highly sensitive to geometric assumptions, making it ideal for verifying the implementation of joint offsets in 3D frame elements.

Model  

We begin by importing necessary libraries and defining a function that computes the analytical solution using a transformation approach:

import xara
import veux

Next, the model is built using the xara API, which mirrors OpenSeesRT’s structured syntax:

def create_diamond():
    m = xara.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))

The rigid offsets are defined by passing an offset dictionary to the geomTransf command:

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
}

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

A vertical unit load is applied at the free end, and a single analysis step is performed:

    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)

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

Validation  

Recall the stiffness $\boldsymbol{K}_e$ of a linear 2D Euler-Bernoulli beam:

This is transformed into the angled configuration using the transformation:

and joint offsets are similarly applied with:

import numpy as np
def diamond_solution(EI, EA, L, off=0):
    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
    return 1 / np.linalg.solve(K, [1, 0, 0])

Visualization  

The resulting deformed shape is plotted using veux, with the displacement field applied:

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)