  We construct a finite element approximation of a strain-limiting elastic
model on a bounded open domain in $\mathbb{R}^d$, $d \in \{2,3\}$. The sequence
of finite element approximations is shown to exhibit strong convergence to the
unique weak solution of the model. Assuming that the material parameters
featuring in the model are Lipschitz-continuous, and assuming that the weak
solution has additional regularity, the sequence of finite element
approximations is shown to converge with a rate. An iterative algorithm is
constructed for the solution of the system of nonlinear algebraic equations
that arises from the finite element approximation. An appealing feature of the
iterative algorithm is that it decouples the monotone and linear elastic parts
of the nonlinearity in the model. In particular, our choice of piecewise
constant approximation for the stress tensor (and continuous piecewise linear
approximation for the displacement) allows us to compute the monotone part of
the nonlinearity by solving an algebraic system with $d(d+1)/2$ unknowns
independently on each element in the subdivision of the computational domain.
The theoretical results are illustrated by numerical experiments.
