The new
c-number calculus for functions of noncommuting operators, developed in Paper I and employed in Paper II to formulate a general phase-space description of boson systems, deals with situations involving equal-time operators only. In the present paper extensions are presented for the treatment of problems involving boson operators at two or more instants of time. The mapping of time-ordered products onto
c-number functions is studied in detail. The results make it possible to evaluate time-ordered products of boson operators by phase-space techniques. The usual Wick theorem for boson systems is obtained as a special case of a much more general theorem on time ordering. Our method of derivation appears to provide the first direct proof of Wick's theorem as well as a clear insight into its true meaning. A closed expression is also obtained for the time-evolution operator in terms of the solution of the
c-number differential equation for the phase-space equivalent of this operator. The new calculus is also applied to the problem of evaluating normally ordered time-ordered, and also the antinormally ordered time-ordered, correlation functions.
The applicability of the so-called truncated Wigner approximation $(\ensuremath{-}W)$ is extended to multitime averages of Heisenberg field operators. This task splits naturally in two. First, what class of multitime averages the $\ensuremath{-}W$ approximates and, second, how to proceed if the average in question does not belong to this class. To answer the first question, we develop a (in principle, exact) path-integral approach in phase space based on the symmetric (Weyl) ordering of creation and annihilation operators. These techniques calculate a new class of averages which we call time-symmetric. The $\ensuremath{-}W$ equations emerge as an approximation within these path-integral techniques. We then show that the answer to the second question is associated with response properties of the system. In fact, for two-time averages, Kubo's renowned formula relating the linear-response function to two-time commutators suffices. The $\ensuremath{-}W$ is directly generalized to the response properties of the system allowing one to calculate approximate time normally ordered two-time correlation functions with surprising ease. The techniques we develop are demonstrated for the Bose-Hubbard model.