Quaternions without unity constraint are used as configuration variables for rotational degrees of freedom of Cosserat rod models, thereby naturally incorporating inflation as well as bending, twisting, extension, and shear deformations of elongate robotic manipulators. The configuration space becomes isomorphic to a subspace of 7-D real-valued functions; thus, an unconstrained local minimizer of total potential energy is a static equilibrium. The ensuing calculus of variations is automated by computer algebra to derive weak-form integral equations that are easily translated to a finite-element package for efficient computation using internal forces. Discontinuities in strain variables are handled in a numerically reliable way. Inextensible, unshearable rod models are derived simply by taking limits of corresponding stiffness parameters. The same procedure facilitates unified software code for both flexible and rigid segments. Simulation experiments with an inflating tube, a helical coil, and a magnetic catheter produce good-quality results and indicate that the computational effort of the proposed method is about two orders of magnitude less than common 3-D finite-element models of large deformation nonlinear elasticity.