Wick's theorem

Wick's theorem is a key result in quantum field theory that expresses the time-ordered product of quantum fields as a sum over all possible normal-ordered products augmented by pairwise contractions of the fields, facilitating the evaluation of correlation functions and scattering amplitudes in perturbative calculations.^[1] Introduced by physicist Gian-Carlo Wick in 1950, the theorem provides an algebraic method to simplify the S-matrix (collision matrix) by reducing complex operator products into manageable terms, building on earlier work by Dyson and Feynman.^[1] In its standard formulation for bosonic fields, the theorem states that for a time-ordered product T(\phi_1 \phi_2 \dots \phi_n), where \phi_i are field operators, it equals the sum over all ways to pair the fields into contractions (denoted by \Delta(x_i - x_j) or the propagator) and the remaining fields in normal order : \dots :, with the number of contractions ranging from 0 to \lfloor n/2 \rfloor.^[2] Contractions represent the difference between time-ordered and normal-ordered two-point functions, effectively capturing vacuum fluctuations or propagators in Feynman diagrams.^[3] For fermionic fields, the theorem includes an additional sign factor of (-1) for each crossing of fermion lines in the contractions.^[3] The theorem's importance lies in its role as a cornerstone of perturbative quantum field theory, enabling the systematic computation of vacuum expectation values of time-ordered products, which vanish for odd numbers of fields in free theories and reduce to sums of fully contracted terms for even numbers.^[2] It directly underpins the derivation of Feynman rules, where each contraction corresponds to a propagator line in diagrams, making it indispensable for calculating processes in quantum electrodynamics, the standard model, and beyond.^[1] Wick's theorem extends naturally to interacting theories via the Dyson series and has applications in statistical mechanics and condensed matter physics for handling operator averages.^[3]

Foundations

Contractions

In quantum field theory, the contraction of two operators A and B is defined as the difference between the vacuum expectation value of their time-ordered product and the vacuum expectation value of their normal-ordered product, denoted as

\langle AB \rangle = \langle T(AB) \rangle - \langle :AB: \rangle,

where T represents the time-ordering operator that arranges the operators in chronological order, and (: \cdot :\ ) denotes normal ordering, which places all creation operators to the left of annihilation operators.^[4]^[1] This definition applies to both bosonic and fermionic operators, with the time-ordering incorporating commutators for bosons and anticommutators for fermions to respect the appropriate statistics.^[4] The contraction operation exhibits several key algebraic properties. It is bilinear, meaning that for scalar coefficients \alpha, \beta and operators A, B, C, D,

\langle \alpha A + \beta B, C \rangle = \alpha \langle A C \rangle + \beta \langle B C \rangle, \quad \langle A, \beta C + \gamma D \rangle = \beta \langle A C \rangle + \gamma \langle A D \rangle.

Contractions vanish for pairs that do not contribute to non-zero vacuum expectations, such as operators from independent fields or those not paired according to the field's commutation relations.^[1] For fermionic operators, the contraction is antisymmetric, satisfying \langle AB \rangle = -\langle BA \rangle, which introduces sign factors in products involving odd permutations of fermion pairs.^[4] The concept of contraction was introduced by Gian-Carlo Wick in 1950 as a foundational tool for simplifying perturbative calculations in quantum field theory, particularly in evaluating matrix elements involving time-ordered products of field operators.^[1]

Normal ordering

In quantum field theory, normal ordering is a prescription for rearranging products of creation and annihilation operators such that all creation operators appear to the left of all annihilation operators, thereby defining a standard form for operator expressions. This procedure was introduced by Wick to handle infinities arising in vacuum expectation values of field operator products. For bosonic operators satisfying the commutation relation [a, a^\dagger] = 1, the normal-ordered product of two operators is given by

: a^\dagger b : = a^\dagger b

if a^\dagger is a creation operator and b is an annihilation operator, or

: b a^\dagger : = a^\dagger b

if the order is reversed, as the commutation relation introduces only a c-number shift without altering the operator structure beyond reordering. The primary purpose of normal ordering is to ensure that the vacuum expectation value of any normal-ordered product vanishes, i.e., \langle 0 | : \mathcal{O} : | 0 \rangle = 0, where |0\rangle is the vacuum state. This property eliminates divergent zero-point contributions in perturbation theory calculations, facilitating the evaluation of correlation functions by isolating the interacting components of operator products. For products involving multiple operators, normal ordering is achieved systematically by successive application of the commutation (or anticommutation) relations to shift all creation operators to the left while preserving the overall structure. In the bosonic case, this process involves no additional phase factors, as interchanges yield only c-number corrections that are discarded in the normal-ordered form. In the fermionic case, operators satisfy anticommutation relations \{c, c^\dagger\} = 1, and normal ordering similarly places creation operators to the left of annihilation operators, but with a sign factor determined by the parity of the permutation required to achieve this order. For two fermionic operators,

: c^\dagger d : = c^\dagger d, \quad : d c^\dagger : = - c^\dagger d,

reflecting the anticommutation that introduces a minus sign upon interchange. For multiple fermions, the extension incorporates a (-1)^P factor, where P is the number of pairwise interchanges (or the parity of the permutation) needed to sort the operators, ensuring consistency with the antisymmetric nature of fermionic statistics. This distinction arises because fermionic normal ordering must account for the intrinsic sign changes under operator exchanges, unlike the bosonic case.

Statement

Finite-dimensional case

Wick's theorem in the finite-dimensional case provides a combinatorial expression for the time-ordered product of a finite number n of free bosonic or fermionic operators A_1, \dots, A_n, typically linear combinations of creation and annihilation operators at specified spacetime points. The theorem equates this product to a sum over all possible normal-ordered products augmented by contractions, where contractions represent pairwise vacuum expectation values. For bosonic fields \phi(x_i), the precise statement is

T\bigl(\prod_{i=1}^n \phi(x_i)\bigr) = \sum : \prod_{i=1}^n \phi(x_i) :

where the sum runs over all possible sets of contractions (including the zero-contraction term, which is the fully normal-ordered product :\prod \phi(x_i):), partial contractions (with some pairs replaced by the Feynman propagator \Delta_F(x_i - x_j) = \langle 0 | T \phi(x_i) \phi(x_j) | 0 \rangle), and, for even n, full contractions consisting of complete pairings of the operators.^[5] In the fermionic case, the structure mirrors the bosonic version, but each term in the sum acquires a sign factor (-1)^\epsilon, where \epsilon counts the parity of the permutation needed to reorder the operators such that the contracted pairs are adjacent in the sequence (arising from the anticommutation relations). Full contractions for even n thus sum over all perfect matchings with these parity signs, while partial contractions leave unpaired operators normal-ordered.^[6]

Infinite-dimensional case

In the infinite-dimensional setting of quantum field theory, Wick's theorem generalizes to time-ordered exponentials of field operators, which is essential for handling continuous degrees of freedom in an infinite-dimensional Hilbert space. The theorem extends the finite-dimensional result by expressing the time-ordered exponential \mathrm{T} \exp\left( \sum_i A_i \right) as a normal-ordered exponential incorporating all possible contractions: \mathrm{T} \exp\left( \sum_i A_i \right) = : \exp\left( \sum_i A_i + \sum \text{contractions} \right) : . This form arises from the cumulant expansion, where the logarithm of the generating functional captures connected contributions, and higher cumulants vanish for free Gaussian theories, reducing the expression to pairwise contractions.^[7] A primary application occurs in the construction of generating functionals for correlation functions. For a free scalar field, the generating functional is given by Z[J] = \langle \mathrm{T} \exp\left( i \int J(x) \phi(x) \, d^4x \right) \rangle = : \exp\left( i \int J(x) \phi(x) \, d^4x + \frac{1}{2} \int J(x) G(x-y) J(y) \, d^4x \, d^4y \right) :, where G(x-y) is the Feynman propagator and the normal ordering ensures the vacuum expectation value is unity. This representation facilitates the extraction of all n-point functions via functional derivatives with respect to J, with contractions providing the propagators in Feynman diagrams.^[8]^[9] Normal ordering plays a crucial role in regulating divergences inherent to infinite-dimensional spaces, such as ultraviolet infinities from zero-point fluctuations. By subtracting the vacuum expectation value, :O: = O - \langle 0 | O | 0 \rangle, the normal-ordered exponential avoids infinite constant shifts, enabling well-defined perturbation theory around the free theory. This regularization is particularly vital in canonical quantization, where the infinite volume of the Hilbert space would otherwise render expressions ill-defined.^[8]^[10] For fermionic fields, the theorem adapts to anticommuting Grassmann variables, replacing bosonic contractions with antisymmetric pairings that introduce signs from permutations. The generating functional involves a pfaffian structure, analogous to the square root of a determinant, yielding Z[\eta, \bar{\eta}] = \exp\left( i \int \bar{\eta} \psi + i \int \eta \bar{\psi} + \int \bar{\eta} S \eta \right) in normal-ordered form, where S is the fermionic propagator; this ensures the correct statistics and handles the finite-dimensional Grassmann integration in the path integral limit.^[11]

Proof

Combinatorial argument

The combinatorial argument provides a proof of Wick's theorem in the finite-dimensional case by expressing the expectation value of a product of Gaussian random variables as a sum over all possible pairings, leveraging the structure of multivariate normal distributions. Consider a zero-mean multivariate Gaussian vector \mathbf{X} = (X_1, \dots, X_n) with covariance matrix C, where n = 2m is even. The theorem states that the moment is

\langle X_1 \cdots X_{2m} \rangle = \sum_{\mathcal{P}} \prod_{(i,j) \in \mathcal{P}} C_{ij},

where the sum runs over all perfect matchings \mathcal{P} (pair partitions) of the $2m indices into m disjoint pairs (i,j) with i < j. For odd n, the moment vanishes due to symmetry. The number of such perfect matchings is the double factorial (2m-1)!! = (2m)! / (2^m m!), which counts the ways to successively pair the indices: choose a partner for the first index in $2m-1 choices, for the next unpaired in $2m-3 choices, and so on.^[12] This result follows from a proof by mathematical induction on the number of variables n. The base cases are straightforward: for n=0, \langle 1 \rangle = 1, corresponding to the empty product with no pairings; for n=2, \langle X_1 X_2 \rangle = C_{12}, the single possible pairing. Assume the statement holds for all moments up to $2m-1 variables. For the inductive step to $2m variables, apply Stein's lemma (Gaussian integration by parts), which states that for a differentiable function f(\mathbf{X}),

\langle X_k f(\mathbf{X}) \rangle = \sum_{l=1}^{2m-1} C_{kl} \langle \partial_{X_l} f(\mathbf{X}) \rangle.

Setting f(\mathbf{X}) = \prod_{i=1}^{2m-1} X_i, the left side is the desired moment \langle X_1 \cdots X_{2m} \rangle (by relabeling), while the right side expands the derivative to sum over cases where \partial_{X_l} acts on one factor, yielding lower-order moments. By the inductive hypothesis, these reduce to sums over pairings of the reduced set, and collecting terms groups them into all possible pairings of the full $2m indices, each contributing the product of covariances with coefficient 1. Odd moments vanish similarly by induction from the base n=1.^[13] In the context of Wick's theorem, this combinatorial structure arises when evaluating time-ordered correlation functions in free field theory, where the time-ordering operator T effectively averages over permutations of the field arguments (assuming generic distinct times, each permutation equally likely with probability $1/(2m)!). The permutations group into equivalence classes corresponding to the pairing structures: for a fixed full contraction (all fields paired), there are $2^m m! permutations that realize it, as each pair can be ordered in 2 ways and the m pairs interleaved in m! ways; this factor cancels with the $1/(2m)! normalization, yielding the unbiased sum over pairings. Partial contractions (with normal-ordered remnants) follow analogously, though their expectations vanish in the free theory.^[12] For fermions, the argument adapts to Grassmann-valued Gaussian variables, where anticommutation introduces signs. The expectation becomes the pfaffian,

\langle \psi_1 \cdots \psi_{2m} \rangle = \sum_{\mathcal{P}} \operatorname{sgn}(\mathcal{P}) \prod_{(i,j) \in \mathcal{P}} \langle \psi_i \psi_j \rangle,

with \operatorname{sgn}(\mathcal{P}) = (-1)^{\pi(\mathcal{P})/2}, where \pi(\mathcal{P}) is the parity of the permutation sorting the indices into consecutive pairs along the pairing \mathcal{P} (equivalent to the number of crossings in the pairing diagram). This sign ensures consistency with antisymmetry: in the time-ordering, the sum over permutations includes \operatorname{sgn}(\sigma) from reordering fermionic fields, and grouping by pairings incorporates the parity via the relative ordering within and between pairs. The number of pairings remains (2m-1)!!, but the signed sum reflects the determinant structure of fermionic two-point functions.

Operator derivation

The operator derivation of Wick's theorem relies on the algebraic structure of creation and annihilation operators satisfying canonical commutation relations, providing a proof that applies uniformly to both finite and infinite collections of operators. For bosonic fields in second quantization, consider a single mode with annihilation operator a and creation operator a^\dagger obeying the commutation relation [a, a^\dagger] = 1, while [a, a] = [a^\dagger, a^\dagger] = 0. The normal-ordered product :a a^\dagger: places all creation operators to the left of annihilation operators, so a a^\dagger = :a a^\dagger: + 1. The time-ordered product for two operators, assuming the standard free-field convention where the time-ordering T enforces chronological ordering, yields T(a a^\dagger) = :a a^\dagger: + \langle 0 | a a^\dagger | 0 \rangle, with the vacuum expectation \langle 0 | a a^\dagger | 0 \rangle = 1 establishing the contraction term.^[14] To generalize to products of multiple operators, employ an inductive argument using recursive commutation relations. Assume the theorem holds for n-1 operators: T(\phi_2 \cdots \phi_n) = :\phi_2 \cdots \phi_n: + sums over all possible pairwise contractions among them, where each field \phi_j = \phi_j^{(+)} + \phi_j^{(-)} decomposes into annihilation (\phi^{(+)}) and creation (\phi^{(-)}) parts, and contractions are c-numbers given by vacuum expectations \langle 0 | T(\phi_i \phi_j) | 0 \rangle. For n operators, insert the nth field \phi_1 and use the time-ordering property: T(\phi_1 \phi_2 \cdots \phi_n) = \phi_1 T(\phi_2 \cdots \phi_n) if \phi_1 is earliest in time, or adjusted by commutation otherwise. Commuting \phi_1^{(+)} through the normal-ordered terms in the induction hypothesis generates additional contractions involving \phi_1, while \phi_1^{(-)} passes freely to the left; the full expansion sums all possible full and partial contractions, with the base case n=2 confirming the pattern. This recursive procedure, leveraging [\phi_i^{(+)}, \phi_j^{(-)}] as the sole non-vanishing commutators, yields the general form T(\phi_1 \cdots \phi_n) = \sum :\phi_1 \cdots \phi_n:_{\text{with contractions}}, where the sum includes all Wick pairings.^[14]^[15] For the infinite-dimensional case, corresponding to continuous fields in quantum field theory, the derivation extends via the generating functional approach, disentangling exponentials of operator sums. The generating functional is Z[J] = \langle 0 | T \exp\left(-i \int J(x) \phi(x) \, d^4x \right) | 0 \rangle, where \phi(x) is the field operator smeared by the source J(x). Split \int J \phi = A + B with A containing creation parts and B annihilation parts; the Baker-Campbell-Hausdorff formula expands \exp(A + B) = \exp(A) \exp(B) \exp\left(-\frac{1}{2} [A, B] + \ higher\ terms\right), but since [A, B] is a c-number (commuting with everything), higher terms vanish, yielding T \exp(A + B) = :\exp(A + B): \exp\left(\frac{1}{2} \langle 0 | [A, B] | 0 \rangle \right). The c-number exponent is -\frac{i}{2} \int J(x) \Delta_F(x-y) J(y) \, d^4x \, d^4y, with \Delta_F the Feynman propagator (contraction); functional derivatives then recover the finite-n time-ordered products as normal-ordered terms plus all pairings of propagators. In the fully infinite case, the Zassenhaus formula provides an alternative disentangling e^{A+B} = e^A e^{Z_1} e^{Z_2} \cdots, where Z_k are nested commutators, converging in the free-field Lie algebra to separate normal ordering from the exponential of pairwise contractions.^[15]^[16] Verification of equivalence to the combinatorial statement follows by expanding the operator sums: the inductive terms match the permutation count of pairings (with signs for fermions via anticommutators), as each commutation step generates unique contraction topologies without overcounting, aligning the algebraic expansion with the symmetric sum over Wick contractions.^[14]

Examples

Two-operator product

Wick's theorem provides a foundational tool for evaluating time-ordered products of operators in quantum field theory and many-body physics, with the two-operator case serving as the simplest illustration beyond trivial single-operator terms. Consider the time-ordered product T(a b) for generic operators a and b, where the time-ordering operator T arranges the operators such that earlier-time operators appear to the right.^[1] The theorem states that this product expands as

T(a b) = :a b: + \langle a b \rangle,

where :a b: is the normal-ordered product, with all creation operators placed to the left of annihilation operators (up to commutation factors for bosons or anticommutation signs for fermions), and \langle a b \rangle denotes the contraction, defined as the difference T(a b) - :a b:, which is a c-number equivalent to the vacuum expectation value \langle 0 | T(a b) | 0 \rangle since normal-ordered products vanish in the vacuum.^[1] This decomposition separates the fully connected (contracted) part from the remaining operator structure, facilitating perturbative calculations. For bosonic fields, a concrete example arises with the annihilation operator a_{\mathbf{k}} and creation operator a_{\mathbf{k}}^\dagger for a single mode \mathbf{k}, obeying the commutation relation [a_{\mathbf{k}}, a_{\mathbf{k}}^\dagger] = 1. In this case, the time-ordered product simplifies to

T(a_{\mathbf{k}} a_{\mathbf{k}}^\dagger) = :a_{\mathbf{k}} a_{\mathbf{k}}^\dagger: + 1,

where the contraction \langle a_{\mathbf{k}} a_{\mathbf{k}}^\dagger \rangle = 1 stems directly from the commutation relation and the fact that the vacuum expectation yields the identity operator's contribution. The normal-ordered term :a_{\mathbf{k}} a_{\mathbf{k}}^\dagger: = a_{\mathbf{k}}^\dagger a_{\mathbf{k}} then represents the operator without the vacuum fluctuation, highlighting how Wick's theorem accounts for the infinite zero-point energy in bosonic systems through the contraction. The fermionic variant follows a similar structure but incorporates anticommutation relations. For fermionic operators c and c^\dagger satisfying \{c, c^\dagger\} = 1, the expansion is

T(c c^\dagger) = :c c^\dagger: + \langle 0 | c c^\dagger | 0 \rangle,

with the vacuum expectation \langle 0 | c c^\dagger | 0 \rangle = 1 for an empty state, reflecting the Pauli exclusion principle's influence on correlations. Here, the normal ordering :c c^\dagger: = -c^\dagger c swaps the operators with a minus sign from anticommutation, but the overall form preserves the theorem's pairwise contraction rule, adapted for the antisymmetric nature of fermions, ensuring no overcounting in higher-order applications. This yields T(c c^\dagger) = 1 - c^\dagger c, consistent with the anticommutation relation.^[1]

Four-operator product

Wick's theorem extends naturally to the time-ordered product of four operators, T(ABCD), where A, B, C, D represent field operators at possibly distinct spacetime points. The theorem decomposes this product into a sum over all possible contractions, encompassing the fully normal-ordered term, terms with partial contractions, and fully contracted terms. For bosonic operators, the possible full contractions correspond to the three distinct ways to pair the four operators: (AB)(CD), (AC)(BD), and (AD)(BC). Partial contractions involve a single pair, such as (AB):CD:, with six such combinations arising from choosing any two operators to contract while normal-ordering the remaining pair. In the context of free bosonic fields, the contractions are evaluated as two-point propagators, denoted \langle AB \rangle or \Delta(x_A - x_B), representing the vacuum expectation value of the time-ordered product minus the normal-ordered part. The explicit expansion is thus

T(ABCD) = :ABCD: + \sum_{i<j} \langle ij \rangle :kl: + \sum_{(ij)(kl)} \langle ij \rangle \langle kl \rangle,

where the first sum runs over the six single-pair contractions (e.g., \langle AB \rangle :CD: + \langle AC \rangle :BD: + \cdots), and the second over the three double contractions (e.g., \langle AB \rangle \langle CD \rangle + \langle AC \rangle \langle BD \rangle + \langle AD \rangle \langle BC \rangle). This yields a total of 1 normal-ordered term, 6 partial contraction terms, and 3 full contraction terms, providing a combinatorial structure essential for perturbative calculations. For fermionic operators, the theorem incorporates anticommutation relations, introducing parity signs based on the number of operator interchanges required to form each contraction. Each term acquires a factor of (-1)^p, where p is the permutation parity of the pairing relative to the original order. This adjustment ensures antisymmetry, effectively reducing the number of independent positive contributions in vacuum expectations to the Pfaffian of the contraction matrix, though the full expansion retains the same 1 + 6 + 3 structural terms with modified signs (e.g., T(\psi_1 \psi_2 \bar{\psi}_3 \bar{\psi}_4) includes -S_F(x_1 - x_4) S_F(x_2 - x_3) for crossed pairings). Contractions here are fermionic propagators S_F(x - y).

Applications

Correlation functions

In quantum field theory, correlation functions are defined as the vacuum expectation values of time-ordered products of field operators, providing essential information about particle interactions and propagation. Wick's theorem is particularly powerful for computing these in the free field limit, where the vacuum expectation value of a time-ordered product of n fields factors into a sum over all possible full contractions. Specifically,

\langle 0 | T(\phi(x_1) \cdots \phi(x_n)) | 0 \rangle = \sum_{\text{all pairings}} \prod \langle 0 | T(\phi(x_i) \phi(x_j)) | 0 \rangle,

with each contraction given by the Feynman propagator \Delta_F(x_i - x_j), the two-point correlation function itself.^[4] In the free field case, this reduction simplifies the computation dramatically, as the theorem equates the n-point function to the sum of all ways to pair the fields into products of two-point propagators. For the two-point function, it directly yields the propagator:

\langle 0 | T(\phi(x) \phi(y)) | 0 \rangle = \Delta_F(x - y) = \int \frac{d^4 p}{(2\pi)^4} \frac{i}{p^2 - m^2 + i\epsilon} e^{-i p \cdot (x - y)},

which encodes the field's propagation. For the four-point function, the result includes three terms corresponding to the distinct pairings:

\langle 0 | T(\phi(x_1) \phi(x_2) \phi(x_3) \phi(x_4)) | 0 \rangle = \Delta_F(x_1 - x_2) \Delta_F(x_3 - x_4) + \Delta_F(x_1 - x_3) \Delta_F(x_2 - x_4) + \Delta_F(x_1 - x_4) \Delta_F(x_2 - x_3).

These terms represent all possible ways to fully contract the operators, yielding a Gaussian structure inherent to free theories.^[4]^[3] The full application of Wick's theorem generates both connected and disconnected contributions to the correlation functions. Disconnected parts arise from pairings that separate the fields into independent clusters, while connected diagrams require all fields to be linked through contractions. In free field theory, the connected n-point functions vanish for n > 2, reflecting the absence of interactions; the complete correlators are instead expressed via cumulants, where the generating functional for connected functions is the logarithm of the full partition function. This distinction is crucial for extracting physically meaningful scattering amplitudes from the total correlations.^[4] Wick introduced his theorem in 1950 motivated by the need to simplify perturbative calculations in quantum field theory, offering an algebraic rule to evaluate these vacuum expectation values without relying solely on diagrammatic methods, at a time when Feynman diagrams were emerging as a complementary tool.^[17]

Time-ordered products

In perturbative quantum field theory, Wick's theorem facilitates the evaluation of the Dyson series for the S-matrix or time-evolution operator in the interaction picture. The Dyson series expands the time-evolution operator as U(t) = T \exp\left( -i \int_{-\infty}^{t} H_{\text{int}}(t') \, dt' \right), where T denotes the time-ordering operator that arranges field operators in chronological order based on their time arguments. Applying Wick's theorem to the time-ordered products of interaction Hamiltonians within this expansion decomposes them into a sum of fully normal-ordered terms plus all possible partial contractions, with contractions representing free-field propagators. In the path integral formulation, Wick's theorem arises naturally when evaluating the generating functional for free fields as a Gaussian integral. The generating functional is given by Z[J] = \int \mathcal{D}\phi \, \exp\left( i S[\phi] + i \int J(x) \phi(x) \, d^4x \right), where S[\phi] is the quadratic action for free fields. Completing the square in this functional integral yields an exact Gaussian result, and the perturbative expansion in powers of the source J produces all connected and disconnected correlation functions as sums over Wick pairings of the sources, corresponding to all possible contractions via the two-point propagator. Wick's theorem applies directly to free fields but extends to interacting theories through perturbative methods combining Dyson series with Wick contractions, often termed the Dyson-Wick formalism. In this approach, the interaction Hamiltonian is treated as a perturbation around the free theory, and each order in the Dyson expansion of time-ordered products is evaluated using Wick's theorem to generate normal-ordered operators and contractions, enabling systematic computation of scattering amplitudes. This formalism offers a significant computational advantage by eliminating the need for explicit enforcement of time-ordering in multi-operator products, instead replacing it with a combinatorial sum over contractions that directly maps to Feynman diagram rules, streamlining higher-order perturbation theory calculations. Time-ordered products under Wick's theorem are essential for deriving correlation functions as vacuum expectations of these expansions.