BEAM2 — 3D 2-node Euler–Bernoulli beam#

Kinematics. Two-node straight beam with six DOFs per node (three translations + three rotations); 12 DOFs per element.

Stiffness. Axial + torsion + two independent bending planes, each assembled in closed form in the local frame and rotated to global.

Mass. Consistent mass matrix with translational + rotational (section inertia) contributions.

Limit. Slender (Euler–Bernoulli) form with shear correction \(\Phi \to 0\); the shear-dominated Timoshenko range (finite \(\Phi\)) is roadmap.

Beam kinematics and the slender limit#

BEAM2 is the classical Timoshenko beam [Hughes1987] [Bathe2014]: a straight reference axis with six kinematic DOFs at every cross- section (three translations + three rotations), plus the shear correction

\[\Phi = \frac{12\, E\, I}{\kappa\, G\, A_s\, L^2},\]

which collapses to zero when the shear area \(A_s\) is infinite or the length-to-depth ratio is large — the slender (Euler–Bernoulli) limit. The stiffness matrix below is the \(\Phi \to 0\) form; femorph-solver’s current BEAM2 kernel implements exactly this limit.

Local DOF ordering (per node): \((u, v, w, \theta_x, \theta_y, \theta_z)\). Global assembly fills rows / cols 0–5 for node I and 6–11 for node J in the canonical \((U_X, U_Y, U_Z, \text{ROT}_X, \text{ROT}_Y, \text{ROT}_Z)\) order.

Decoupled axial / torsion / bending blocks#

With the local \(\hat{x}\) along I → J, the four elastic contributions to \(K_\text{local}\) are decoupled:

  • Axial (\(u_I, u_J\)):

    \[\begin{split}K_a = \frac{E A}{L} \begin{bmatrix} \,\phantom{-}1 & -1 \\ -1 & \phantom{-}1 \end{bmatrix}.\end{split}\]

  • Torsion (\(\theta_{xI}, \theta_{xJ}\)):

    \[\begin{split}K_t = \frac{G J}{L} \begin{bmatrix} \,\phantom{-}1 & -1 \\ -1 & \phantom{-}1 \end{bmatrix}, \qquad G = \frac{E}{2 (1 + \nu)}.\end{split}\]

  • Bending about local z (\(v_I, \theta_{zI}, v_J, \theta_{zJ}\)):

    \[\begin{split}K_{bz} = \frac{E I_{zz}}{L^3} \begin{bmatrix} \phantom{-}12 & \phantom{-}6L & -12 & \phantom{-}6L \\ \phantom{-}6L & \phantom{-}4L^2 & -6L & \phantom{-}2L^2 \\ -12 & -6L & \phantom{-}12 & -6L \\ \phantom{-}6L & \phantom{-}2L^2 & -6L & \phantom{-}4L^2 \end{bmatrix}.\end{split}\]

  • Bending about local y (\(w_I, \theta_{yI}, w_J, \theta_{yJ}\)):

    \[\begin{split}K_{by} = \frac{E I_{yy}}{L^3} \begin{bmatrix} \phantom{-}12 & -6L & -12 & -6L \\ -6L & \phantom{-}4L^2 & \phantom{-}6L & \phantom{-}2L^2 \\ -12 & \phantom{-}6L & \phantom{-}12 & \phantom{-}6L \\ -6L & \phantom{-}2L^2 & \phantom{-}6L & \phantom{-}4L^2 \end{bmatrix}.\end{split}\]

The sign convention on \(K_{by}\) differs from \(K_{bz}\) because deflection along local \(\hat{z}\) is driven by rotation \(\theta_y\) through the opposite sense of the right-hand rule [Cook2001].

Orientation#

Without an explicit orientation node, femorph-solver picks the local \(\hat{y}\) so that local \(\hat{z}\) stays close to world \(+Z\) whenever the beam is not parallel to the global Z-axis, falling back to world \(+Y\) for vertical beams. For axisymmetric sections (\(I_{yy} = I_{zz}\)) the orientation is immaterial.

Consistent mass#

The 12 × 12 local consistent mass [Cook2001] carries a translational part \(\rho A L\) and a rotational / transverse-inertia part coupling deflection to rotation. The closed-form blocks for each plane are the classical Euler–Bernoulli consistent mass; femorph-solver implements them directly — no quadrature.

Real constants#

Slot

Symbol

Meaning

real[0]

\(A\)

Cross-sectional area.

real[1]

\(I_{zz}\)

Bending inertia about local z (resists y-deflection).

real[2]

\(I_{yy}\)

Bending inertia about local y (resists z-deflection).

real[3]

\(J\)

Saint-Venant torsion constant.

Material#

EX (Young’s modulus), PRXY (Poisson’s ratio; only used to derive \(G = E / (2 (1 + \nu))\) for torsion), DENS (mass density for \(M_e\)).

Validation#

Cantilever first mode. Euler–Bernoulli theory gives the first natural frequency of a clamped–free beam [Timoshenko1937] [Cook2001] as

\[f_1 = \frac{\beta_1^2}{2 \pi L^2}\, \sqrt{\frac{E I}{\rho A}}, \qquad \beta_1 L \approx 1.8751.\]

The BEAM2 cantilever benchmark under tests.elements.beam188 reproduces this to sub-percent on a 20-element mesh.

Axial natural frequency. A free-free bar of length \(L\) has the classical wave modes \(f_n = (n / 2) \sqrt{E / \rho} / L\); the BEAM2 axial block reproduces these exactly (one DOF per node per mode family).

API reference#

class femorph_solver.elements.beam2.Beam2[source]#

Bases: ElementBase

static distributed_load(coords: ndarray, real: ndarray, material: dict[str, float], *, face: int, value: float) → ndarray[source]#

Consistent nodal-force vector for a uniform distributed load.

Implements the Hermite-cubic equivalent nodal forces for a uniform line load q = value (force per unit length) applied along the element axis, perpendicular to the requested local face. Standard textbook result (Cook §5.7 Table 5.7-1) on a Hermite-cubic Euler-Bernoulli beam:

F_v_each_end = q · L / 2 M_each_end = ± q · L² / 12 (signs flip end-to-end)

Face conventions (BEAM188 / BEAM2 local frame, X = element axis):

  • face 1 → +Y (load acts in local +Y direction)

  • face 2 → +Z

  • face 3 → -Y

  • face 4 → -Z

Returns a 12-element local-frame load vector laid out as [ux_I, uy_I, uz_I, θx_I, θy_I, θz_I, ux_J, uy_J, uz_J, θx_J, θy_J, θz_J] already transformed to the global frame via the same rotation matrix R _k_local() uses, so the caller can scatter it directly into F.

References

  • Cook et al., Concepts and Applications of Finite Element Analysis, 4th ed., Table 5.7-1 (consistent equivalent nodal forces for uniform / triangular loads on Hermite beams).

  • Logan, D.L., A First Course in the Finite Element Method, 5th ed., Cengage, 2012, §5.6 (work-equivalent loads).

static fiber_stress(coords: ndarray, material: dict[str, float], real: ndarray, u_e: ndarray, *, fiber_y: float, fiber_z: float = 0.0, station_s: float = 0.0) → float[source]#

Axial stress σ_xx at a section fiber along a Beam2 element.

Combines axial, strong-axis bending, and weak-axis bending:

σ_xx(s, y, z) = N(s)/A − M_z(s)·y/Iz + M_y(s)·z/Iy

where s is the natural station along the element axis (s = 0 at node I, s = L at node J), (y, z) are the fiber’s location relative to the section centroid in the beam’s local frame, and (N, M_z, M_y) are the internal forces recovered from the element’s nodal displacement vector u_e via the same shape functions used in Beam2.ke():

  • Linear axial: N(s) = E·A · (u_J − u_I)/L (constant)

  • Hermite-cubic bending about local z (v(s)): M_z(s) = E·Iz · v''(s) — linear in s.

  • Hermite-cubic bending about local y (w(s)): with the dw/dx = −θ_y sign convention used by _k_local, M_y(s) = −E·Iy · w''(s).

Parameters:

  • coords ((2, 3)) – Element node coordinates in global frame.

  • material – Same as ke(). real carries (A, Iz, Iy, J) in slots 0..3.

  • real – Same as ke(). real carries (A, Iz, Iy, J) in slots 0..3.

  • u_e ((12,)) – Per-element global-frame DOF vector laid out as [ux_I, uy_I, uz_I, θx_I, θy_I, θz_I, ux_J, …].

  • fiber_y – Fiber location relative to the section centroid in the beam’s local frame. y is the depth direction (the orientation node defines local +y), z is the cross-bending direction. Defaults to z = 0.

  • fiber_z – Fiber location relative to the section centroid in the beam’s local frame. y is the depth direction (the orientation node defines local +y), z is the cross-bending direction. Defaults to z = 0.

  • station_s – Natural station s ∈ [0, L] (default 0 — node I).

References

  • Cook, R.D., Malkus, D.S., Plesha, M.E., Witt, R.J., Concepts and Applications of FEA, 4th ed., Wiley, 2002, §2.5 (Hermite-cubic Euler-Bernoulli beam, internal moment recovery from second derivative of the deflection field).

  • Logan, D.L., A First Course in the Finite Element Method, 5th ed., Cengage, 2012, §5.4 (axial-bending superposition for σ_xx at a section fiber).

static initial_static_nonlinear_state(coords: ndarray, material: dict[str, float], real: ndarray, plastic: ElastoplasticMaterial | None, fiber_grid: FiberGrid | None = None)[source]#

Initial plastic state — shape depends on plasticity layout.

  • No fibre grid: a single IntegrationPointState (axial- only, matches the Phase 3 part 1 path).

  • Fibre grid: a list-of-lists indexed [axial_station] [fiber_idx] — each fibre at each axial Gauss point carries its own (ε_p, α, ε_p_cum).

static ke(coords: ndarray, material: dict[str, float], real: ndarray) → ndarray[source]#

Return the element stiffness matrix in global DOF ordering.

Parameters:

  • coords ((n_nodes, 3) float64)

  • material (mapping with property names as keys ("EX", "PRXY", "DENS", …))

  • real (1-D array of real constants)

static me(coords: ndarray, material: dict[str, float], real: ndarray, lumped: bool = False) → ndarray[source]#

Return the element mass matrix in global DOF ordering.

static strain_batch(coords: ndarray, u_e: ndarray, material: dict[str, float] | None = None) → ndarray | None[source]#

Per-element axial-only 3D Voigt strain at the two end nodes.

Centroid (neutral-axis) recovery only - bending fibre stress is a per-section follow-up that needs y_f / z_f positions from the registered BeamSectionProperties. The recovered strain corresponds to the pure-uniaxial centroid stress state sigma = E * (du_axial / L) * (e^e) where e is the element-axis unit vector and du_axial = (u_J - u_I) . e.

Mirrors Truss2.strain_batch() so compute_nodal_stress() recovers the right global Cartesian stress projection via the standard isotropic C @ eps contraction. Voigt order is [exx, eyy, ezz, gxy, gyz, gxz] (engineering shears).

static stress_update_static_nonlinear(coords: ndarray, material: dict[str, float], real: ndarray, u_e_global: ndarray, state, plastic: ElastoplasticMaterial | None, fiber_grid: FiberGrid | None = None) → tuple[ndarray, ndarray, object][source]#

Phase-3 plasticity hook on a Beam2 — axial or fibre.

Two integration regimes are supported, keyed on fiber_grid:

  • Axial-only (fiber_grid is None): the section’s axial fibre is treated as a single integration point. The consistent algorithmic D_t · A / L axial block sits in the element tangent; bending and torsion remain linear- elastic. This matches the Phase 3 part 1 path.

  • Fibre integration (fiber_grid is a FiberGrid): full distributed-plasticity beam-column at 2 axial Gauss stations × n_fiber cross-section fibres. Each fibre runs stress_update_uniaxial() on its uniaxial strain ε_f = ε_axial − y_f · κ_z(s) − z_f · κ_y(s). Section resultants N, M_z, M_y come from fibre sums; the consistent algorithmic D_section couples axial- bending DOFs through Σ_f D_t · y_f · dA etc. Torsion stays linear-elastic.

Parameters:

  • coords ((2, 3)) – Element node coordinates in the global frame.

  • material – Same as ke(). real carries (A, Iz, Iy, J).

  • real – Same as ke(). real carries (A, Iz, Iy, J).

  • u_e_global ((12,)) – Per-element global-frame DOF vector.

  • state (IntegrationPointState) – Trial state at the start of the Newton iteration (committed state from the previous load step until convergence).

  • plastic (ElastoplasticMaterial or None) – The bilinear (BKIN / BISO) law deposited by TB,PLAS,…. None means a fully elastic element — return-map is skipped and the elastic axial stiffness is used.

Returns:

  • K_t_global ((12, 12)) – Consistent algorithmic tangent in the global frame.

  • F_int_global ((12,)) – Internal-force vector in the global frame.

  • new_state (IntegrationPointState) – Trial state after this Newton step’s return-map (axial).

References#

[Hughes1987]

Hughes, T. J. R. (1987 / Dover 2000). The Finite Element Method: Linear Static and Dynamic Finite Element Analysis. Ch. 5 on beams and thin-walled structures. https://store.doverpublications.com/0486411818.html

[Bathe2014]

Bathe, K.-J. (2014). Finite Element Procedures, 2nd ed. Ch. 5 (beams). https://www.klausjurgenbathe.com/fepbook/

[Cook2001] (1,2,3)

Cook, R. D., Malkus, D. S., Plesha, M. E., and Witt, R. J. (2001). Concepts and Applications of Finite Element Analysis, 4th ed. Wiley. Ch. 4 (beams), Ch. 11 (consistent mass). https://www.wiley.com/en-us/Concepts+and+Applications+of+Finite+Element+Analysis%2C+4th+Edition-p-9780471356059

[Timoshenko1937]

Timoshenko, S. (1937). Vibration Problems in Engineering, 2nd ed. D. Van Nostrand Company (reprinted by Wiley). The derivation of the first-mode Euler–Bernoulli cantilever frequency. https://archive.org/details/vibrationproblem00timo