  Spherical Mat\'ern--Whittle Gaussian random fields are considered as
solutions to fractional elliptic stochastic partial differential equations on
the sphere. Approximation is done with surface finite elements. While the
non-fractional part of the operator is solved by a recursive scheme, a
quadrature of the Dunford--Taylor integral representation is employed for the
fractional part. Strong error analysis is performed, obtaining polynomial
convergence in the white noise approximation, exponential convergence in the
quadrature, and quadratic convergence in the mesh width of the discretization
of the sphere. Numerical experiments for different choices of parameters
confirm the theoretical findings.
