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.

4 thoughts on “Most Solvers Can Be Marplots

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Google photo

You are commenting using your Google account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s

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