A novel approach to perturbative calculations for a large class of interacting boson theories

TL;DR

The paper introduces a new method for perturbative calculations of the S-matrix in a broad class of boson theories, leveraging solutions to Diophantine equations, bypassing traditional Feynman diagrams, and applicable to various quantum field models.

Contribution

It presents a novel algebraic approach linking Green's functions to Diophantine equations for perturbative expansions, extending calculability to arbitrary orders and models.

Findings

Method computes S-matrix to any order for many boson theories.

Avoids traditional Feynman diagram techniques.

Demonstrated on Unruh-DeWitt detectors and connected to Hafnian enumeration.

Abstract

We present a method of calculating the interacting S-matrix to an arbitrary perturbative order for a large class of boson interaction Lagrangians. The method takes advantage of a previously unexplored link between the $n$-point Green's function and a certain system of linear Diophantine equations. By finding all nonnegative solutions of the system, the task of perturbatively expanding an interacting $S$-matrix becomes elementary for any number of interacting fields, to an arbitrary perturbative order (irrespective of whether it makes physical sense) and for a large class of scalar boson theories. The method does not rely on the position-based Feynman diagrams and promises to be extended to many perturbative models typically studied in quantum field theory. Aside from interaction field calculations we showcase our approach by expanding a pair of Unruh-DeWitt detectors coupled to Minkowski vacuum to an arbitrary perturbative order in the coupling constant. We also link our result to Hafnian as introduced by Caianiello and present a method to list all (2n-1)!! perfect matchings of a complete graph on 2n vertices.

Equations

$$⟨\Omega |\mathsf{T}\{\mathbf{\varphi }(x_1)\dots \mathbf{\varphi }(x_k)\}|\Omega ⟩,$$

$$⟨\Omega |\mathsf{T}\{\mathbf{\varphi }(x_1)\dots \mathbf{\varphi }(x_k)\}|\Omega ⟩=\frac{⟨0_M|\mathsf{T}\{\varphi (x_1)\dots \varphi (x_k)\exp [i\int {\mathrm{d}}^{d}z\text{pazocal}{L}_{\mathrm{int}}]\}|0_M⟩}{⟨0_M|\mathsf{T}\{\exp [i\int {\mathrm{d}}^{d}z\text{pazocal}{L}_{\mathrm{int}}]\}|0_M⟩},$$

$$\text{pazocal}{L}_0=\sum _{k=1}^N\frac{1}{2}{\partial }^{\nu }{\varphi }_k{\partial }_{\nu }{\varphi }_k-\frac{1}{2}{\mu }_k^2{\varphi }_k^2$$

$$\text{pazocal}{L}_{\mathrm{int}}$$

$$S^{(m)}={i^{m}\over m!}\int\mathrm{d}^{d}\boldsymbol{z}\mathop{\left<0_{M}\,\right|}\nolimits\mathsf{T}\big{\{}\phi(x_{1})\dots\phi(x_{k})\pazocal{L}_{\mathrm{int}}(z_{1})\dots\pazocal{L}_{\mathrm{int}}(z_{m})\big{\}}\mathop{|0_{M}\rangle}\nolimits$$

$$⟨0|\mathsf{T}\{{\partial }_{\mu }\varphi (x){\partial }_{\nu }\varphi (y)\}|0⟩={\partial }_{\mu }{\partial }_{\nu }⟨0_M|\mathsf{T}\{\varphi (x)\varphi (y)\}|0_M⟩-g_{\mu 0}\delta (x^{\prime }-y^{\prime })⟨0_M|\mathsf{T}\{\varphi (x){\partial }_{\nu }\varphi (y)\}|0_M⟩$$

$$⟨0|\mathsf{T}\{{\partial }_{\mu }\varphi (x){\partial }_{\nu }\varphi (y)\}|0⟩\text{mbox}‘{=}^{\prime }{\partial }_{\mu }{\partial }_{\nu }⟨0_M|\mathsf{T}\{\varphi (x)\varphi (y)\}|0_M⟩$$

$$\mathbb{E}\Big{[}\prod_{i=1}^{2m}x_{i}\Big{]}=\sum_{r=1}^{(2m-1)!!}\prod_{\genfrac{}{}{0.0pt}{2}{j,k=1}{j<k}}^{m}\mathbb{E}_{r}[x_{j}x_{k}],$$

$$\mathbb{E}\Big{[}\prod_{i=1}^{2m}x_{i}\Big{]}=\mathop{\left<0_{M}\,\right|}\nolimits\mathsf{T}\big{\{}\prod_{i=1}^{2m}\phi_{i}\big{\}}\mathop{|0_{M}\rangle}\nolimits\equiv\langle\prod_{i=1}^{2m}i\rangle_{0},$$

$$\mathop{\left<0_{M}\,\right|}\nolimits\mathsf{T}\big{\{}\prod_{i=1}^{\ell_{1}}\phi_{1}\prod_{i=1}^{\ell_{2}}\phi_{2}\dots\prod_{i=1}^{\ell_{f}}\phi_{f}\big{\}}\mathop{|0_{M}\rangle}\nolimits\equiv\langle 1^{\ell_{1}}2^{\ell_{2}}\dots f^{\ell_{f}}\rangle_{0}=\langle\underbrace{1\dots 1}_{\ell_{1}}\underbrace{2\dots 2}_{\ell_{2}}\dots\underbrace{f\dots f}_{\ell_{f}}\rangle_{0}$$

$$⟨1^{{\ell }_1}{2}^{{\ell }_2}\dots f^{{\ell }_f}{⟩}_0=\sum _{\mathbf{\alpha }}{\mu }_{\mathbf{\alpha }}\prod _{\genfrac{}{}{0ex}{}{i,j=1}{i<j}}^f⟨ij{⟩}_0^{{\alpha }_{ij}},$$

$$2{\alpha }_{ii}+\sum _{\genfrac{}{}{0ex}{}{j=1}{j\ne i}}^f{\alpha }_{ij}={\ell }_i,$$

$${\alpha }_{12}+{\alpha }_{13}+\cdots +{\alpha }_{1f}=\ell -2{\alpha }_{11}.$$

$$\sum _{k=1}^{f-1}\left(\genfrac{}{}{0ex}{}{f-1}{k}\right)\left(\genfrac{}{}{0ex}{}{\ell -1}{k-1}\right)=\left(\genfrac{}{}{0ex}{}{f-2+\ell }{f-2}\right).$$

$$\sum _{\tilde{\ell }=0}^{\ell }\left(\genfrac{}{}{0ex}{}{f-2+\tilde{\ell }}{f-2}\right)=\frac{1+\ell }{f-1}\left(\genfrac{}{}{0ex}{}{f-1+\ell }{f-2}\right)=\frac{1}{(f-1)!}\prod _{j=1}^{f-1}(\ell +j)=\frac{1}{(f-1)!}p(\ell ),$$

$$\sum _{\genfrac{}{}{0ex}{}{j=1}{j\ne i}}^f{\alpha }_{ij}\le {\ell }_i$$

$$\left(\genfrac{}{}{0ex}{}{s}{n}\right)t\times \dots \times (t-n+1)$$

$${\mu }_{\mathbf{\alpha }}=\Gamma \prod _{\genfrac{}{}{0ex}{}{i,j=1}{i<j}}^f{\gamma }_{ij},$$

$$\Gamma$$

$${\gamma }_{ij}$$

$$n$$

$$s$$

$$t$$

$${\gamma }_{ii}=\left(\genfrac{}{}{0ex}{}{{\ell }_i}{2{\alpha }_{ii}}\right)(2{\alpha }_{ii}-1)!!,$$

$$l(G)=\frac{p!}{|\mathsf{Aut}[G]|}.$$

$$V_{\mathrm{ext}}=\bigcup _{j=1}^J{V}_{\mathrm{ext}}^{(j)},$$

$${\xi }_{\mathbf{\alpha }}=\prod _{i=1}^J\left(\genfrac{}{}{0ex}{}{E-{\sum }_{j=1}^{i-1}{E}^{j-1}}{E^i}\right)\frac{I!}{|\mathsf{Aut}[G_{\mathbf{\alpha },\mathrm{int}}]|}$$

$$\int {\mathrm{d}}^4\mathbf{z}⟨0_M|\mathsf{T}\{\varphi (x_1)\varphi (x_2)\varphi (z_1{)}^4\varphi (z_2{)}^4\}|0_M⟩=\int {\mathrm{d}}^4\mathbf{z}⟨1^1{2}^1{3}^4{4}^4{⟩}_0.$$

$$⟨1^{{\ell }_1}{2}^{{\ell }_2}{3}^{{\ell }_3}{4}^{{\ell }_4}{⟩}_0\equiv ⟨\underset{{\ell }_1}{\underbrace{1\dots 1}}\underset{{\ell }_2}{\underbrace{2\dots 2}}\underset{{\ell }_3}{\underbrace{3\dots 3}}\underset{{\ell }_4}{\underbrace{4\dots 4}}{⟩}_0,$$

$$2{\alpha }_{11}+{\alpha }_{12}+{\alpha }_{13}+{\alpha }_{14}$$

$${\alpha }_{12}+2{\alpha }_{22}+{\alpha }_{23}+{\alpha }_{24}$$

Full text

1. Introduction

The calculation of an interacting $S$-matrix is one of the first problems encountered in quantum field theory (QFT) [1] . The majority of physical models must be calculated perturbatively in the coupling constant and there exists a whole calculational industry with one sole purpose: to make the increasingly tedious calculations manageable for a large variety of interaction Lagrangians [2, 3, 4, 5, 6, 7, 8]. Selected theories (such as the $\phi^{4}$ model) were studied even in more detail providing many-loop expansions in terms of the Feynman diagrams including closed expressions for their multiplicity factors [9, 10, 11]. These works fall into a broader effort of the Feynman diagram enumeration [12, 13, 14, 15, 16, 17], including generalizations to various field theories [18, 19, 20, 21, 22, 23].

The purpose of this paper is to develop a different method to construct perturbative contributions which is based on (what seems to be) an unexplored link between Wick’s theorem [24] (or its statistics equivalent due to Isserlis [25]) and a certain system of linear Diophantine equations. Based on this insight, our approach is readily applicable to a large class of boson scalar theories and have a great chance of being generalized to a much larger class of theories typically considered in QFT. It is not our goal to compete with the aforementioned multipurpose packages to perform perturbative calculations. Rather, we hope that our method will offer a new conceptual insight into perturbative calculations where a speedy calculation of scattering amplitudes and multiplicities would be a nice bonus. Indeed, from the practical point of view, the method provides an extremely streamlined and versatile way of calculating the $S$-matrix contributions including multiplicities to a high perturbative order, for any number of interacting fields and for a large class of real scalar bosonic theories without the need to summon Feynman’s diagrams (all position-based Feynman diagrams are automatically obtained in the Diophantine approach). A strategy to go beyond scalar boson models is outlined as well.

The effortless nature of the presented method is quite surprising due to the fact that solving Diophantine equations has a well-earned reputation of being hard. In more detail, it turns out that the number of ways to simplify an $n$-point Green’s function is given by all nonnegative solutions of a linear Diophantine system. This, on the other hand, is a problem equivalent to counting the number of interior points of convex polyhedra – a major topic in algebraic geometry [26, 27].

There are circumstances where the increase of perturbative components does not make sense from the physical point of view. It is well known that some theories are (super)renormalizable and others are not (depending of the spacetime dimension $d$) [1]. Additionally, infinities of a different type creep even into those theories which are renormalizable because their perturbative expansion is asymptotic [28, 29, 30]. A typical example is the quantum electrodynamics Lagrangian [1] whose perturbative contributions start to increase at the order of the inverse of the fine-structure constant. Nonetheless, we feel that the generality of the presented method and its potential to go beyond boson theories together with so far unnoticed connection (to the best author’s knowledge) of perturbative calculations to algebraic geometry/number theory makes the method interesting and we return to the potentially fruitful link between QFT perturbative calculations and lattice polyhedra in the last section.

It may seem that the number of real scalar boson theories in QFT that are worth of exploring to any order and for any number of interacting fields is limited. Even if it was true, there exists a number of physical processes where such expansions are desirable. One of them is a model of a real scalar field linearly coupled to a two-level quantum system known as the Unruh-DeWitt (UDW) particle detector. It was conceived in [31], improved [32] and studied in many physical situations [33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44]. The starting point for perturbative calculations is the Hamiltonian formalism (equivalent to the Lagrangian approach in the absence of a derivative coupling). But one quickly notices a qualitative difference compared to the perturbative expansions in QFT since the problem does not seem to suffer from the convergence problems [28, 29]. Contrary to the typical situation in QFT, it is shown that the number of perturbative terms grows polynomially with the perturbative order. It is not clear whether an actual non-perturbative solution exists (we leave it as an interesting open problem) but the main result of our calculation based on a Diophantine system of linear equations is the next best thing: an efficiently calculable expansion to an arbitrary order in the coupling constant.

Finally, we explore the connection to the so-called Hafnian introduced in [45] in the context of Dyson’s series expansion. This insight enables to apply our Diophantine algorithm to list all $(2n-1)!!$ perfect matchings of a complete graph on $2n$ vertices. Hafnians have recently attracted a lot of attention and their applications go well beyond perturbative QFT [46, 47, 48, 49].

The paper is organized as follows. After a brief recollection of the scalar boson theory in Section 2 through its Lagrangian formulation we present the main findings of this paper in Section 3. This is accompanied by several examples in Section 4 of the perturbative contributions of the $\phi^{4}$ and $\phi^{3}$ theory (for an easy comparison with the published results) and their multiplicities and one major application which is the perturbative calculation of a pair of UDW detectors to an arbitrary order in the coupling constant. In Section 5 we discuss the necessary steps in order to generalize the current formalism to more complicated perturbative models in QFT and Section 6 concludes with several open problems.

2. Boson scalar theories

An important object of study in interacting QFT is the S-matrix, which is proportional to the $n$-point interacting Green’s function

[TABLE]

where $\mathop{|\Omega\rangle}\nolimits$ is the interacting vacuum for the given theory, $\boldsymbol{\phi}(x_{i})$ are the interacting fields and $\mathsf{T}$ stands for the time-ordering operator. In the Lagrangian formulation of the theory the way of calculating Green’s function (1) is through the formula [1]

[TABLE]

where the interacting part of the complete Lagrangian $\pazocal{L}=\pazocal{L}_{0}+\pazocal{L}_{\mathrm{int}}$ appears, $\phi(x_{i})$ are free fields and $\mathop{|0_{M}\rangle}\nolimits$ is a free ground state in $d$-dimensional Minkowski spacetime (Minkowski vacuum). The type of theories we will study in this work are the $N$-mode boson scalar theories whose free Lagrangian reads

[TABLE]

and the interaction part can take many different forms. Instead of trying to write the most general $\pazocal{L}_{\mathrm{int}}$, we will list several types of interaction (omitting the multiplicative terms making the Lagrangian densities of dimension $m^{d}$):

[TABLE]

In the first line we set $N=1$, in the second item we skipped the numerical prefactors and in the third line we have $\mu_{k}=\mu$ with $(\phi_{1},\dots,\phi_{N})$ forming a $N$-component boson field exhibiting the $O(N)$ symmetry (the so-called $O(N)$ sigma model). The fields in (4) are real but the current analysis is equally applicable to complex fields with only minor modifications (and for suitable interaction Lagrangians – we discuss this and other possible generalizations in Section 5). This is due to the fact that $\phi$ and $\phi^{\dagger}$ must be considered independent. After all, the $O(2)$ sigma model can be rewritten as a one-mode complex theory whose interaction term reads $\pazocal{L}_{\mathrm{int}}=-g/4(\phi^{\dagger}\phi)^{2}$.

By looking at the RHS of (2), the solution is given by expanding the numerator around the coupling parameter. The denominator (the free $S$-matrix) is known to factor out and this removes the disconnected contributions of the scattering amplitudes in the numerator. Here comes our main contribution. We present an efficient way of calculating the perturbative expansion of the numerator of (2) for arbitrary number $k$ of interacting fields $(\phi(x_{i}))_{i=1}^{k}$, for an extensive class of interaction Lagrangians such as those in Eqs. (4) and to any perturbative order in the coupling constant. A generic expansion element in the numerator reads

[TABLE]

and the expectation value (the $k+m$-point Green’s function) is the focus of this work. Here we make the perturbative calculations extremely straightforward for any $m$ and $k$ without relying on Feynman rules or constructing Feynman position-space diagrams. We show how to efficiently factorize the $k+m$-point Green’s function for large $k,m$ into a sum of products of two-point Green’s functions (Feynman propagators) with no effort whatsoever. This is typically the most tedious step when dealing with perturbation techniques.

Lagrangians may also contain derivative interactions $\pazocal{L}_{\mathrm{int}}\equiv\pazocal{L}_{\mathrm{int}}(\phi_{i},\partial_{\nu}{\phi_{i}})$. This is more or less a technical issue despite the fact that the identity

[TABLE]

complicates the extraction of the derivatives out of the propagator. It was shown in [50] that the delta function contributions cancel out in the perturbative expansion of the whole $S$-matrix. Hence, in principle, we could use the ‘essentially’ correct identity

[TABLE]

and perform the same analysis we will present here.

3. $n$-point Green’s functions solved through a system of linear Diophantine equations

Theorem 1** (Isserlis’ [25]).**

Let $x_{i}$ be a Gaussian random variable satisfying $\mathbb{E}[\prod_{i=1}^{2m+1}x_{i}]=0$. Then

[TABLE]

where the sum goes over the products of bivariate expectation values $\mathbb{E}_{r}$.

Remark* (notational).*

In our case the Gaussian random variable will be a product of time-ordered free scalar fields $\phi(x_{i})$ and $\phi(z_{j})$ and the expectation value will be taken w.r.t. to Minkowski vacuum $\mathop{|0_{M}\rangle}\nolimits$:

[TABLE]

where $\phi_{i}\equiv\phi(x_{i})$ (or $\phi(z_{i})$) and on the RHS is the minimalist notation we will be mostly using.

Remark*.*

The theorem is also known as Wick’s theorem [24]. Wick’s theorem transforms time-ordered operator expressions into a normal form [1] and Isserlis’ result is recovered upon taking a (free) vacuum expectation value. As a matter of fact, Isserlis’ theorem is more general since naturally there is no notion of time-ordering in Eq. (8) and so the theorem applies even for ‘ordinary’ products of free scalar fields. Of course, we are actually calculating the time ordered version as it appears in the definition of the sought Green’s function.

For $2m$ different random variables $x_{i}$ (scalar fields $\phi_{i}$) in Eq. (9) there is nothing else to say but Isserlis’ theorem can be refined if the number of different variables $x_{i}$ in (8) is less than $2m$.

Definition 1* ([27]).*

A $k$-composition of $n\in\mathbb{Z}_{>0}$ is $n$ written as an ordered sum of $k$ strictly positive integers.

Remark*.*

Order matters unlike for integer partitions. There are $\binom{n-1}{k-1}$ $k$-compositions of $n$. For example, for $n=4$, there are three 2-compositions: $1+3,2+2$ and $3+1$.

Definition 2*.*

Given a finite set of positive integers $S=(1,2,\dots,f)$ we define the lexicographic ordering on the subset of two elements of $i,j\in S$ as $ij\leq i^{\prime}j^{\prime}$ iff $i\leq i^{\prime}$ or $i=i^{\prime}$ together with $j\leq j^{\prime}$.

Theorem 2**.**

Let $\ell_{i}\in\mathbb{Z}_{\geq 0}$ such that $\sum_{i=1}^{f}\ell_{i}=2m$. Then, Green’s function

[TABLE]

can be written as a product of two-point correlation functions

[TABLE]

where $\mu_{\boldsymbol{\alpha}}\in\mathbb{Z}_{>0}$ is the product multiplicity, $\alpha_{ij}\in\mathbb{Z}_{\geq 0}$ and $\boldsymbol{\alpha}\overset{\mathrm{df}}{=}(\alpha_{ij})_{1\leq i\leq j\leq f}$. The number of products is a polynomial function of $\ell_{i}$ for a fixed number of fields $f$.

Definition 3*.*

A product $\prod_{\genfrac{}{}{0.0pt}{2}{i,j=1}{i<j}}^{f}\langle ij\rangle^{\alpha_{ij}}_{0}$ is called a conformation and $\boldsymbol{\alpha}$ a conformation exponent.

Remark*.*

The conformation exponent $\boldsymbol{\alpha}$ will always be lexicographically ordered but it is not necessary for the proof of the above theorem.

Proof of Theorem 2.

To solve Eq. (11) means to find all nonnegative solutions $\alpha_{ij}$ of the system

[TABLE]

where $1\leq i\leq f$. It is a system of linear Diophantine equations and they have zero, one or infinitely many solutions (counting both positive and negative ones). In this case, the numerical coefficients for the variables $\alpha_{ij}$ are simple so we can easily guess a solution for $i=1$ in (12) (say $\alpha_{12}=\alpha_{13}=\ldots=\alpha_{1\ f-1}=0$ and so $\alpha_{1f}=\ell_{1}-2\alpha_{11}$ which has infinitely many solutions). It also means that there are finitely many nonnegative solutions $\alpha_{ij}$ we are looking for. Let us prove the second statement first and polynomially upper bound the number of non-negative solutions. By setting $\ell=\max_{i}{\ell_{i}}$, the first row of (12) becomes

[TABLE]

Let us also set $a_{11}=0$ for the moment. The task of finding the number of nonnegative solutions resembles the problem of finding the number of $k$-compositions of $\ell$ (see Definition 1). Indeed, if we take $1\leq k\leq f-1$ then we only need to know how to distribute $\binom{\ell-1}{k-1}$ $k$-composition in $(f-1)$ ‘slots’. Or, put differently, how many ways we can pad a $k$-composition with zeros. This equals to $\binom{f-1}{k}\binom{\ell-1}{k-1}$ and by summing over all $k$-compositions we get

[TABLE]

We have to add all possible $\ell-2\alpha_{11}$ on the RHS of (13). Since $\alpha_{11}\in\mathbb{Z}_{\leq 0}$, the parameter $\ell$ decreases by the multiples of two so we will need to distinguish between odd and even $\ell$ and sum only the appropriate ones. But we are interested in an upper bound so let’s pretend $\ell$ can be any positive integer and just calculate

[TABLE]

where $p(\ell)$ is a polynomial of degree $f-1$. There is $f$ linear equations in system (12). Each consecutive equation has one independent variable less than the previous one but if we assume for a while that all $f$ equations contain $f$ independent variables then the number of solutions is upper bounded by a polynomial of degree $(f-1)^{f}$.

How do we obtain all nonnegative solutions in a systematic way? Since $0\leq\alpha_{ii}\leq\lfloor\ell_{i}/2\rfloor$ and $0\leq\alpha_{ij}\leq\min{\{\ell_{i},\ell_{j}\}}$ for $i\neq j$ we start in the first equation of system (12) ($i=1$) by fixing the lowest possible values of $\alpha_{11},\alpha_{12}$ up to $\alpha_{1\ f-2}$ ($\alpha_{11}=\alpha_{12}=\ldots\alpha_{1\ f-2}=0$) and simply list all the admissible $\alpha_{1j}$’s for $f-2<j\leq f$. We continue by increasing $\alpha_{1\ f-2}$ by one and repeat the procedure until we hit $\min{[\ell_{1},\ell_{2}]}$. Then we repeat the whole process for $\alpha_{1\ f-3}$ all the way to $\alpha_{11}$. In the next step we move to (12) for $i=2$ by inserting all found $f$-tuples $(\alpha_{1j})_{1\leq j\leq f}$ and repeat the procedure until we find all admissible solutions in the second row. As a result we obtain an $(2f-1)$-tuple $(\alpha_{1j},\alpha_{2j^{\prime}})_{1\leq j\leq f,2\leq j^{\prime}\leq f}$. As it is clear from (12) the second row ‘inherits’ one variable ($\alpha_{12}$) from the first row. We continue in a similar fashion for all $i$. The result is a complete list of $f(f+1)/2$-tuples $\boldsymbol{\alpha}$ thus solving Eqs. (12). ∎

Remark*.*

Note that the number of nonnegative solutions is gruesomely overestimated. One can convert the second part of the proof of Theorem 2 into a program that systematically finds all the solutions. Essentially, if it takes one time step to find the first solution then all the $t$ solutions can be find in $t$ time steps where $t$ increases polynomially. There is really no need for solving a Diophantine system by some sophisticated number-theoretic methods since due to the simple form of (12) we get all nonnegative solutions by inspection and only need to list them (i.e., save them to the memory).

A related problem is how to get a closed form for the number of nonnegative solutions for given $\ell_{i}$ and $f$ without actually listing and counting the solutions. This is a nontrivial task attempted by the author in a very special case of $\ell_{i}=\ell$ and $f=4$ [51] (the chosen values have no relevance to perturbative QFT). It is a problem studied by the Ehrhart theory [52] where the interior points of convex lattice polyhedra are counted [26]. The Ehrhart theory does not provide a closed expression for the number of lattice points, however, and, typically, computer algebra systems are used to count the interior points. The connection to convex geometry comes from rewriting system (12) as a system of inequalities

[TABLE]

which is an algebraic definition of a convex polyhedron.

To get the conformation multiplicity factor $\mu_{\boldsymbol{\alpha}}$ in Theorem 2 we present an auxiliary lemma.

Lemma 3**.**

Let $S$ and $T$ be discrete sets where $|S|=s,|T|=t$. Then there is

[TABLE]

bijective functions $f:X\mapsto Y$ where $X\subset S,Y\subset T$ and $|X|=|Y|=n$ where $0<n\leq\min{\{s,t\}}$. For $n=0$ we set Eq. (16) to one.

Remark*.*

If $s=t=n$ then (16) becomes $n!$ which is known to be the number of bijections from a set to itself (the number of permutations). To make the notation more concise in following text we will use the definition of the falling factorial: $[m]_{n}\overset{\mathrm{df}}{=}m\times\ldots\times(m-n+1)$. So Eq. (16) becomes $\binom{s}{n}[t]_{n}$.

Proof.

The coefficient $\binom{s}{n}$ is simply a number of all possible domains $X\subset S$ whose cardinality is $n$. Then for every domain $X$ there is $\binom{t}{n}n!$ codomains $Y\subset T$ of the same cardinality. The coefficient $n!$ comes from the number of permutations within each of $\binom{t}{n}$ codomains $Y$. ∎

Theorem 4**.**

The multiplicity factor $\mu_{\boldsymbol{\alpha}}$ for a given conformation exponent $\boldsymbol{\alpha}$ reads

[TABLE]

where

[TABLE]

where $[\bullet]_{\alpha_{ij}}$ denotes the falling factorial introduced in the remark below Lemma 3.

Proof.

We will use repeatedly Lemma 3 by identifying $\langle ij\rangle_{0}^{\alpha_{ij}}$ for $i\neq j$ from Theorem 2 by setting:

[TABLE]

Set $S$ is the set of all $i$’s and $T$ is the set of all $j$’s. However, if we calculate more than one two-point Green’s function (like in Theorem 2) we have to take into account the fact that the sets $S$ and $T$ might have shrunk. This depends on whether in the preceding Green’s function $\langle kl\rangle_{0}^{\alpha_{kl}}$ we had $k=i$ or $l=j$. The strategy we will follow here is to first count the multiplicity of $\langle ii\rangle_{0}^{\alpha_{ii}}$ as they are independent (meaning non-overlapping for different $i$’s) and then the multiplicities of $\langle ij\rangle_{0}^{\alpha_{ij}}$. In that case the derivation follows the lexicographic ordering, Def. 2, which is a handy tool here.

The multiplicity of $\langle ii\rangle_{0}^{\alpha_{ii}}$ immediately follows from Isserlis’ theorem. Looking at Eq. (8) we see that if $x_{i}=x_{j}$ there will be $(2m-1)!!$ identical products. Hence ($m$ is $\alpha_{ii}$ here)

[TABLE]

as follows from Eq. (16)111Indeed, the double factorial follows either from Theorem 1 as stated or from the counting done in Lemma 3.. The factor of two in the binomial ‘denominator’ accounts for the two $i$’s in $\langle ii\rangle_{0}$. Repeating this procedure for all $i$ we get $\Gamma$ in Eq. (18a). To get $\gamma_{ij}$ we then have to keep track of the set cardinality and this is the point where the used lexicographic order becomes important. We follow the ordered set of $ij$ where $1\leq i<j\leq f$ (note the sharp inequality). Hence, $ij=(12,\dots,1f,23,\dots,ff)$ and we find all $\gamma_{ij}$ where $i\neq j$. The first one is $\gamma_{12}$ and so we set $n=\alpha_{12}$ (using the notation of Lemma 3) and notice that the cardinality of the discrete set of ‘ones’ is decreased by those contributing to the previously analyzed case $i=j$, namely $\gamma_{11}$. So $s=\ell_{1}-2\alpha_{11}$. Similarly, the cardinality of the ‘twos’ is decreased by the number of elements contributing to $\gamma_{22}$. So $t=\ell_{2}-2\alpha_{22}$ and $\gamma_{12}=\binom{\ell_{1}-2\alpha_{11}}{\alpha_{12}}[\ell_{2}-2\alpha_{22}]_{\alpha_{12}}$ follows. At this point we can proceed and construct a generic $\gamma_{ij}$. The ‘source pool’ of $i$’s will be depleted by three contributions: $2\alpha_{ii},\sum_{m=1}^{i-1}\nolimits\alpha_{mi}$ and $\sum_{n=i+1}^{j-1}\nolimits\alpha_{in}$. The bounds of the sums are chosen such that only the preceding contributions to the sought $\gamma_{ij}$ are accounted for. Similarly, the ‘target pool’ of $j$’s is depleted by $2\alpha_{jj}$ and $\sum_{m=1}^{i-1}\alpha_{mj}$, but the sum bounds only choose those contributions preceding $\gamma_{ij}$ (using the lexicographic ordering). We notice an interesting asymmetry, where in the falling factorial part of (18b) there is only one sum. This is because all (negative) contributions $\gamma_{jl}$ come after any $\gamma_{ij}$ given the lexicographic ordering and $\gamma_{jj}$ was taken care of separately. ∎

Remark*.*

We again emphasize the importance of the lexicographic ordering by which we followed the construction of $\gamma_{ij}$’s. Other orderings are certainly plausible and have to provide the same $\mu_{\boldsymbol{\alpha}}$, Eq. (17). But we believe the lexicographic ordering and the construction based on it is the most intuitive one.

Having a decomposition into a product of two-point Green’s functions obtained through Theorem 2 we can proceed as usual. The position-based Feynman diagrams can be immediately read off from the products of Green’s functions. However, even different Diophantine solutions (different Green’s products) lead to the same Feynman diagram, that is, to the same physical processes. This is typical for real scalar theories since when we permute the internal degrees of freedom it is often the same physical process (the same integral in the $S$-matrix). However, when charges are involved, like for complex fields, the internal lines are directed and we may need to distinguish among the Diophantine solutions. For more, see Sections 4 and 5. At this point it is advantageous to express our results in the language of graph theory which leads to Feynman diagrams. This is just an advantageous way of presenting the results and not a starting point as it typically is in QFT. The connection between Feynman diagrams and graph theory is well-documented [53] and here we recall a few basic concepts and results [54, 55].

Definition 4*.*

A labeled undirected graph $G=(V,E)$ is a set of vertices $V$ and a set $E$ of unordered pairs of elements of $V$ called edges. A graph is called simple if no loops and no multiple edges connecting the same pair of vertices are allowed. Graphs allowing both loops and multiple edges are called pseudographs.

The following definition of graph isomorphism is split into two: the standard one [54] and a version upgraded by a requirement on the labeling, see [55, Chapter 1]. We will need both definitions.

Definition 5*.*

A type-1 graph isomorphism ${}_{1}(G_{1},G_{2})$ is a bijection ${}_{1}:V_{1}\mapsto V_{2}$ such that ${}_{1}(E_{1})=E_{2}$. A type-2 graph isomorphism ${}_{2}(G_{1},G_{2})$ is a bijection ${}_{2}:V_{1}\mapsto V_{2}$ such that ${}_{2}(E_{1})=E_{2}$ while also preserving the labeling. An automorphism of a labeled graph $G$ is defined as $\euA(G)={}_{2}(G,G)$. The set of all automorphisms of $G$ forms the group $\mathsf{Aut}{[G]}$.

If we restrict ourselves to boson scalar theories then a position-based Feynman diagram is a labeled undirected pseudograph $G=(V,E)$ (we will omit the prefix pseudo- and all other adjectives from now on unless we need them to avoid confusion). The main observation here is that the conformation exponent $\boldsymbol{\alpha}=(\alpha_{ij})_{1\leq i\leq j\leq f}$ from Theorem 2 generates $G$ where $|V|=f$ and $0<\alpha_{ij}\in E$. Let us summarize in

Definition 6*.*

We call $G_{\boldsymbol{\alpha}}$ an $\boldsymbol{\alpha}$-generated labeled undirected graph if every $\alpha_{ij}\neq 0$ (from Theorem 2) denotes $\alpha_{ij}$ edges connecting two vertices $i$ and $j$. We split the edges and vertices into external and internal ones according to the details of the studied boson theory, calculated perturbative order and the number of interacting fields and write $G_{\boldsymbol{\alpha}}=(V_{\mathrm{int}}\cup V_{\mathrm{ext}},E_{\mathrm{int}}\cup E_{\mathrm{ext}})$. We further define the internal subgraph of $G_{\boldsymbol{\alpha}}$ as $G_{\boldsymbol{\alpha},\mathrm{int}}\overset{\mathrm{df}}{=}(V_{\mathrm{int}},E_{\mathrm{int}})$.

The internal vertices are the internal degrees of freedom that are integrated over in Lagrangian (2) ($\int\mathrm{d}^{d}\boldsymbol{z}$). The external vertices represent the interacting fields ($(x_{i})_{i=1}^{k}$) and are connected to the internal vertices by external edges. The internal edges connect the internal vertices. Finally, only connected graphs are considered.

Remark*.*

A symmetric adjacency matrix is a way of encoding an undirected graph. If we take an upper-triangular part of the matrix and flatten it row by row we get our conformation exponent $\boldsymbol{\alpha}$.

Theorem 5**.**

[56]** The number of labelings of a given graph $G$ with $p$ vertices is

[TABLE]

Theorem 5 answers the question of how many labeled type-1 isomorphic graphs to $G$ there are.

Lemma 6**.**

Let $G_{\boldsymbol{\alpha}}$ be an $\boldsymbol{\alpha}$-generated labeled undirected graph where $I\overset{\mathrm{df}}{=}|V_{\mathrm{int}}|,E\overset{\mathrm{df}}{=}|V_{\mathrm{ext}}|$ and $(z_{i})_{i=1}^{I}\in V_{\mathrm{int}}$ together with $(x_{i})_{i=1}^{E}\in V_{\mathrm{ext}}$. Given a partition of $V_{\mathrm{ext}}$ into $J\leq I$ disjoint sets $V_{\mathrm{ext}}^{(j)}$ of cardinality $E^{j}\overset{\mathrm{df}}{=}|V_{\mathrm{ext}}^{(j)}|$

[TABLE]

there is

[TABLE]

Feynman diagrams representing the same physical interaction. We set $E^{0}=0$.

Remark* (important).*

Let us emphasize what we do and do not count in this lemma. In Fig. 1 (a), there is an example of a simple Feynman diagram with 4! ways of permuting the external vertices. But at the level of a system of Diophantine equations, all 4! possibilities corresponds to a single nonnegative solution and its multiplicities were calculated in Theorem 4. In Fig. 1 (b), there is a generic $\boldsymbol{\alpha}$-generated graph where the grey circle represents an arbitrary relation among the six internal vertices $z_{i}$. There can be multiple edges and loops and internal vertices need not be connected to any external vertex. The external vertices are labeled by $x_{i}$ and a swap of $x_{i}$ with different $z_{i}$ corresponds to a different solution of a Diophantine system. All such possibilities are then enumerated by $\xi_{\boldsymbol{\alpha}}$ in (24).

Proof of Lemma 6.

The graph $G_{\boldsymbol{\alpha}}$ is $\boldsymbol{\alpha}$-generated and so we assume that if an internal vertex connects to one or more external vertices their number is fixed. Then, the product of binomials $\prod_{i=1}^{J}\binom{E-\sum_{j=1}^{i-1}E^{j-1}}{E^{i}}$ counts the number of ways the external vertices can be connected to the internal ones. The $i$-th binomial ‘numerator’ is depleted by the external vertices already connected to the internal ones and there are naturally many ways (depending on the order of $E^{j}$) leading to the same result. The coefficient $I!$ in (24) comes from Theorem 5 and it is a permutation of all internal vertices. We now introduce two procedures named amputation (not to be confused with the equally called procedure for removing infinities in the momentum-based Feynman diagrams!) and grafting in order to determine $|\mathsf{Aut}{[G_{\boldsymbol{\alpha},\mathrm{int}}]}|$. For a given $G_{\boldsymbol{\alpha}}$ we first remove all the external vertices but keep the external edges ‘freely floating’. This is called amputation and at this point it is not a graph. We make it a graph again by creating loops of out the amputated external edges. This step will be called grafting and the whole procedure is illustrated in Fig. 2. The newly created loops hold the information about the internal subgraph symmetry and thus its automorphism group (using Definition 5 type-2 isomorphism). Then, ${I!\over|\mathsf{Aut}{[G_{\boldsymbol{\alpha},\mathrm{int}}]}|}$ is the total number of internal graphs with different $\boldsymbol{\alpha}$ corresponding to different solutions of the Diophantine system leading to the same graph. ∎

Remark*.*

Note that a grafted graph is not a vacuum Feynman diagram. This can be seen in Fig. 2 on the right where some vertices have different vertex orders (the number of edges meeting there). It is just a pseudograph we have created to help us calculate the overall multiplicity factor of the graph we started with. Even though the graph automorphism order is always studied when dealing with Feynman diagrams (explicitly like in [57] or implicitly like everywhere else), our approach invokes it at a different stage of the calculations and for grafted graphs compared to other studies.

How do we calculate $\xi_{\boldsymbol{\alpha}}$? There are two different routes. We can follow Theorem 5 and that includes the calculation of $|\mathsf{Aut}{[G_{\boldsymbol{\alpha},\mathrm{int}}]}|$. For moderate graphs, it can be done by realizing that the whole list of Diophantine solutions $\boldsymbol{\alpha}$ splits into several classes of physically indistinguishable processes (the same Feynman diagrams). We take any representative of the class we are interested in, process it according to the previous lemma and calculate $|\mathsf{Aut}{[G_{\boldsymbol{\alpha},\mathrm{int}}]}|$. The problem of graph automorphism is not known to be tractable and is closely related to the well-known graph isomorphism membership problem [58, 54]. However, practically efficient algorithms are known and implemented (for pseudographs and discrete mathematics in general, SAGE [59] is a great tool). Even though almost all finite graphs possess no symmetry (that it, a trivial $\mathsf{Aut}{[G_{\boldsymbol{\alpha},\mathrm{int}}]}$, as proved in [60] for simple graphs), we need the exact counting.

To check the graph automorphism order we can take the second route to get $\xi_{\boldsymbol{\alpha}}$. By solving a Diophantine system we get all admissible conformation exponents $\boldsymbol{\alpha}$ and within each class there is obviously never the same $\boldsymbol{\alpha}$ more than once. In other words, the graph automorphism order is already accounted for and $\xi_{\boldsymbol{\alpha}}$ is nothing else than the cardinality of the class we are interested in. We did not get away with the hardness of the problem, though. In order to compute the cardinality we have to take a representative and check how many graphs from the list of all $\boldsymbol{\alpha}$’s it is isomorphic to (using Definition 5, type-1 isomorphism). This approach is most likely harder from the computational complexity viewpoint because the graph isomorphism membership algorithm has to be invoked many times and it is computationally equivalent to the graph automorphism order [58]. On the other hand, recall that the Diophantine system of solutions grows polynomially and, mainly, many solutions are computationally easy to disprove to be isomorphic to the studied representative by checking whether the graph is (dis)connected, by counting the edges or checking any other graph invariant [54] that is computationally cheap to implement on the level of the conformation exponent $\boldsymbol{\alpha}$ (to name a few more it is the number of loops, the ways the external edges are connected to the internal vertices etc.). But to the author’s knowledge there is no guarantee that two different classes always differ at least in one invariant (for any bosonic theory). After finding $\xi_{\boldsymbol{\alpha}}$, we can immediately obtain the graph automorphism order by counting $E,E^{j}$ and $I$ and inserting it into Eq. (24). Both routes will be illustrated in the next section.

Our hope is that the structure of a convex lattice polyhedron where all nonnegative solutions of a Diophantine system live will provide yet another way of counting the graph automorphism order, this time without the need to run computationally costly algorithms. We will briefly return to this open problem in Section 6.

4. Applications and examples

Several explicit examples are worked out in the next section. One can promptly (that is, in linear time in the number of Diophantine solutions) processed the solutions of the Diophantine system and keep only the diagrams of interest such as all connected diagrams or all vacuum diagrams and so on. It is well-known that the vacuum contributions in the numerator of Eq. (2) cancel with the free $S$-matrix in the denominator [1] and the cancellation can be done at the level of Diophantine solutions. Again, we illustrate it in the next section.

Example 1* ($\phi^{4}$ for $k=2$, second perturbative order in detail).*

Let’s calculate the third term from the expanded numerator of Eq. (2). It is proportional to

[TABLE]

Upon setting $f=4$, Eq. (10) becomes

[TABLE]

Next, we solve (12) that took the following form:

[TABLE]

where $\ell_{1}=\ell_{3}=1$ and $\ell_{3}=\ell_{4}=4$. The result is eleven conformation exponents

[TABLE]

whose lexicographic ordering is

[TABLE]

The solutions and the values of $\ell_{i}$ are inserted into

[TABLE]

to get Eq. (17) whose explicit form reads

[TABLE]

We obtain

[TABLE]

and perform a highly non-trivial check of both Theorems 2 and 4:

[TABLE]

We obtained eleven diagrams but not all are physically distinguishable. To account for them we use Theorem 4 (valid for connected diagrams which are physically most interesting). The unique diagrams and their conformations are depicted in Fig. 3. In the first row we put all three connected diagrams (the first three rows in (28)) with multiplicity two which we will now verify. $E=2$ and $I=2$ hold for all diagrams. For the left upper diagram we find $E^{1}=2$ while for the other two we get $E^{1}=E^{2}=1$. Then, following the procedure of amputation and grafting, we create new graphs where the upper left one has no symmetry (and so $|\mathsf{Aut}{[G_{\boldsymbol{\alpha}_{1},\mathrm{int}}]}|=1$) while the other two are identical upon swapping $3\leftrightharpoons 4$ (hence $|\mathsf{Aut}{[G_{\boldsymbol{\alpha}_{2,3},\mathrm{int}}]}|=2$). So we get

[TABLE]

and indeed $\boldsymbol{\alpha}_{4,5,6}$ in (28) are the other class members for $\boldsymbol{\alpha}_{1,2,3}$ (in this order). We just explored the two routes of calculating $\xi_{\boldsymbol{\alpha}}$ as described below Lemma 6. A final check is the overall multiplicity of a given Feynman diagram: $\nu_{\boldsymbol{\alpha}_{i}}=\mu_{\boldsymbol{\alpha}_{i}}\xi_{\boldsymbol{\alpha}_{i}}=(288,192,288)$ for $i=1,2,3$ in accord with [9, Table II]. ∎

Remark*.*

By solving (27) for $\ell_{i}=2$ we get 17 conformations $\langle 1^{2}2^{2}3^{2}4^{2}\rangle_{0}$ and find

[TABLE]

from (30). A calculation performed in [51] leads to a closed expression for the number of non-negative solutions of (27) (for $\ell_{i}=\ell$ even) to be

[TABLE]

For $\ell=2$ we get $\mathsf{e}(2)=17$ in accord with the above result. In theory, one can get closed expressions for the number of nonnegative Diophantine solutions following the methods of Ehrhart theory [52, 26].

Example 2* ($\phi^{4}$ for $k=4$, third perturbative order).*

We study

[TABLE]

By solving the corresponding Diophantine system we obtain 960 nonnegative solutions whose conformation exponents will not be listed. After filtering out the disconnect graphs we find in total 8 classes of connected graphs just like in [9, Table III, $p=3$]. Our results completely agree with the multiplicities of all 8 Feynman diagrams. Let’s provide two class examples whose representatives are

[TABLE]

where $\mu_{\boldsymbol{\alpha}_{1}}=6912,\mu_{\boldsymbol{\alpha}_{2}}=3456$. For $\boldsymbol{\alpha}_{1}$ we find $E=4,E^{1}=E^{2}=2$ and $I=3$. Following Lemma 6 we get a grafted graph $G_{\boldsymbol{\alpha}_{1},\mathrm{int}}$ where $|\mathsf{Aut}{[G_{\boldsymbol{\alpha}_{1},\mathrm{int}}]}|=2$ leading to

[TABLE]

Thus $\nu_{\boldsymbol{\alpha}_{1}}=\mu_{\boldsymbol{\alpha}_{1}}\xi_{\boldsymbol{\alpha}_{1}}=124416$ as it should be. In the second case we have $E=4,E^{1}=3,E^{2}=1$ and $I=3,|\mathsf{Aut}{[G_{\boldsymbol{\alpha}_{2},\mathrm{int}}]}|=1$ and so

[TABLE]

Therefore, $\nu_{\boldsymbol{\alpha}_{2}}=\mu_{\boldsymbol{\alpha}_{2}}\xi_{\boldsymbol{\alpha}_{2}}=82944$ again in agreement with [9]. ∎

Example 3* ($\phi^{3}$ for $k=2$, second perturbative order plus vacuum amplitudes).*

Here we analyze

[TABLE]

Similarly, the second expansion term of the denominator of (2) (the vacuum digrams) is proportional to

[TABLE]

Using Theorem 2 we find 8 conformations for Eq. (40)and Theorem 4 yields their multiplicities

[TABLE]

We check the solution by inspecting

[TABLE]

where the RHS follows from $\ell_{1}=\ell_{2}=1$ and $\ell_{3}=\ell_{4}=3$. By considering just connected diagrams, $(\boldsymbol{\alpha}_{1},\boldsymbol{\alpha}_{8})$ and $(\boldsymbol{\alpha}_{4},\boldsymbol{\alpha}_{6})$ represent the same processes as confirmed by Lemma 6 after we do the proper counting, grafting and graph automorphism counting:

[TABLE]

Hence $\nu_{\boldsymbol{\alpha}_{1}}=\nu_{\boldsymbol{\alpha}_{4}}=36$ as we can find in [1, Chapter 7].

For (41) we get only two solutions

[TABLE]

with multiplicities $\mu_{\boldsymbol{\alpha}_{c}}=(6,9)$. As these are vacuum graphs, there is no need to remove the external edges and graft and so following Lemma 6 we find that $\xi_{\boldsymbol{\alpha}_{1}}=\xi_{\boldsymbol{\alpha}_{2}}=2!/2=1$. The final answer is then $\nu_{\boldsymbol{\alpha}_{1}}=6,\nu_{\boldsymbol{\alpha}_{2}}=9$ in agreement with [1, Chapter 7]. ∎

A pair of Unruh-DeWitt detectors

A realistic Unruh-DeWitt detector [31, 32] is described by the following interaction Hamiltonian

[TABLE]

and

[TABLE]

is a real scalar field where $(\omega_{k},k)$ is a four-momentum, $w(\tau)$ is a detector’s window function ($\tau$ is its proper time and $t(\tau),x(\tau)$ are Minkowski coordinates), $f(x)$ a smearing function, $\sigma^{\pm}$ are the detector ladder operators, $\delta$ stands for the energy gap and $\lambda$ for a coupling constant. The evolution operator for observer $A$ reads

[TABLE]

Out task is to efficiently calculate the following unitary (tensor product and time-ordering implied)

[TABLE]

for any $n$, no matter how large, and for any input atomic state. We denoted

[TABLE]

and so $A_{-}=A_{+}^{\dagger},B_{-}=B_{+}^{\dagger}$ holds. Note that we set $f(x)$ to be a delta function for convenience. The result of this calculation is oblivious to the details of smearing, window functions or trajectories. The expression that matters is Eq. (4) together with the fact fact that $A_{\pm},B_{\pm}$ are boson scalar fields. Our final goal is to express $2n$-point Green’s functions in terms of a polynomial number of two-point correlators and not their calculation where these details are relevant.

The matrix elements of (4) will be written as

[TABLE]

where $i,j,k,l\in\{0,1\}$ denote the canonical basis of the detectors’ initial and final state and a tensor product is implied. In all previous examples we calculated a transition amplitude $\mathop{\left<f\,\right|}\nolimits S^{(m)}\mathop{|i\rangle}\nolimits$ for a perturbatively evolved $S$-matrix. This is a typical task in QFT. Here we are after a different (but, of course, related) quantity: the transition probability. More precisely, we want to perturbatively expand density matrix (52) where the probabilities lie on the diagonal. We will derive $\omega_{AB}$ for $i,j=0$ (that is, assuming the detectors are in their ground states) and later argue that a trivial modification of our analysis will enable us to calculate Eq. (52) for any canonical basis state $\mathop{|ij\rangle}\nolimits$. This, on the other hand, will lead to a ‘standalone’ unitary matrix $U_{A}\otimes U_{B}$ expressed perturbatively which then can be used to find $\omega_{AB}$ for any detector input state and thus completely solving the problem.

We start by taming the exponential number (in $n$) of the summands of the core expression, Eq. (4). By plugging the result into Eq. 52 we get again a polynomial number of $2n$-point Green’s functions and apply the results obtained in this paper to express them in terms of the two-point Green’s functions. Let us recall the algebra of Pauli operators representing the two-level detector(s). They satisfy $(\sigma^{+})^{2}=(\sigma^{-})^{2}=0$ and also $\sigma_{A}^{\pm}\sigma_{B}^{\pm}=\sigma_{B}^{\pm}\sigma_{A}^{\pm}$ together with $\sigma_{A}^{\pm}\sigma_{B}^{\mp}=\sigma_{B}^{\mp}\sigma_{A}^{\pm}$. That is, the atomic operators for the two detectors commute in the Unruh-DeWitt model.

Proposition 7**.**

For $n=2m$ we obtain from Eq. (4)

[TABLE]

and for $n=2m+1$ we get

[TABLE]

Proof.

The initial state $\mathop{|0\rangle}\nolimits$ is a ground state: $\sigma^{-}\mathop{|0\rangle}\nolimits=0$. There are two possibilities for the output states: $\mathop{|0\rangle}\nolimits$ and $\mathop{|1\rangle}\nolimits$ but taking into account the nilpotency property of $\sigma^{\pm}$ we can also enumerate all possible ways the output states can be obtained. It is simply

[TABLE]

for all $p>0$. All other operator sequences result in zero. In our case we have two sets of ladder operators – for the $A$ and $B$ Hilbert spaces and so

[TABLE]

In the first two lines $n$ is even and in the last two lines $n$ is odd. The sums should be understood as counting all possible strings of the sigma operators leading to a given output state on the left and we will now find the explicit expressions. Let’s call them nonzero strings. Given the desired output state, the main observation is that there are always only two sigma operators that do not destroy the previous nonzero string. This is because every nonzero string acting on $\mathop{|00\rangle}\nolimits$ can land in only one of the four possible states. For instance, if the state is $\mathop{|01\rangle}\nolimits$ the next sigma operator can be $\sigma^{+}_{A}$ or $\sigma_{B}^{-}$ to get another nonzero string. Given the offspring node, the choice of either $+$ or $-$ in the branching is unambiguous in order to get a nonzero string. Similarly for the other three states occupying the node and so we can represent all nonzero strings as a complete binary tree of the depth $n$ where we adopt the following convention: The root node is the initial state $\mathop{|00\rangle}\nolimits$ where the left offspring corresponds to the action of $\sigma_{A}^{+}$ and the right offspring corresponds to $\sigma_{B}^{+}$. The nodes are labeled by the resulting state $\mathop{|ij\rangle}\nolimits$, in this case $\mathop{|01\rangle}\nolimits$ or $\mathop{|10\rangle}\nolimits$, see Fig. 4. In all the branchings that follow the left offsprings correspond to the action of $\sigma_{A}^{\pm}$ and the right offsprings correspond to $\sigma_{B}^{\pm}$.

But this immediately provides the desired counting because we have just shown a bijection between the number of branches of a perfect binary tree of the depth $n$ and the number of non-zero strings of the same length. This is because each branch is uniquely represented by a nonzero string of the length $n$ and no string out of $2^{n}$ of them is missing. Consider Eq. (57a). We set the total number of branchings $n=2m$ and $p=\ell$ giving us $q=m-\ell$. The number of paths with $2\ell$ left branchings (one $\ell$ for $\sigma_{A}^{+}$ and one $\ell$ for $\sigma_{A}^{-}$ making it even since we have to end up in $\mathop{|0\rangle}\nolimits_{A}$) out of the total $n=2m$ branchings is $\binom{2m}{2\ell}$ and so is the number of nonzero strings in Eq. (57a). Identically, for (57b) we again set $n=2m$ and $p=\ell$ and obtain $\binom{2m}{2\ell-1}$ of nonzero strings ($2\ell-1$ left branchings and $m-\ell+1$ right branchings). For Eq. (57c) we set $n=2m+1$ and $p=\ell$. Hence $\binom{2m+1}{2\ell}$ is the total number of nonzero strings. Finally, for Eq. (57d) we have $n=2m+1$ and we now set $q=\ell$. In this way we obtain the same number of nonzero strings $\binom{2m+1}{2\ell}$ as in the previous case. We derived the binomial coefficients in Eqs. (7) and (7). The sums collect all Pauli operator products leading to the same target state $\mathop{|ij\rangle}\nolimits$. The proposition statement is then obtained by inspecting the LHS of Eqs. (7) and (7) where the boson operators $A_{\pm},B_{\pm}$ accompany the corresponding ladder operators and the RHS consist of the products of the boson operators from nonzero strings. ∎

Remark*.*

Fig. 5 illustrates the proof. It also illustrates the fact that the same analysis can be done for any initial state $\mathop{|ij\rangle}\nolimits$. Depending on the values of $i$ and $j$ the arrows in the binary tree will represent different ladder operators leading to the same counting and different boson operators on the RHS of Eqs. (7) and (7). In this way we are able to find all the unitary matrix elements and therefore we can act on any input state (pure or mixed) from a space spanned by $\mathop{|ij\rangle}\nolimits_{AB}$.

Example 4*.*

Setting $m=0,1,2$ in Eq. (7) and $m=0,1$ in Eq. (7) (corresponding to the perturbative order $n=4$) we quickly reproduce the fourth order calculation in [44] and the result can be applied to a variety of situations [35, 36, 42, 38]. With the same ease we can obtain an output state for $m\gg 0$.

A non-trivial check is the trace of (52) being equal to one irrespective of the perturbative order. Case $n=4$ was verified in [44] and a straightforward calculation shows that it is true for $n=6$ as well222Note that the fourth order is an absolute must if one wants to calculate any entropic quantity to properly assess the importance of entanglement for quantum communication. If only the second order is calculated (almost always the case with a notable exception of [44]) then the $|11\rangle\!\langle 11|$ component of (52) is zero and two of the eigenvalues of $\omega_{AB}$ are negative.. ∎

We use Proposition 7 to insert Eqs. (7) and (7) to Eq. (52) and get the output state components to the $n$-th order. Let’s take a look at the $|00\rangle\!\langle 00|$ component:

[TABLE]

We identified $A_{+}\rightleftharpoons 1,A_{-}\rightleftharpoons 2,B_{+}\rightleftharpoons 3$ and $B_{-}\rightleftharpoons 4$ and in this form it is ready for Theorem 2 and onto Theorem 4. The number of terms in Eq. (4) is polynomial in $m$ and therefore in $n$ as well.

Hafnians and the number of perfect matchings of a graph

Following Def. 4 we introduce some additional concepts from graph theory.

Definition 7*.*

A simple undirected graph $G=(V,E)$ on $m=|V|$ vertices is called complete if all vertices are connected with each other. That is, $|E|=m(m-1)/2$ is the total number of edges. A graph which is not complete is called incomplete.

Definition 8*.*

A perfect matching (or 1-factor) of a simple undirected graph $G=(V,E)$ is a graph $G_{PM}=(V^{\prime},E^{\prime})\subset G$ where every vertex is paired with exactly one other vertex. Hence $G$ contains a perfect matching if $|V|=m$ is even in which case $V^{\prime}=V$ and $|E^{\prime}|=|V|/2$.

The number of perfect matchings of a graph varies, depending on the graph type. The majority of graphs are incomplete and the calculation of the number of perfect matchings not known to be tractable. One of the very few cases where the analytical form is known is $K_{m}$ – a complete graph on $m=2n$ vertices where $\#(PM_{K_{2n}})=(2n)!/(2^{n}n!)=(2n-1)!!$.

Example 5*.*

The following picture shows all the perfect matching for $K_{4}$.

∎

The picture resembles the way a four-point correlation function splits into as a result of Wick’s (Isserlis’) theorem. Caianiello formalized this connection [45] by introducing the hafnian of a $2n\times 2n$ matrix $A$ (here we use the definition from [61])

[TABLE]

where the sum goes over all $(2n-1)!!$ unordered disjoint pairs of the set $\{1,\dots,2n\}$. For the above example of a complete graph on four vertices we calculate the hafnian of its adjacency matrix

[TABLE]

and find $\operatorname*{Haf}{[A(K_{4})]}=1+1+1=3$ agreeing with Example 5. Indeed, the hafnian of any complete graph is precisely the number of summands of the products of two-point correlation functions in Wick’s (Isserlis’) formula, Eq. 8, as envisaged by Caianiello.

Our contribution in this section is the following. While the number of perfect matchings for a complete graph is known we may face the problem of systematically listing all of them. Following the proof of Theorem 2, namely the set of Diophantine equations in (12), we set $f=2n$ and $\ell_{i}=1,\forall i$ and by solving the system obtain the desired perfect matchings.

5. Generalization beyond boson scalar theories

Wick’s theorem is at the core of any perturbative calculation and so the method based on finding all nonnegative solutions of a Diophantine system can in principle be used everywhere. But there are several, mainly technical, obstacles in extending the current formalism to more general theories. One of the was mentioned in Section 2 regarding the presence of derivative interactions in the Lagrangian.

Another obstacle to overcome in order to generalize the developed scheme is to consider Wick’s theorem for anticommuting fields. This is a standard tool in QFT [1] where one has to keep track of the parity and change the sign whenever two fermions swap their place in a correlation function. There seems to be no equivalent of Isserlis’ theorem for anticommuting variables but one can envisage a generalization of random Gaussian variables to a Gaussian integral over (anticommuting) random Grassmann variables. After all, Grassmann variables are routinely used to study the properties of fermions in the path integral formulation of QFT using similar algebraic properties of the CAR and Grassmann algebra. Then, a generalization of our refinement of Wick’s/Isserlis’ theorem (Theorem 2) is straightforward as well but due to the nilpotency of Grassmann variables there is really not much to generalize as the higher powers are missing.

Once we know how treat anticommuting fields we proceed in the same way by solving a system of Diophantine equations. The difference is (and this is true for any complex boson theory as well) that we have to disregard the solutions where a contraction of two fields is zero. For example, for the second order expansion of the Yukawa interaction term $\pazocal{L}_{\mathrm{int}}=-g\psi^{\dagger}(x)\psi(x)\phi(x)$, where $\psi$ is anticommuting, the contribution given by contracting $\psi(x)\psi(y)$ is zero. The Diophantine system, of course, treats all fields equally and so we have to exclude such contributions manually. But this can be easily done just by going through the list of Diophantine solutions and removing all those containing at least one zero propagator such as $\mathop{\left<0_{M}\,\right|}\nolimits\mathsf{T}\{\psi^{\dagger}(x)\psi^{\dagger}(y)\}\mathop{|0_{M}\rangle}\nolimits,\mathop{\left<0_{M}\,\right|}\nolimits\mathsf{T}\{\psi(x)\phi(x)\}\mathop{|0_{M}\rangle}\nolimits$ or any other kinematically forbidden process. By easily we mean without an additional computational complexity cost.

Finally, one has to face the calculation of the graph automorphism order to get the right multiplicities. Already starting with complex boson theories the graphs are in general directed, distinguishing between scattering of particles and antiparticles. This will result in a rather minor modification of Lemma 6. In particular, directed edges or edges of different types (different fields) will decrease the grafted graph symmetry compared to the real scalar case.

6. Discussion and open questions

Among several open problems stands out the following one. We have seen that in order to get the right multiplicity factor of a perturbative contribution one has to use a potentially computationally expensive procedure of calculating the graph automorphism order or graph isomorphism membership problem. This procedure bundles the graphs from the same class of physically indistinguishable scattering processes. Recall that they are represented by different solutions of a Diophantine system and every such a solution is an interior point of a certain convex polyhedron. Could it be that solutions from the same class form an exceptional subset in the polyhedron (perhaps even with nice geometric properties) or are they scattered randomly? In the former case we could use such knowledge to avoid the potentially costly calculation of the graph isomorphisms (costly for large graphs in extremely high perturbative orders).

A closely related question is whether the study of the contributions from one perturbative order (for a given Lagrangian and a fixed number of interacting fields) tells us something about the structure of the contributions from higher perturbative orders. All nonnegative solutions come from a Diophantine system whose size (the number of interior lattice points) grows but its structure is very similar across the perturbative orders. Does geometry of the corresponding convex polyhedra have anything in common across perturbative orders? If the answer is affirmative, could a new insight be gained into how to sum the related perturbative contributions?

Another outstanding problem is a generalization of the major application of our method: an efficient perturbative calculation of the interaction Hamiltonian for a pair of Unruh-DeWitt (UDW) detectors in Minkowski spacetime. The result is independent on the detector details and works for any smearing, envelope function or trajectory by efficiently decomposing a generic $2n$-point correlation function to a product of two-point Green’s functions where the detector details play a role. An efficient calculation means that the number of perturbative terms increases polynomially with the perturbative order of the coupling constant. It would be interesting to generalize the proofs for any number of Unruh-DeWitt detectors to study the multipartite structure of Minkowski and other spacetimes.

A related open problem is to clarify whether the UDW perturbative expansion suffers from the convergence issues. Everything suggests that the series is not of an asymptotic character. It is unlikely but not entirely impossible that the density matrix components or even the unitary entries themselves are actually analytically summable. To answer these questions one could perhaps first look at single UDW detector coupled to Minkowski vacuum as the simplest case. Another interesting question is what happens with the divergencies in the momentum space after all perturbative components are added. Is renormalization still necessary as claimed in [62] or do the divergencies cancel out?

Acknowledgements

This material is based upon work supported by the Air Force Office of Scientific Research under award number FA9550-17-1-0083.