In this paper, a finite element model for efficient nonlinear analysis of the mechanical response of viscoelastic beams is presented. The principle of virtual work is utilized in conjunction with the third-order beam theory to develop displacement-based, weak-form Galerkin finite element model for both quasi-static and fully-transient analysis. The displacement field is assumed such that the third-order beam theory admits C0Lagrange interpolation of all dependent variables and the constitutive equation can be that of an isotropic material. Also, higher-order interpolation functions of spectral/hp type are employed to efficiently eliminate numerical locking. The mechanical properties are considered to be linear viscoelastic while the beam may undergo von Krmn nonlinear geometric deformations. The constitutive equations are modeled using Prony exponential series with general n-parameter Kelvin chain as its mechanical analogy for quasi-static cases and a simple two-element Maxwell model for dynamic cases. The fully discretized finite element equations are obtained by approximating the convolution integrals from the viscous part of the constitutive relations using a trapezoidal rule. A two-point recurrence scheme is developed that uses the approximation of relaxation moduli with Prony series. This necessitates the data storage for only the last time step and not for the entire deformation history.