Discontinuities in moment, deflection diagrams

[connorferster]

Describe the bug There seems to be inconsistent results on a beam model.

1. My checks on the sum of the forces does not seem to check out with the statics check.
2. There seems to discontinuties (steps) in both the moment and deflection diagram which are unexpected. These steps increase in magnitude as the point load is increased.

To Reproduce I have created a gist that shows the model that I am using: https://gist.github.com/connorferster/e50aff86e552bb6249fde9427af1c095

Expected behavior No steps in moment diagram. No steps in deflection diagram.

Screenshots Images are in the gist.

Additional context I am not sure if I am missing something but there definitely seems to be some strange behaviour going on.

[JWock82]

Yikes. I’ll have a look this week. I think I might already know what’s going on. It might be a precision error discrepancy between the user defined load position and the calculated beam length.

[JWock82]

I've started looking into this. I've confirmed it's a precision error for loads at the end of the beam. Move the load one millimeter to the left and the diagrams look correct. I've got some more digging to do to get to the bottom of it, but at least I know what I'm looking for a little better. I'll keep working at it - should have a fix soon. In the meantime, if you really want the load at the end of the member, I suggest a nodal load instead of a member load.

[JWock82]

Just found it! That was easier than I thought. I had introduced this error when I created physical members late last year. Bug squished!!! Look for the fix in v0.0.72. I noticed in your code that you're still defining E and G for each member. Note that you'll need to start using materials to use v0.0.72. The updated examples show how to do that. It's pretty painless.

[connorferster]

Thanks so much! Is it possible that this same bug could have caused a singular matrix? @aedificatio found this bug during the PfSE course. He also reported another issue that I think could be from the same bug: in one beam, of length 4800, it analyzes properly and everything checks out. However, if he changes the length of the beam to 4900 (with nodes and everything being updated as required), he gets a singular matrix. I tried reproducing it but, for some reason, I could not trigger the singular matrix using his code on my machine. The weird part was that both of our global K matrices were identical.

Thoughts?

[connorferster]

By the way, the physical members are an excellent feature. So convenient!

[JWock82]

A singular matrix is typically an instability in the model. I also noticed in your model that global rotation about the Y-axis was not restrained:

beam_model.def_support("N1", 1, 1, 1, 1, 0, 0) # @ 1000
beam_model.def_support("N2", 0, 1, 0, 0, 0, 0) # @ 3600

There are 6 global degrees of freedom that must be constrained to prevent rigid body motion: 3 translations (X, Y, Z) and 3 rotations (RX, RY, RZ). In the support definitions, translation in 3 directions is constrained, rotation about the global X is constrained at N1 by the rotationally fixed support, and rotation about the global Z is constrained by the additional vertical support at N2 that provides a reaction "couple". There is nothing to constrain the 6th degree of freedom: rotation about global Y (or any axis parallel to it). What's to stop the beam from rotating about the global Y-axis at N1? Think if someone were to push laterally on your beam at N2, what would prevent rotation? Nothing. It's static, but it's unstable. My guess is that you're getting inconsistent behavior because of tiny precision errors in unstable models that differ across machines. I was surprised the support definitions above worked on my machine, but they did. I think the more correct way to model this would be as follows:

beam_model.def_support("N1", 1, 1, 1, 1, 0, 0) # @ 1000
beam_model.def_support("N2", 0, 1, 1, 0, 0, 0) # @ 3600

The extra support in the Z-direction will have a zero (or near zero - precision errors) reaction, since all your loads are in the XY plane, but it ensures stability if you were to apply a load in the Z direction. Keep in mind that a precision error induced load of 0.000000000001 is a load in an unstable direction nonetheless. Think to yourself, how would this thing respond in the real world if I were to apply a load in any direction? That can help you identify where to add supports and what types to use.

A lot of commercial programs will introduce weak springs to address instabilities, which I think can become a crutch if we're not careful. If a program is doing this, it's because there's an error (even if only a minor one) in the model, and possibly an oversight by the engineer. PyNite is not so forgiving. It requires the user to ensure their model is stable, not just static.

[aedificatio]

Glad to see this issue being resolved so quickly. Furthermore I can confirm that adding the DZ supports for a 2D geometry makes the singular matrix error disappear. Thanks @JWock82 and @connorferster

[JWock82]

One other feature I forgot to about: if you ever need help identifying where an instability is in your model, you can set the check_stability flag to True when you run the analyze method, and Pynite will isolate where the instabilities are for you. It does add a little to the solve time as Pynite has to traverse the diagonal of the stiffness matrix in search of instabilities, so it defaults to False.