# Are linear beam element solutions ever exact and when are they used in practice?

In qualifying the results of a Finite Element (FE) beam model, I'm realizing that I dont understand FE beam elements as well as I should.

Is there a procedure to determine what order beam element is required?

Specifically, for the 2D tip-loaded cantilever beam shown below: The analytical solution given by Euler-Bernoulli Beam Theory is used to qualify my FE results, where:

• $y(x) = -\frac{F}{6EI_{zz}}(3Lx^{2}-x^{3}) \; \rightarrow \; \therefore y(x=L) = -\frac{FL^{3}}{3EI_{zz}}$
• $\sigma_{xx} = \pm \frac{M(x)c}{I_{zz}}$
• $\tau_{xy} = \frac{F}{2I_{zz}}(\frac{h^{2}}{2I_{zz}}-y^{2}) \; \rightarrow \; \therefore \tau_{xy}(y) = \frac{3}{2}\frac{V}{A}$

I am testing Beam 188 elements in ANSYS using Linear (2 Node) and Quadratic (3 Node) beam elements. For both linear and quadratic elements, I increment the number of elements required until the FE model convergences with the analytical (Euler-Bernoulli) solution.

For the single element, linear Beam 188 solution: • The maximum displacement is ~7% less than the analytic solution!
• The maximum bending stress is 50% less than the analytic solution!

I thought that the stiffness matrix for a linear beam element was formulated from the differential beam equations, making it an exact solution for a beam subjected to point loads. As I understand it, higher order elements are required for distributed loading,curved beams, etc.

I understand the stress gradient cannot be captured without additional elements, but is it possible to calculate stress along the beam from the beam equations and interpolated (shape function) values? Still, are the nodal displacements accurate for a single element?

Are linear beam element solutions ever exact and when are they used in practice?

• Nothing is ever exact in engineering. Even your analytical calculation is not exact with respect to reality as it makes simplifying assumptions that do not hold true in reality. In this case looking at your FEA soution one of these assumptions is different in the FEA model, your FEA model accounts for rotation, partially anyway, while your analytic solution does not. – joojaa Sep 22 '16 at 5:54
• @joojaa Thanks for your interest. It is true that ANSYS uses Timoshenko or Mindlin beam elements. While this does account for some error in comparison to Euler-Bernoulli, Im wondering why the result doesnt better reflect solutions obtained analytically. BTW, all beam theories do account for rotation (defined as $ROTZ=\frac{dy}{dx}$ etc.) . – OnStrike Sep 22 '16 at 6:04
• not the rotation of the beam cross section. – joojaa Sep 22 '16 at 6:16
• @joojaa As I understand it: Euler-Bernoulli assumes the cross section is perpendicular to the N.A. (best for slender beams); Timoshenko assumes cross sections remain planar but does not assume they are perpendicular to the N.A.; My knowledge of Mindlin elements is minimal but I dont think they require the cross section be planar. Regardless, the difference between theories is negligible compared to the FE error. – OnStrike Sep 22 '16 at 6:22
• Euler-bernulli assumes the cross section is constant and stays in same orientation over the entire movement, in other words does not rotate. That is why you can just roll the Izz into the equation (all of these therefore are not exact). But given that all of these are approximation what exactly is exact in your mind, exact to the model? Or exact to reality. – joojaa Sep 22 '16 at 6:23

I'm not familiar with ANSYS, but the basic reason would seem to be that the ANSYS beam is not exactly the same formulation as the OP's version of Euler beam theory.

Many FE systems have beam elements formulated using Euler-Timoshenko beam theory which includes some additional terms.

The OP's link to the documentation implies this, in the quote

These elements support an elastic relationship between transverse shear forces and transverse shear strains. Users may specify transverse shear stiffnesses using real constants.

If your input generator allows you to define the section as a geometric shape (i.e. a rectangle) the relevant factors will probably be calulcated automatically. There will probably be some way to switch the extra terms off, if you dig deeper into the documentation.

But, the extra terms should make the beam more flexible, i.e. the tip deflection of the FE model should have been bigger than the Euler beam theory value, not smaller as the OP reported. Maybe there are other differences in the formulation, or the OP did a non-linear analysis instead of a linear one, or whatever....

As for the "50% discrepancy between the maximum bending stress", the OP's plot seems to show constant stress along the length of the beam which is obviously wrong. This may be an artefact of the way ANSYS presents the output, and the plotting program displays it. The stress at the mid-point along the length of the beam should of course be 50% of the stress at the fixed end. You always need to take care to understand exactly what data your "pretty colored and/or animated output pictures" are based on, before you actually believe them!

You have to understand that discretization with both linear and quadratic elements converge in the limit. The accuracy you acheive depends on both the order and size of the element. More elements should result in greater accuracy. You will never acheive "exact", analytical results. The only way to tell if your solution makes sense is to keep refining your mesh and verifying that the error keeps getting smaller (watch out, and make sure you don't refine the mesh too much though... this can lead to loss of precision due to numerical subtraction of nearly equal values from assembling an fem system from mesh points that are too close to each other).