As Glen_b pointed out, the specifications of the predict.lm() package in R state the following, thus giving exactly what I wanted:
"The prediction intervals are for a single observation at each case in newdata (or by default, the data used for the fit) with error variance(s) pred.var. This can be a multiple of res.var, the estimated value of σ^2: the default is to assume that future observations have the same error variance as those used for fitting. If weights is supplied, the inverse of this is used as a scale factor. For a weighted fit, if the prediction is for the original data frame, weights defaults to the weights used for the model fit, with a warning since it might not be the intended result. If the fit was weighted and newdata is given, the default is to assume constant prediction variance, with a warning."
Also, the difference between what Faraway gives and does not give is merely an addition of 1 underneath the square root:
For a general answer from Faraway (2002, 39-41) on the (1) prediction interval for a single future response when given $x_0$:
$$
\hat{y}_0 \pm t^{(\alpha/2)}_{n-p} \hat\sigma \sqrt{1 + x^T_0 (X^TX)^{-1} x_0}
$$
For a general answer from Faraway on the (2) prediction interval for the average of responses when given $x_0$:
$$
\hat{y}_0 \pm t^{(\alpha/2)}_{n-p} \hat\sigma \sqrt{x^T_0 (X^TX)^{-1} x_0}
$$
I think that the above formulas look like the following in the simplest case of one predictor variable [please correct me if I have made an error in supposing this special case is equivalent to Faraway's general case above].
First (1):
$$
\hat{y}_i \pm t^{(\alpha/2)}_{{crit}_{(n-p)}}*(s_{y*x} \sqrt{1 + {1 \over n} + {(x_i - \bar x)^2 \over SS_x}})
$$
then (2):
$$
\hat{y}_i \pm t^{(\alpha/2)}_{{crit}_{(n-p)}}*(s_{y*x} \sqrt{ {1 \over n} + {(x_i - \bar x)^2 \over SS_x}})
$$
?predict.lm
, which explains what argument you need to change $\endgroup$