# Most Solvers Can Be Marplots

Have you ever tried to replicate the familiar beam stiffness coefficients $12EI/L^3$, $6EI/L^2$, $4EI/L$, and $2EI/L$ (there’s a poem about them here) by imposing unit displacements and rotations at fixed supports? It should be one of the first sanity checks you make when using or developing new structural analysis software.

You can perform this check in OpenSees by defining a single elastic beam element between two nodes, fixing each node, then using the sp command to impose a unit displacement or rotation at one end. The ensuing nodal reactions should equal the stiffness coefficients.

``````L = 100
E = 29000
A = 20
I = 1400

ops.node(1,0,0); ops.fix(1,1,1,1)
ops.node(2,L,0); ops.fix(2,1,1,1)

ops.geomTransf('Linear',1)
ops.element('elasticBeamColumn',1,1,2,A,E,I,1)

ops.timeSeries('Constant',1)
ops.pattern('Plain',1,1)
ops.sp(2,2,1.0) # Unit displacement at node 2, DOF 2

ops.constraints('Transformation')
ops.system('ProfileSPD') # THIS MATTERS!
ops.analysis('Static')
ops.analyze(1)
ops.reactions()``````

Note that you have to use a constraint handler for the imposed displacement. Although not shown here, you can also apply the same approach to finding fixed-end shears and moments for various member loads.

Unfortunately, working with fixed-fixed beams is easier said than done in OpenSees because many of the equation solvers just can’t handle the fact that there are no equations to solve, i.e., there are no unconstrained degrees of freedom. (It’s not the solver routines that cause problems, it’s the OpenSees code that wraps those routines that causes problems.)

At any rate, I think you would agree that an analysis should not fail because there are no equations to solve. Don’t call the cops on a party to which you were not invited!

Some equation solvers in OpenSees do not check if the number of equations is zero before calling the solver routine. You’ll see an odd error message, e.g., an illegal parameter value going into the `BandSPD` solver routine (LAPACK‘s DPBSV function).

Other solvers check for zero equations and bypass the solver routines, but when the integrator asks the system of equations for the $X$ vector, the poop hits the fan because null pointers will be dereferenced. You’ll see an error message like this when using `ProfileSPD` and many other solvers in OpenSees.

You can give the solver something to work with by using two elements with an interior node (and you’d be half way to passing a patch test). However, the only equation solver I’ve seen in OpenSees that can handle a fixed-fixed beam analysis is `UmfPack`. It checks for zero equations and doesn’t use a pointer to store $X$.

Until we fix these issues with the solvers in OpenSees, party on with `UmfPack` when analyzing fixed-fixed beams, or any other model with all degrees of freedom constrained.

## 6 thoughts on “Most Solvers Can Be Marplots”

1. MSB says:

Hi PD,
Interesting.
I love your blog post .

Liked by 1 person

2. Art says:

Interesting!

Liked by 1 person

3. ZHU says:

This must be a chance to travel deep into OpenSees code. So exciting…..

Liked by 1 person

1. Positive Definite says:

Hello ZHU,
The fixed-fixed beam thing is something where I’ve always just used a workaround like defining an extra node. But, such workarounds should not be necessary!
PD

Like

This site uses Akismet to reduce spam. Learn how your comment data is processed.