  For Euclidean quantum field theories, Holland and Hollands have shown
operator product expansion (OPE) coefficients satisfy "flow equations": For
interaction parameter $\lambda$, the partial derivative of any OPE coefficient
with respect to $\lambda$ is given by an integral over Euclidean space of a sum
of products of other OPE coefficients. In this paper, we generalize these
results for flat Euclidean space to curved Lorentzian spacetimes in the context
of the solvable "toy model" of massive Klein-Gordon scalar field theory, with
$m^2$ viewed as the "self-interaction parameter". Even in Minkowski spacetime,
a serious difficulty arises from the fact that all integrals must be taken over
a compact spacetime region to ensure convergence but any integration cutoff
necessarily breaks Lorentz covariance. We show how covariant flow relations can
be obtained by adding compensating "counterterms" in a manner similar to that
of the Epstein-Glaser renormalization scheme. We also show how to eliminate
dependence on the "infrared-cutoff scale" $L$, thereby yielding flow relations
compatible with almost homogeneous scaling of the fields. In curved spacetime,
the spacetime integration will cause the OPE coefficients to depend non-locally
on the spacetime metric, in violation of the requirement that quantum fields
should depend locally and covariantly on the metric. We show how this
potentially serious difficulty can be overcome by replacing the metric with a
suitable local polynomial approximation about the OPE expansion point. We
thereby obtain local and covariant flow relations for the OPE coefficients of
Klein-Gordon theory in curved Lorentzian spacetimes. In an appendix, we develop
an algorithm for constructing local and covariant flow relations beyond our
"toy model" based on the associativity properties of OPE coefficients, and we
apply our method to $\lambda\phi^4$-theory.
