  We apply polynomial approximation methods -- known in the numerical PDEs
context as spectral methods -- to approximate the vector-valued function that
satisfies a linear system of equations where the matrix and the right hand side
depend on a parameter. We derive both an interpolatory pseudospectral method
and a residual-minimizing Galerkin method, and we show how each can be
interpreted as solving a truncated infinite system of equations; the difference
between the two methods lies in where the truncation occurs. Using classical
theory, we derive asymptotic error estimates related to the region of
analyticity of the solution, and we present a practical residual error
estimate. We verify the results with two numerical examples.
