ADG: Automated generation and evaluation of many-body diagrams II. Particle-number projected Bogoliubov many-body perturbation theory
aa r X i v : . [ nu c l - t h ] J u l ADG : Automated generation and evaluation of many-body diagramsII. Particle-number projected Bogoliubov many-body perturbation theory
P. Arthuis a,b,c , A. Tichai d,e,f,g , J. Ripoche c , T. Duguet b,h a Department of Physics, University of Surrey, Guildford GU2 7XH, United Kingdom b IRFU, CEA, Université Paris-Saclay, F-91191 Gif-sur-Yvette, France c CEA, DAM, DIF, F-91297 Arpajon, France d Max-Planck-Institut für Kernphysik, Heidelberg, Germany e Institut für Kernphysik, Technische Universität Darmstadt, Darmstadt, Germany f ExtreMe Matter Institute EMMI, GSI Helmholtzzentrum für Schwerionenforschung GmbH, Darmstadt, Germany g ESNT, IRFU, CEA, Université Paris-Saclay, F-91191 Gif-sur-Yvette, France h KU Leuven, Instituut voor Kern- en Stralingsfysica, 3001 Leuven, Belgium
Abstract
We describe the second version (v2.0.0) of the code
ADG that automatically (1) generates all valid off-diagonal
Bogoliubovmany-body perturbation theory diagrams at play in particle-number projected Bogoliubov many-body perturbationtheory (PNP-BMBPT) and (2) evaluates their algebraic expression to be implemented for numerical applications. Thisis achieved at any perturbative order p for a Hamiltonian containing both two-body (four-legs) and three-body (six-legs)interactions (vertices). All valid off-diagonal BMBPT diagrams of order p are systematically generated from the setof diagonal , i.e., unprojected, BMBPT diagrams. The production of the latter were described at length in Ref. [1]dealing with the first version of ADG . The automated evaluation of off-diagonal BMBPT diagrams relies both on theapplication of algebraic Feynman’s rules and on the identification of a powerful diagrammatic rule providing the resultof the remaining p -tuple time integral. The new diagrammatic rule generalizes the one already identified in Ref. [1] toevaluate diagonal BMBPT diagrams independently of their perturbative order and topology. The code ADG is written in
Python3 and uses the graph manipulation package
NetworkX . The code is kept flexible enough to be further expandedthroughout the years to tackle the diagrammatics at play in various many-body formalisms that already exist or are yetto be formulated.
Keywords: many-body theory, ab initio , perturbation theory, diagrammatic expansion
PACS:
NEW VERSION PROGRAM SUMMARY
Program Title:
ADG
Licensing provisions:
GPLv3
Programming language:
Python3
Supplementary material:
Yes.
Journal reference of previous version:
P. Arthuis, T. Duguet,A. Tichai, R.-D. Lasseri and J.-P. Ebran, "ADG: Automatedgeneration and evaluation of many-body diagrams I. Bogoli-ubov many-body perturbation theory", Computer PhysicsCommunications 240 (2019), pp. 202-227.
Does the new version supersede the previous version?:
Yes.
Reasons for the new version:
Incorporation of a new formalisminto the program.
Summary of revisions:
Addition of off-diagonal BMBPT tothe formalisms for which diagrams can be generated, fix of awrong symmetry factor, move of the codebase from Python2to Python3 while maintaining support for Python2, variousoptimizations to reduce the time and memory necessary to the
Email addresses: [email protected] (P. Arthuis), [email protected] (A. Tichai), [email protected] (J. Ripoche), [email protected] (T. Duguet) program.
Nature of problem:
As formal and numerical developments inmany-body-perturbation-theory-based ab initio methods makehigher orders reachable, manually producing and evaluatingall the diagrams becomes rapidly untractable as both theirnumber and complexity grow quickly, making it prone tomistakes and oversights.
Solution method:
Diagonal BMBPT diagrams are encodedas square matrices known as oriented adjacency matrices ingraph theory, and then turned into graph objects using the
NetworkX package. Off-diagonal BMBPT diagrams can thenbe generated from those graphs. Checks on the diagrams andevaluation of their time-integrated expression are eventuallydone on a purely diagrammatic basis. HF-MBPT diagramsare produced and evaluated as well using the same principle.
1. Introduction
In Ref. [1], the first version of the code
ADG was de-scribed. The code was designed to automatically (1) gen-erate all valid Bogoliubov many-body perturbation the-
Preprint submitted to Computer Physics Communications July 6, 2020 ry (BMBPT) diagrams and (2) evaluate their algebraicexpression to be implemented for numerical applications.This was achieved at any perturbative order p for a Hamil-tonian containing both two-body (four-legs) and three-body (six-legs) interactions (vertices). This code devel-opment took place in the context of the rapidly evolvingfield of quantum many-body calculations for fermionic sys-tems, i.e., atomic nuclei, molecules, atoms or solids. In nu-clear physics in particular, the past decade has witnessedthe development and/or the application of formalisms [2–13] among which is Bogoliubov many-body perturbationtheory [12, 14]. This profusion of diagrammatic methods,along with the rapid progress of computational power al-lowing for high-order implementations, welcomes the de-velopment of a versatile code capable of both generatingand evaluating many-body diagrams. This was the goalachieved in Ref. [1] for BMBPT.The aim of BMBPT is to tackle (near) degenerateFermi systems, e.g., open-shell nuclei displaying a super-fluid character, by perturbatively expanding the exact so-lution of the Schrödinger equation around a so-called Bo-goliubov reference state. This tremendous benefit is ob-tained at the price of using a reference state that breaks U (1) global-gauge symmetry associated with the conser-vation of particle number in the system. Given that thebreaking of a symmetry cannot actually be realized in afinite quantum system, BMBPT calculations come witha (hopefully mild) contamination of physical observables.In this context, BMBPT must actually be seen as a firststep towards a more general diagrammatic method, theparticle-number projected Bogoliubov many-body pertur-bation theory (PNP-BMBPT) [9], in which the brokensymmetry is exactly restored at any truncation order.It is thus the goal of the present paper to extend theformal and code developments performed in Ref. [1] toPNP-BMBPT. For a reason that will become clear lateron, the diagrammatics at play in PNP-BMBPT is coinedas the off-diagonal BMBPT diagrammatics from which the diagonal
BMBPT diagrammatics encountered in straightBMBPT is recovered in a particular limit, i.e., diagonalBMBPT diagrams characterize the subset of off-diagonalBMBPT diagrams that remains non-zero in that limit. Inthis context, the new version of the code
ADG is capable ofautomatically (1) generating all valid off-diagonal
BMBPTdiagrams and of (2) evaluating their algebraic expressionto be implemented for numerical applications. This isachieved at any perturbative order p for a Hamiltoniancontaining both two-body (four-legs) and three-body (six-legs) interactions (vertices). The numerical tool remainsflexible enough to be expanded throughout the years totackle the diagrammatics at play in yet other many-bodyformalisms (already existing or yet to be formulated).The paper is organized as follows. Section 2 introducesthe basic ingredients that are needed to elaborate on thePNP-BMBPT formalism in Sec. 3. The master equations,the associated diagrammatics and the difficulties to over-come in order to achieve an automated generation and eval- uation of diagrams of arbitrary orders are detailed. Build-ing on this, Secs. 5 and 6 detail the method developped toreach such an objective. While Sec. 9 details how the ADG code operates, conclusions are given in Sec. 10. Three ap-pendices follow to provide details regarding the formalismand the structure of the program.
2. Basic ingredients
Expansion many-body methods consist of expandingexact A-body quantities, e.g., the ground-state energy,around a reference state. By default, the reference state isnaturally chosen to carry the symmetry quantum numbersof the exact target state dictacted by the Hamiltonian. Asfor U (1) global gauge symmetry associated with particle-number conservation, this typically leads to using a Slaterdeterminant as a reference state. In open-shell nuclei, how-ever, the reference Slater determinant is degenerate withrespect to elementary, i.e., particle-hole, excitations, whicheventually makes the many-body expansion singular. Thissingularity relates to Cooper pair’s instability associatedwith the superfluid character of singly (doubly) open-shellnuclei.One option to address Cooper pair’s instability is to al-low the reference state to break U (1) symmetry, i.e., tobe chosen under the form of a more general Bogoliubovproduct state. Doing so, the degeneracy of the Slater de-terminant with respect to particle-hole excitations is liftedand commuted into a degeneracy with respect to symme-try transformations of the group. As a consequence, the ill-defined (i.e., singular) expansion of exact quantities is re-placed by a well-behaved one. This task is indeed achievedby straight BMBPT [1, 12] in a perturbative setting.Eventually, however, the degeneracy with respect to U (1) transformations must also be lifted via the restora-tion of the symmetry. Indeed, the breaking of a symme-try is fictitious in finite quantum systems and can at bestbe characterized as being emergent [15, 16]. In this re-gards, straight BMBPT only restores the symmetry in thelimit of an all-orders resummation, and, thus, includes asymmetry contamination at any finite order. Thus, anexplicit extension of the formalism is necessary in prac-tice to restore good particle-number at any truncation or-der. Within a perturbative setting, this task is achieved byPNP-BMBPT [9]. Before detailing this many-body formal-ism and its associated diagrammatics, the present sectionintroduces the necessary ingredients. U (1) group Let us consider the abelian compact Lie group U (1) ≡{ S ( ϕ ) , ϕ ∈ [0 , π ] } associated with the global rotation ofan A-body fermion system in gauge space. As U (1) isconsidered to be a symmetry group of the Hamiltonian H ,commutation relations[ H, S ( ϕ )] = [ A, S ( ϕ )] = 0 , (1)2old for any ϕ ∈ [0 , π ], where A denotes the particle-number operator. We utilize the unitary representation of U (1) on Fock space F given by S ( ϕ ) = e iAϕ . (2)One is presently interested in evaluating a ground-stateobservable O A0 whose associated operator O commuteswith H and A . Generically speaking, this means that theeigenequations of interest in PNP-BMBPT are given by O | Ψ A0 i = O A0 | Ψ A0 i , (3)where O ≡ H or A and where Eq. (3) explicitly stipulatesthat energy eigenstates | Ψ A µ i carry good symmetry quan-tum number A.Given these many-body eigenstates, matrix elements ofthe irreducible representations (IRREPs) of U (1) are givenby h Ψ A µ | S ( ϕ ) | Ψ A ′ µ ′ i ≡ e i A ϕ δ AA ′ δ µµ ′ . (4)The volume of the group is given by v U (1) ≡ Z π dϕ = 2 π , while the orthogonality of IRREPs reads as12 π Z π dϕ e − i A ϕ e + i A ′ ϕ = δ AA ′ . (5)Whenever a many-body state is not an eigenstate ofthe particle-number operator, it is easy to demonstratethat the eigen-component of A associated with particle-number A can be extracted from it via the application ofthe projection operator P A defined through P A ≡ π Z π dϕ e − i A ϕ S ( ϕ ) . (6) The set up of the many-body formalism starts with theintroduction of the Bogoliubov reference state | Φ i ≡ C Y k β k | i , (7)where C is a complex normalization constant and | i de-notes the physical vacuum. The Bogoliubov state is a vac-uum for the set of quasi-particle operators obtained fromthose associated with a basis of the one-body Hilbert spacevia a unitary linear transformation of the form [17] β k ≡ X p U ∗ pk c p + V ∗ pk c † p , (8a) While the focus is presently on ground-state quantities, exten-sions of the formalism to excited states or transitions of non-scalaroperators can be envisioned. β † k ≡ X p U pk c † p + V pk c p , (8b)i.e., β k | Φ i = 0 for all k . One possiblity to specify theBogoliubov reference state | Φ i is to require that it solvesthe Hartree-Fock-Bogoliubov (HFB) variational problem.This fixes the transformation matrices ( U, V ) [17] and de-livers the set of quasi-particle energies { E k > } defin-ing the unperturbed part of the Hamiltonian later on (seeEqs. (27)-(28)). We do not impose this choice here suchthat the reference state and the associated unperturbedHamiltonian can be defined more generally.The Bogoliubov reference state is not an eigenstate ofthe particle-number operator A . Although the presentobjective is to perturbatively correct | Φ i while exactlyrestoring the particle number at any truncation order,the fact that the reference state breaks the symmetry inthe first place requires to work with the grand potentialΩ ≡ H − λA rather than with H itself in the set up of themany-body formalism [9], where λ denotes the chemicalpotential. The operator O of interest typically contains one-body,two-body and three-body terms . The operator in theSchrödinger representation is expressed in an arbitrary ba-sis of the one-body Hilbert space as O ≡ o [2] + o [4] + o [6] (9) ≡ o + o + o ≡ X p p o p p c † p c p + 1(2!) X p p p p o p p p p c † p c † p c p c p + 1(3!) X p p p p p p o p p p p p p c † p c † p c † p c p c p c p . Each term o kk of the particle-number conserving opera-tor O is characterized by the equal number k of particlecreation and annihilation operators. The class o [2 k ] is noth-ing but the term o kk of k -body character. The maximumvalue deg_max ≡ Max 2 k defines the rank of the operator O . Matrix elements are fully antisymmetric, i.e., o kkp ...p k p k +1 ...p k = ( − σ ( P ) o kkP ( p ...p k | p k +1 ...p k ) , (10)where σ ( P ) refers to the signature of the permutation P .The notation P ( . . . | . . . ) denotes a separation into the k particle-creation operators and the k particle-annihilation Higher-body operators can be employed as well. From the for-mal point of view, it poses no fundamental difficulty but furthercomplicates the diagrammatic formalism. As for the automated gen-eration of diagrams, it poses no fundamental difficulty but requires tohandle the memory needed to deal with the increased combinatorialcomplexity. O with respectto the Bogoliubov vacuum | Φ i , thus, obtaining O ≡ O [0] + O [2] + O [4] + O [6] (11) ≡ O + h O + { O + O } i + h O + { O + O } + { O + O } i + h O + { O + O } + { O + O } + { O + O } i = O + 1(1!) X k k O k k β † k β k + 12! X k k n O k k β † k β † k + O k k β k β k o + 1(2!) X k k k k O k k k k β † k β † k β k β k + 13!1! X k k k k n O k k k k β † k β † k β † k β k + O k k k k β † k β k β k β k o + 14! X k k k k n O k k k k β † k β † k β † k β † k + O k k k k β k β k β k β k o + 1(3!) X k k k k k k O k k k k k k β † k β † k β † k β k β k β k + 12! 4! X k k k k k k n O k k k k k k β † k β † k β † k β † k β k β k + O k k k k k k β † k β † k β k β k β k β k o + 15!1! X k k k k k k n O k k k k k k β † k β † k β † k β † k β † k β k + O k k k k k k β † k β k β k β k β k β k o + 16! X k k k k k k n O k k k k k k β † k β † k β † k β † k β † k β † k + O k k k k k k β k β k β k β k β k β k o , (12)where the expressions of the matrix elements of each oper-ator O ij in terms of those of the operators o kk and of the( U, V ) matrices can be found in Ref. [18] for terms up to O [4] and in Ref. [19] for O [6] . Each term O ij is character-ized by its number i ( j ) of quasi-particle creation (annihila-tion) operators. Because O has been normal-ordered with respect to | Φ i , all quasi-particle creation operators (if any)are located to the left of all quasi-particle annihilation op-erators (if any). The class O [2 k ] groups all the terms O ij of effective k -body character, i.e., with i + j = 2 k . Theoperator being overall unchanged by the normal-orderingprocedure, its rank deg_max ≡ Max 2 k remains itself un-changed. Matrix elements are fully antisymmetric, i.e., O ijk ...k i k i +1 ...k i + j = ( − σ ( P ) O ijP ( k ...k i | k i +1 ...k i + j ) . (13)More details and properties can be found in Refs. [9, 18,19].State-of-the-art many-body calculations are typicallyperformed within the normal-ordered two-body (NO2B)approximation [20–23]. However, the naive adaptation ofthe NO2B approximation to many-body formalisms basedon a particle-number breaking reference state, which wouldresults in neglecting the residual three-body part O [6] , hasbeen shown to be fundamentally inappropriate. As aresult, a particle-number conserving normal-ordered two-body (PNO2B) approximation was designed [19]. The neteffect of the PNO2B approximation is to modify in a spe-cific way the matrix elements at play in O [4] in additionto fully neglecting O [6] . In the present work, however, thediagrammatic is anyway worked out in presence of the ef-fective three-body part, i.e., in presence of six-legs vertices(see below), which significantly increases the number ofpossible diagrams at a given order and the complexity oftheir topology. Correspondingly, the code can eventuallybe run with or without including the effective three-bodypart of the operators.
3. Many-body formalism
Taking the Hermitian conjugate of Eq. (3) and right-multiplying by an arbitrary auxiliary state | Θ i (such that h Ψ A0 | Θ i , A0 = h Ψ A0 | O | Θ ih Ψ A0 | Θ i . (14)Choosing | Θ i ≡ | Φ i and expanding h Ψ A0 | around it leadsto straight, i.e., symmetry-breaking, BMBPT [1]. In thepresent work, the auxiliary state is taken as | Θ i ≡ P A | Φ i such that the symmetry is explicitly restored by the pres-ence of the projection operator P A even after expand-ing and truncating h Ψ A0 | around the Bogoliubov referencestate.In this context, Eq. (14) becomesO A0 = h Ψ A0 | OP A | Φ ih Ψ A0 | P A | Φ i , (15)such that inserting Eq. (6) as well as introducing the so-called off-diagonal norm and operator kernels N ( ϕ ) ≡ h Ψ A0 | Φ( ϕ ) i , (16a)4 ( ϕ ) ≡ h Ψ A0 | O | Φ( ϕ ) i , (16b)where | Φ( ϕ ) i ≡ S ( ϕ ) | Φ i denotes the gauge-rotated refer-ence state, leads to the working formO A0 = R π dϕ e − i A ϕ O ( ϕ ) R π dϕ e − i A ϕ N ( ϕ ) . (17)Equation (17) constitutes the master equation on whichthe PNP-BMBPT formalism is built. In absence of theprojection operator, one recovers BMBPT’s master equa-tion under the form O A0 = O (0) , (18)where intermediate normalization N (0) = h Ψ A0 | Φ i = 1with the unrotated Bogoliubov reference state has beenused. Equations (17) and (18) are obviously equivalentin the exact limit but differ as soon as h Ψ A0 | is expandedaround h Φ | and truncated. Introducing the evolution operator in imaginary time U ( τ ) ≡ e − τ Ω , (19)with τ real, allows one to write the ground state as | Ψ A0 i = lim τ →∞ | Ψ( τ ) i ≡ lim τ →∞ U ( τ ) | Φ ih Φ |U ( τ ) | Φ i , (20)where h Φ | Ψ( τ ) i = 1 for all τ . With this definition at hand,the off-diagonal kernels entering Eq. (17) read as N ( ϕ ) ≡ N ( ϕ ) N (0) = lim τ →∞ h Φ |U ( τ ) | Φ( ϕ ) ih Φ |U ( τ ) | Φ i , (21a) O ( ϕ ) ≡ O ( ϕ ) N (0) = lim τ →∞ h Φ |U ( τ ) O | Φ( ϕ ) ih Φ |U ( τ ) | Φ i . (21b)The off-diagonal kernels N ( ϕ ) and O ( ϕ ) are the many-body quantities to be approximated via a viable expansionmethod from which N (0) and O (0) can be obtained as aparticular case [1, 12]. The time is given in units of MeV − . The result is obtained by inserting a complete set of energy eigen-states in both the numerator and the denominator. The chemical potential λ is fixed such that Ω A for the targetedparticle number A is the lowest value of all Ω A µ over Fock space,i.e., it penalizes systems with larger number of particles such thatΩ A < Ω A µ for all A > A while maintaining at the same time thatΩ A < Ω A µ for all A < A . This is practically achievable only if E A0 is strictly convex in the neighborhood of A , which is generally butnot always true for atomic nuclei. The off-diagonal norm kernel plays a particular rolegiven that it does not actually involve a non-trivial op-erator, which makes its perturbative expansion differentfrom the expansion of an operator kernel [9]. In fact, itcan be trivially related to the particle-number operatorkernel through ddϕ N ( ϕ ) = i A ( ϕ ) . (22)Accessing N ( ϕ ) via the integration of Eq. (22) ensures thatEq. (3) applied to O ≡ A and rewritten as Eq. (17) deliversthe expected result A A0 = A even when A ( ϕ ) is computedapproximately through, e.g., perturbation theory. Indeed,as long as Eq. (22) is enforced, one has R π dϕ e − i A ϕ A ( ϕ ) R π dϕ e − i A ϕ N ( ϕ ) = − i R π dϕ e − i A ϕ ddϕ N ( ϕ ) R π dϕ e − i A ϕ N ( ϕ )= + i R π dϕ [ ddϕ e − i A ϕ ] N ( ϕ ) R π dϕ e − i A ϕ N ( ϕ )= A , which is the required result for the symmetry of presentinterest to be exactly restored. Further introducing thefactorization of an arbitrary operator kernel O ( ϕ ) ≡ o ( ϕ ) N ( ϕ ) (23)where o ( ϕ ) denotes the connected/linked part of the op-erator kernel [9], one arrives at the first-order ordinarydifferential equation (ODE) fulfilled by the norm kernel ddϕ N ( ϕ ) − i a ( ϕ ) N ( ϕ ) = 0 , (24)whose closed-form solution reads as N ( ϕ ) = e i R ϕ dφ a ( φ ) . (25)From the computation of a ( ϕ ), the off-diagonal normkernel is consistently obtained. Eventually, the con-nected/linked part o ( ϕ ) of an operator kernel O ( ϕ ) is thesole quantity one needs to effectively focus on in order toimplement the complete PNP-BMBPT formalism. As amatter of fact, it is not by chance given that, as will bediscussed below, o ( ϕ ) is size-extensive and properly scaleswith system size, which translates into the fact that it ef-fectively displays a connected expansion.
4. Perturbation theory
The grand potential is split into an unperturbed partΩ and a residual part Ω Ω = Ω + Ω . (26)5or a given number of interacting fermions, the key is tochoose Ω with a low-enough symmetry for its ground state | Φ i to be non-degenerate with respect to elementary ex-citations. For open-shell superfluid nuclei, this leads tochoosing an operator Ω that breaks particle-number con-servation, i.e., while Ω commutes with U (1) transforma-tions, we are interested in the case where Ω , and thus Ω ,do not. More specifically, the partioning is defined throughΩ ≡ Ω + ¯Ω , (27a)Ω ≡ Ω + ˘Ω + Ω + Ω + Ω + Ω + Ω + Ω + Ω + Ω + Ω + Ω + Ω + Ω + Ω , (27b)with ˘Ω ≡ Ω − ¯Ω and where the normal-ordered one-body part of Ω is diagonal, i.e., ¯Ω ≡ X k E k β † k β k , (28)with E k > for all k .Introducing many-body states generated via an evennumber of quasi-particle excitations of the Bogoliubov vac-uum | Φ k k ... i ≡ β † k β † k . . . | Φ i , (29)the unperturbed grand potential Ω is fully characterizedby its complete set of orthonormal eigenstates in Fockspace Ω | Φ i = Ω | Φ i , (30a) Ω | Φ k k ... i = (cid:2) Ω + ǫ k k ... (cid:3) | Φ k k ... i , (30b)where the strict positivity of unperturbed excitations ǫ k k ... ≡ E k + E k + . . . characterizes the lifting of theparticle-hole degeneracy authorized by the spontaneousbreaking of U (1) symmetry in open-shell nuclei at themean-field level.In the particular case where | Φ i solves the HFB varia-tional problem, one has that Ω = ˘Ω = Ω = 0 suchthat Ω reduces to Ω [4] + Ω [6] . This choice defines the canonical version of (PNP-)BMBPT and reduces signifi-cantly the number of non-zero diagrams to be considered.However, we do not make this a priori hypothesis suchthat the reference state | Φ i and the corresponding unper-turbed grand potential Ω can be defined more generally,eventually leading to the appearance of non-canonical di-agrams involving Ω , ˘Ω and Ω vertices.On the basis of the above splitting of Ω , one introducesthe interaction representation of operators in the quasi-particle basis, e.g., O ( τ ) ≡ e + τ Ω O e − τ Ω (31) = 13! X k k k k O k k k k β † k ( τ ) β † k ( τ ) β † k ( τ ) β k ( τ ) , where β k ( τ ) ≡ e + τ Ω β k e − τ Ω = e − τE k β k , (32a) β † k ( τ ) ≡ e + τ Ω β † k e − τ Ω = e + τE k β † k . (32b) Expanding the evolution operator in powers of Ω [24] U ( τ ) ≡ e − τ Ω = e − τ Ω T e − R τ d τ Ω ( τ ) , (33)where T denotes the time-ordering operator , one ob-tains [9] the expansion of interest o ( ϕ ) ≡ lim τ →∞ h Φ |U ( τ ) O | Φ( ϕ ) ih Φ |U ( τ ) | Φ( ϕ ) i = lim τ →∞ h Φ | T e − R τ dt Ω ( t ) O | Φ( ϕ ) i c = h Φ | O | Φ( ϕ ) i− Z + ∞ dτ h Φ | T [Ω ( τ ) O (0)] | Φ( ϕ ) i c + 12! Z + ∞ dτ dτ h Φ | T [Ω ( τ ) Ω ( τ ) O (0)] | Φ( ϕ ) i c − ... , (34)where the lower index c refers to the restriction to con-nected terms, i.e., contributions arising from the appli-cation of Wick’s theorem in which the associated stringof contractions necessarily involves all the operators atplay in the many-body matrix element under considera-tion. The time-independent operator O could be insertedat no cost within the time-ordering by providing it with afictitious and harmless time dependence t = 0 . Indeed, all Ω ( τ k ) operators appear to the left of O and occur at alarger time given that their corresponding time variablesare positive.Invoking perturbation theory consists of truncatingthe Taylor expansion of the time-evolution operator inEq. (34). Gathering all terms up to order p , o ( ϕ ) sumsmatrix elements of products of up to p + 1 time-dependentoperators . The running time variables are integrated overfrom to τ → + ∞ whereas the time label attributed to theoperator O itself remains fixed at t = 0 , i.e., contributionsof order p contain a p -tuple time integral that needs to beperformed to generate the end result under the requiredform.Given the off-diagonal character of the kernels, each ma-trix element in Eq. (34) is computed via the application The time-ordering operator orders a product of operators in de-creasing order according to their time labels (i.e., larger times to theleft) and multiplies the result with the signature of the permutationused to achieve the corresponding reordering. In agreement with Eq. (18), straight BMBPT is recovered fromEq. (34) for ϕ = 0 given that O A0 = O (0) = o (0) in this formalism. The expansion starts at order p = 0 that corresponds to the termcontaining no Ω operator and no time integral in Eq. (34). off-diagonal Wick’s theorem [25], which is applicableto matrix elements of operators between any two (non-orthogonal) left and right product states. As a result, dia-grams at play invoke a set of four off-diagonal unperturbedpropagators defined in the quasi-particle basis { β k ; β † k } as G + − (0) k k ( τ , τ ; ϕ ) ≡ h Φ | T [ β † k ( τ ) β k ( τ )] | Φ( ϕ ) ih Φ | Φ( ϕ ) i , (35a) G −− (0) k k ( τ , τ ; ϕ ) ≡ h Φ | T [ β k ( τ ) β k ( τ )] | Φ( ϕ ) ih Φ | Φ( ϕ ) i , (35b) G ++(0) k k ( τ , τ ; ϕ ) ≡ h Φ | T [ β † k ( τ ) β † k ( τ )] | Φ( ϕ ) ih Φ | Φ( ϕ ) i , (35c) G − +(0) k k ( τ , τ ; ϕ ) ≡ h Φ | T [ β k ( τ ) β † k ( τ )] | Φ( ϕ ) ih Φ | Φ( ϕ ) i . (35d)By virtue of the off-diagonal elementary contractionsworked out in Appendix A.4, the four off-diagonal prop-agators are equal to G + − (0) k k ( τ , τ ; ϕ ) = − e − ( τ − τ ) E k θ ( τ − τ ) δ k k , (36a) G −− (0) k k ( τ , τ ; ϕ ) = + e − τ E k e − τ E k R −− k k ( ϕ ) , (36b) G ++(0) k k ( τ , τ ; ϕ ) = 0 , (36c) G − +(0) k k ( τ , τ ; ϕ ) = + e − ( τ − τ ) E k θ ( τ − τ ) δ k k , (36d)where both normal propagators are actually related via an-tisymmetry under the exchange of time and quasi-particlelabels. The higher generality and complexity of the off-diagonal BMBPT diagrammatics of present interest com-pared to the straight BMBPT diagrammatics discussed inRef. [1] is due to the presence of the anomalous propagator G −− (0) ( ϕ ) that carries the full gauge-angle dependence. Inparticular, the possibility to form anomalous propagatorssignificantly increases the combinatorics [9]. Eventually,the two diagrammatics coincide in the limit ϕ = 0 giventhat G −− (0) (0) = 0 . All in all, the present extension of the ADG code amounts to dealing with this higher generalityand complexity, which itself originates from the presenceof different left and right vacua in the off-diagonal kernel o ( ϕ ) (see Eq. (34)).Equal-time propagators can solely arise from contract-ing two quasi-particle operators belonging to the samenormal-ordered operator displaying creation operators tothe left of annihilation ones. As a result, one finds that [9] G + − (0) k k ( τ, τ ; ϕ ) ≡ , (37a) G −− (0) k k ( τ, τ ; ϕ ) ≡ + e − τ ( E k + E k ) R −− k k ( ϕ ) , (37b) G ++(0) k k ( τ, τ ; ϕ ) ≡ , (37c) G − +(0) k k ( τ, τ ; ϕ ) ≡ , (37d)such that the sole non-zero equal-time contraction, andthus the sole contraction of an interaction vertex onto it- O [0] = O O [2] = O + O + O O [4] = O + O + O + O + O O [6] = O + O + O + O + O + O + O Figure 1: Canonical diagrammatic representation of normal-orderedcontributions to the operator O in the Schrödinger representation. self, is of anomalous character. Correspondingly, no con-traction of an interaction vertex onto itself can occur inthe diagonal case, i.e., for ϕ = 0 . The pedestrian application of the off-diagonal Wick’stheorem becomes quickly cumbersome as the order p in-creases. Furthermore, it leads to computing independentlymany contributions that are in fact identical. By iden-tifying the corresponding pattern, one can design a dia-grammatic representation of the various contributions andevaluate their algebraic expressions such that a single dia-gram captures all identical contributions at once. In orderto achieve this goal, one must first introduce the diagram-matic representation of the building blocks.The operator O expressed in the quasi-particle basis isdisplayed in the Schrödinger representation in Fig. 1 asa sum of Hugenholtz vertices denoting its various normal-ordered contributions O ij . The antisymmetrized matrixelement O ijk ...k i k i +1 ...k i + j must be assigned to the corre-sponding square vertex, where i ( j ) denotes the numberof lines traveling out of (into) the vertex and representing7 k k k + O k k k k = k k k k + O k k k k = k k k k − O k k k k Figure 2: Rules to apply when departing from the canonical dia-grammatic representation of a normal-ordered operator. Orientedlines can be rotated through the dashed line but not through the fullline. quasi-particle creation (annihilation) operators. The oper-ator O ( τ ) in the interaction representation possesses thesame diagrammatic except that a time τ is attributed toeach of the vertices, i.e., to each of the lines coming in orout of them.In the canonical representation used in Fig. 1, all ori-ented lines go up, i.e., lines representing quasi-particle cre-ation (annihilation) operators appear above (below) thevertex. Accordingly, indices k . . . k i must be assigned con-secutively from the leftmost to the rightmost line abovethe vertex, while k i +1 . . . k i + j must be similarly assignedconsecutively to lines below the vertex. In the diagram-matic representation of the observable O A0 , it is howeverpossible for a line to propagate downwards. This can beobtained unambiguously starting from the canonical repre-sentation of Fig. 1 at the price of adding a specific rule. Asillustrated in Fig. 2 for the diagram representing O , linesmust only be rotated through the right of the diagram, i.e.,going through the dashed line, while it is forbidden to ro-tate them through the full line. Additionally, a minus signmust be added to the amplitude O ijk ...k i k i +1 ...k i + j associ-ated with the canonical diagram each time two lines crossas illustrated in Fig. 2.Since the grand canonical potential Ω is involved in theevaluation of any observable O A0 , its own diagrammaticrepresentation is needed and displayed in Fig. 3. The onlydifference with Fig. 1 relates to the use of dots rather thansquare symbols to represent the vertices. The same is eas-ily done for other operators of interest, i.e., H and A . Itis to be noted that Ω has the same diagrammatic repre-sentation as Ω except that Ω must be omitted and Ω replaced by ˘Ω , which requires to use a different symbolfor that particular vertex .As the off-diagonal Wick theorem contracts pairs ofquasi-particle operators together, the lines entering the di-agrammatic representation of operators are eventually con-nected in the computation of the kernel o ( ϕ ) , thus, form-ing elementary contractions. Consequently, the four un-perturbed propagators at play also need to be representeddiagrammatically, which is done in Fig. 4. Here, the con-vention is that the left-to-right reading of a matrix element We omit to use a different symbol for ˘Ω in the following al-though it must be clear that the vertex with one line coming in andone line coming out does represent ˘Ω whenever it originates fromthe perturbative expansion of the evolution operator. This may beconfusing whenever O = Ω since in this case there can also be avertex Ω at fixed time t = 0. Ω [0] = Ω Ω [2] = Ω + Ω + Ω Ω [4] = Ω + Ω + Ω + Ω + Ω Ω [6] = Ω + Ω + Ω + Ω + Ω + Ω + Ω Figure 3: Canonical diagrammatic representation of normal-orderedcontributions to the grand potential operator Ω in the Schrödingerrepresentation. corresponds to the up-down reading of the diagram.
With the building blocks at hand, off-diagonal BMBPTFeynman diagrams representing the contributions to o ( ϕ ) are generated by assembling them according to a set oftopological rules [9]1. A Feynman diagram of order p consists of p vertices Ω i k j k ( τ k ) , i k + j k = 2 , or , along with one vertex O mn (0) , m + n = 0 , , or , that are connected byfermionic quasi-particle lines, i.e., via non-zero prop-agators G + − (0) , G − +(0) or G −− (0) .2. Each vertex is labeled by a time variable while eachline is labeled by two time labels associated with thetwo vertices the line is attached to.3. Generating all contributions to Eq. (34) requires toform all possible diagrams, i.e., contract quasi-particlelines attached to the vertices in all possible ways whilefulfilling the following restrictions.8 τ k τ k τ k τ k τ k τ k τ k τ G + − (0) k k ( τ , τ ; ϕ ) G −− (0) k k ( τ , τ ; ϕ ) G ++(0) k k ( τ , τ ; ϕ ) G − +(0) k k ( τ , τ ; ϕ ) Figure 4: Diagrammatic representation of the four unperturbed elementary one-body propagators G gg ′ (0) ( ϕ ). The convention is that theleft-to-right reading of a matrix element corresponds to the up-down reading of the diagram. (a) Restrict equal-time propagators starting andending at the same vertex to anomalous prop-agators. In the diagonal case, i.e., for ϕ = 0 , nosuch self-contraction may occur.(b) Restrict the set to connected diagrams, i.e., omitdiagrams containing parts that are not connectedto each other by either propagators or vertices.This implies in particular that the vertex O with no line can only appear at order p = 0 .(c) Because of the time-ordering relations carriedby the propagators (see Eq. (36)), normal lineslinking a set of vertices must not form an ori-ented loop. For two given vertices Ω i k j k ( τ k ) and Ω i k ′ j k ′ ( τ k ′ ) , it means that normal lines mustpropagate between them in the same direction.Correspondingly, normal lines connected to thegeneric operator O at fixed time must go outof it, i.e., upwards in time. Anomalous lines donot carry time-ordering relations and are, thus,not concerned by these restrictions. In the di-agonal limit, i.e., ϕ = 0 , where no anomalousline may be formed, the above constraint imposesthat contributing vertices O mn (0) can only havelines going out, i.e., one necessarily has n = 0 .(d) Restrict the set to vacuum-to-vacuum diagramsforming a set of closed loops with no external,i.e., unpaired, lines. This condition, togetherwith the fact that G ++(0) ( ϕ ) is identically zero,strongly constrains which normal-ordered parts Ω i k j k ( τ k ) and O mn (0) of the p + 1 involved oper-ators can be combined, i.e., the condition n a ≡ p X k =1 ( j k − i k ) + n − m ≥ , must be fulfilled. The number n a corresponds tothe number of anomalous propagators G −− (0) ( ϕ ) in the diagram. In the diagonal limit for which G −− (0) (0) = 0 , the set of combined operators arefurther reduced to n a = 0 .(e) Restrict the set to topologically distinct time-unlabelled diagrams, i.e., time-unlabelled dia-grams that cannot be obtained from one anothervia a mere displacement, i.e., translation, of thevertices. k k k k G −− (0) k k ( τ, τ ; ϕ ) O k k k k τ Figure 5: Convention to draw and read anomalous self-contractions.The example is given for a vertex Ω displaying a self-contraction. The way to translate off-diagonal BMBPT Feynman di-agrams into their mathematical expressions follows a setof algebraic rules1. Each of the p + 1 vertices contributes a factor, e.g., Ω ijk ...k i k i +1 ...k i + j with the sign convention detailed inSec. 4.3.2. Each of the n b ≡ p X k =1 ( j k + i k ) + n + m ! / , lines contributes a factor G gg ′ (0) k k ( τ k , τ k ′ ) , where g and g ′ characterize the type of elementary propagator theline corresponds to. According to Eq. (36), each ofthe n a anomalous propagators carries an exponentialfunction of the two time labels and an anomalous con-traction R −− k k ( ϕ ) while each of the n b − n a normalpropagators carries an exponential function and a stepfunction of the two time labels.3. A normal line can be interpreted as G − +(0) or G + − (0) depending on the ascendant or descendant reading ofthe diagram. Similarly, the ordering of quasi-particleand time labels of a propagator depends on the as-cendant or descendant reading of the diagram. All the lines involved in a given diagram must be inter-preted in the same way, i.e., sticking to an ascen-dant or descendant way of reading the diagram allthroughout. In the present work, the chosen conven-tion corresponds to reading diagrams from top to bot-tom, which further relates to reading the many-body9atrix element it originates from in Eq. (34) in a left-right fashion. It is the convention employed to repre-sent the four propagators in Fig. 4.4. The reading of an anomalous line linking two differ-ent vertices is unambiguous as long as one stick tothe up-down convention displayed in Fig. 4. How-ever, the up-down reading of a self-contraction is po-tentially ambiguous depending on the way the lineis actually drawn. As illustrated in Fig. 5, one mustfurther fix a convention based on the insertion of a fic-titious semi-straight, e.g., horizontal line originatingfrom the vertex that the self-contraction is forbiddento cross. Taking the semi-straight line as a referencepoint, the quasi-particle indices must be attributed tothe equal-time propagator in the order the lines arecrossed when going around the vertex in a clockwisefashion.5. All quasi-particle labels must be summed over whileall running time variables must be integrated overfrom to τ → + ∞ .6. A sign factor ( − p + n c , where p denotes the order ofthe diagram and n c denotes the number of crossinglines in the diagram, must be considered . The over-all sign results from multiplying this factor with thesign associated with each matrix element.7. Each diagram comes with a numerical prefactor ob-tained from the following combination(a) A factor / ( n e )! must be considered for each group of n e equivalent lines. Equivalent linesmust all begin and end at the same vertices (orvertex, for anomalous propagators starting andending at the same vertex), and must correspondto the same type of contractions, i.e., they mustall correspond to propagators characterized bythe same superscripts g and g ′ in addition tohaving identical time labels.(b) Given the previous rule, an extra factor / mustbe considered for each anomalous propagatorthat starts and ends at the same vertex.(c) A symmetry factor /n s must be considered inconnection with exchanging the time labels ofthe vertices in all possible ways, counting theidentity as one. The factor n s corresponds tothe number of ways exchanging the time labelsprovides a time-labelled diagram that is topolog-ically equivalent to the original one.In order to illustrate the typical expression of off-diagonal Feynman BMBPT diagrams and to anticipateseveral key characteristics, let us compute the three second-order diagrams displayed in Fig. 6, i.e.,PO2.2.1 = 14 X k i Ω k k Ω k k k k O k k In case a line is drawn such that it crosses itself, the crossing(s)must be omitted when evaluating p . × lim τ →∞ Z τ d τ d τ θ ( τ − τ ) (38a) × e − τ ǫ k k e − τ ǫ k k k k PO2.2.2 = 12 X k i Ω k k Ω k k k k O k k R −− k k ( ϕ ) × lim τ →∞ Z τ d τ d τ θ ( τ − τ ) (38b) × e − τ ǫ k k e − τ ǫ k k k k PO2.2.3 = 14 X k i Ω k k Ω k k k k O k k R −− k k ( ϕ ) R −− k k ( ϕ ) × lim τ →∞ Z τ d τ d τ (38c) × e − τ ǫ k k e − τ ǫ k k k k , where the extended notation ǫ k a k b ...k i k j ... ≡ E k i + E k j + . . . − E k a − E k b − . . . , (39)was introduced. In each case, the sign, the combinatorialfactors and the three matrix elements directly reflect Feyn-man’s algebraic rules listed above and are easy to interpret.Eventually, the final form of the integrand originates fromexpliciting the n b = 4 propagators via Eq. (36), whichinduces the presence of one off-diagonal elementary con-tractions per anomalous propagator.The three chosen diagrams display the same overalltopology , i.e., while the vertex O is at fixed time ,the vertex Ω / Ω / Ω is at running time τ and the Ω vertex is at running time τ . However, the three diagramsdiffer in their number of anomalous lines and, as such,clearly illustrate key consequences of going from diagonalto off-diagonal BMBPT. The first diagram, PO2.2.1, con-tains no anomalous line ( n a = 0 ) and already occurs instraight, i.e., diagonal, BMBPT . By turning the secondvertex Ω into Ω ( Ω ), PO2.2.2 (PO2.2.3) contains n a = 1 ( n a = 2 ) anomalous line(s) between the secondand the third vertices. As a consequence, the integrandsdisplay typical structures that need to be scrutinized forthe following. • The fact that the two running variables τ and τ arepositive is directly encoded into the boundary of thedouble integral. • In PO2.2.1, the explicit step function characterizesthe time ordering induced between Ω and Ω ver-tices by the two normal propagators connecting them.This step function, i.e., time ordering, remains at playin PO2.2.2 given than one normal line still connectsthe second and third vertices. Contrarily, the absence The number of quasi-particle indices on which summation is per-formed increases by one per anomalous propagator due to the factthat the matrix R −− ( ϕ ) is not diagonal in quasi-particle space. This diagram is the one denoted as PO2.2 in Fig. 6 of Ref. [1]. τ τ PO2.2.1 k k k k +Ω k k + O k k +Ω k k k k k k τ τ PO2.2.2 k k k +Ω k k + O k k +Ω k k k k k k k k τ τ PO2.2.3 k k +Ω k k + O k k +Ω k k k k Figure 6: Selected second-order off-diagonal Feynman BMBPT diagrams. of step function in PO2.2.3 characterizes the fact that Ω and Ω are solely connected via anomalous prop-agators that do not induce any time-ordering relationbetween them. While in the first two cases the in-tegral over τ depends on the integral over τ , bothintegrals are independent from each other in PO2.2.3. • Grouping appropriately the exponential functionscoming from the four propagators, the integrand dis-plays one exponential factor per running time, i.e., per Ω i k j k ( τ k ) vertex. The relevant energy factor ǫ k a k b ...k i k j ... multiplying the variable τ k in this exponential func-tion denotes the sum/difference of quasi-particle en-ergies associated with the lines entering/leaving thecorresponding vertex.The three diagrams exemplify the fact that the off-diagonalBMBPT diagrammatics differentiates itself by the pres-ence of anomalous lines that, depending on the situation,may change the time-ordering structure between the ver-tices compared to the diagonal BMBPT diagram display-ing the same overall topology. The expression obtained via the application of Feyn-man’s algebraic rules does not yet constitute the formneeded for the numerical implementation of the formalism.While the sign, the combinatorial factor and the matrixelements will remain untouched, the p -tuple time integralmust be performed in order to obtain the needed expres-sion.A major part of Ref. [1] was dedicated to the automatedcomputation of the p -tuple time integrals via the introduc-tion of the so-called time-structure diagram (TSD) under-lying any given BMBPT diagram of arbitrary order andtopology. We thus refer to Ref. [1] for the general theory ofTSDs and will come back later on to its implementation formore general off-diagonal BMBPT diagrams. For now, itis sufficient to focus on the main consequence of the aboveanalysis, i.e., while the presence of one anomalous line inPO2.2.2 does not change its time structure compared toPO2.2.1, turning the other propagator connecting the sec-ond and third vertices into an anomalous line does modifyit. Consequently, while the TSD associated to diagrams PO2.2.1 and PO2.2.2 is T2.1 (see Fig. 12), it becomes T2.2for PO2.2.3. Generically denoting as a k the energy factormultiplying the time label τ k in the integrand, the inte-grals associated with the examples given in Eq. (38) areT2.1 = lim τ →∞ Z τ dτ dτ θ ( τ − τ ) e − a τ e − a τ = 1 a ( a + a ) , (40a)T2.2 = lim τ →∞ Z τ dτ dτ e − a τ e − a τ = 1 a a , (40b)the first (second) of which applies to PO2.2.1 and PO2.2.2(PO2.2.3).In order to obtain the final, i.e. time-integrated, expres-sion of each of the three diagrams, the factors a and a must be expressed back in terms of quasi-particle energies.As discussed in Ref. [1] for diagonal BMBPT diagrams,and as generalized to off-diagonal BMBPT diagrams be-low, the specific combinations of these factors emergingfrom the TSDs correspond necessarily to positive sumsof quasi-particle energies that can be straightforwardlyextracted from the diagram itself. Combining Eqs. (38)and (40) before inserting the appropriate combinations ofquasi-particle energies, one eventually obtains the desiredexpressions under the formPO2.2.1 = 14 X k i Ω k k Ω k k k k O k k ǫ k k ǫ k k , PO2.2.2 = 12 X k i Ω k k Ω k k k k O k k ǫ k k k k ǫ k k R −− k k ( ϕ ) , PO2.2.3 = 14 X k i Ω k k Ω k k k k O k k ǫ k k k k ǫ k k R −− k k ( ϕ ) R −− k k ( ϕ ) . Off-diagonal BMBPT diagrams of order p = 0 and have been generated and evaluated manually for two-body11O0.1.1 (1) n a = O PO0.1.1 (2) n a = O PO0.1.1 (3) n a = O n a = a = PO1.1.1 (1) τ Ω O PO1.1.1 (2) τ Ω O PO1.1.1 (3) τ Ω O PO1.1.1 (4) τ Ω O PO1.1.1 (5) τ Ω O PO1.1.1 (6) τ Ω O PO1.1.2 (1) τ Ω O PO1.1.2 (2) τ Ω O PO1.1.2 (3) τ Ω O PO1.1.2 (4) τ Ω O PO1.1.2 (5) τ Ω O PO1.1.2 (6) τ Ω O PO1.2.1 (1) τ Ω O PO1.2.1 (2) τ Ω O PO1.2.1 (3) τ Ω O PO1.2.1 (4) τ Ω O PO1.2.1 (5) τ Ω O Figure 7: Zero- and first-order off-diagonal
Feynman BMBPT diagrams. Diagrams are grouped vertically according to the number n a ofanomalous lines they contain. operators, i.e., operators summing O [ k ] with k ≤ [9]. Thetwenty corresponding diagrams are displayed in Fig. 7 forillustration. Among these twenty diagrams, only the threediagrams in the first column ( n a = 0 ) remain in straight,i.e., diagonal, BMBPT that has been dealt with in the firstversion of the ADG code [1].While it was already challenging to automatically gen-erate and evaluate diagonal BMBPT diagrams of arbi-trary orders and topologies, off-diagonal BMBPT obvi-ously reaches yet another level of complexity related to the proliferation of diagrams, itself increasing with theperturbative order, associated with the possibility to formoff-diagonal propagators. In this context, the step accom-plished in Ref. [1] will however happen to be of tremendoushelp to achieve the automatic generation and evaluation ofthe off-diagonal BMBPT diagrams.
5. Generation of off-diagonal BMBPT diagrams
The automated generation of diagonal BMBPT Feyn-man diagrams via elements of graph theory was explained12t length in Ref. [1]. We do not repeat it here and refer thereader to Ref. [1] for details. As a matter of fact, the strat-egy presently employed is not to follow a similar methodto generate off-diagonal diagrams from scratch but ratherto take advantage of having already done so for the diago-nal ones, i.e., to start from the order p diagonal BMBPTdiagrams to systematically produce their off-diagonal part-ners. Given that diagonal BMBPT diagrams constitute thebase line for generating the off-diagonal ones, the elevenzero-, first- and second-order diagonal BMBPT diagramsgenerated from operator vertices containing four legs atmost are displayed in Fig. 8 for reference. One recog-nizes in particular the three zero- and first-order diagramsPO0.1, PO1.1 and PO1.2 already appearing in Fig. 7 witha slightly different denomination whose aim is to group alldiagrams originating from the same diagonal diagram.Diagonal and off-diagonal diagrams of order p derivefrom the same many-body matrix element in Eq. (34), theemergence of both categories at once being authorized bythe fact that the ket is gauge rotated. The latter featureleads to the necessity to consider contractions betweenpairs of quasi-particle annihilation operators in additionto only contracting one creation and one annihilation op-erators in the strict diagonal limit. Starting from diago-nal BMBPT diagrams of order p , the complete set of off-diagonal diagrams can thus be obtained via the applicationof two basic operations1. adding self-contractions to each vertex, while chang-ing the nature of the vertex accordingly, until the sumof lines entering and leaving the vertex is equal to therank deg_max of the associated operator.2. incrementally changing normal propagators linkingtwo vertices into anomalous ones. This is achievedby turning the arrow associated with an outgoing linein the original propagator into an incoming line, thusmodifying the concerned vertex accordingly.Let us now exemplified the two above operations that musteventually be applied systematically.Considering the zero-order diagonal diagram PO0.1.1 (1) in Fig. 7 (i.e., PO0.1 in Fig. 8), and a two-body operator O ( deg_max = 4 ), the vertex O has no line entering orleaving it. Replacing it by O and O , one generates twovalid off-diagonal diagrams containing one and two anoma-lous contractions denoted as PO0.1.1 (2) and PO0.1.1 (3) in Fig. 7. One can proceed similarly starting from thefirst-order diagram denoted as PO1.1.1 (1) in Fig. 7 (i.e.,PO1.1 in Fig. 8). Adding a self-contraction to each ofthe two vertices provides three additional off-diagonal dia-grams containing one or two anomalous lines and denotedas PO1.1.1 (4) , PO1.1.2 (1) and PO1.1.2 (4) in Fig. 7.To illustrate the second operation, let us consider thefirst-order diagonal diagram PO1.2.1 (1) in Fig. 7 (i.e.,PO1.2 in Fig. 8). This diagram contains four normal O PO0.1 O Ω τ O Ω τ PO1.1 PO1.2 τ τ Ω O ˘Ω τ τ Ω O Ω τ τ Ω O Ω PO2.1 PO2.2 PO2.3 τ τ Ω O Ω τ τ Ω O Ω τ τ ˘Ω O Ω PO2.4 PO2.5 PO2.6 τ τ Ω O Ω τ τ Ω O Ω PO2.7 PO2.8
Figure 8: Zero-, first- and second-order diagonal
Feynman BMBPTdiagrams generated from operator vertices containing four legs atmost, i.e., with deg_max = 4. lines between the two vertices, each of which can be trans-formed into an anomalous line. Doing so generates fouradditional topologically-distinct off-diagonal diagrams de-noted as PO1.2.1 (2) , PO1.2.1 (3) , PO1.2.1 (4) and PO1.2.1 (5) in Fig. 7.Combining the transformation of normal lines intoanomalous lines and the addition of self-contractions, oneobtains the other topologically distinct off-diagonal dia-grams displayed in Fig. 7.13O0.1.1 ˜ O ( ϕ ) PO1.1.1 Ω ˜ O ( ϕ ) PO1.1.2 Ω ˜ O ( ϕ ) PO1.2.1 Ω ˜ O ( ϕ ) Figure 9: Zero- and first-order effective off-diagonal BMBPT dia-grams recasting the twenty displayed in Fig. 7.
An important feature is that the bottom vertex O m appearing in diagonal BMBPT diagrams is always at fixedtime zero. Consequently, off-diagonal diagrams generatedby adding self-contractions to it and/or by transforminga normal line leaving it into an anomalous line entering itpossess the same time structure as the diagonal diagramit derives from. Indeed, a self-contraction carries no timedependence and thus cannot impact the time structure ofthe diagram. Furthermore, the fact that all Ω ij vertices areat higher times than O m remains true even if all the linesattached to the bottom vertex are changed into anomalousones such that the time structure is invariant under thistransformation.A key consequence of the above observations is that alldiagrams differing only by the number of self-contractionsonto the bottom vertex and/or the number of anomalouspropagators connected to it can be grouped into a singlediagram in which the bottom vertex is replaced by its sim-ilarity transformed partner at gauge angle ϕ [9] ˜ O ( ϕ ) ≡ e − Z ( ϕ ) Oe Z ( ϕ ) , (42)where Z ( ϕ ) is the Thouless operator, see Appendix A.3.As explained at length in Appendix B, the transformedoperator ˜ O ( ϕ ) possesses the same formal structure as theinitial operator O . As such, it is decomposed as a sum ofterms ˜ O mn ( ϕ ) with the same overall rank as O , i.e., m + n ≤ deg_max . The only difference relates to the definitionof the (gauge-dependent) matrix elements entering eachterm ˜ O mn ( ϕ ) . The expression of these matrix elementsin terms of the original ones were provided in Ref. [9] for deg_max = 4 .Exploiting this key observation, one can reduce drasti-cally the number of effective diagrams at play. Employingthe transformed operator for the bottom vertex and for-bidding any anomalous line to connect to it, the twentyoff-diagonal diagrams displayed in Fig. 7 are recasted intothe four effective off-diagonal diagrams displayed in Fig. 9.This feature being generic, the recasting procedure extendsto any order p . The analysis provided above puts us in position tostate the systematic rules used to generate all effectiveoff-diagonal BMBPT diagrams of order p from the diag-onal ones. Starting from a diagonal BMBPT diagram oforder p
1. replace the bottom vertex O m by its transformedpartner ˜ O m ( ϕ ) ,2. for each energy vertex Ω ij (a) transform l outgoing arrows into incoming ar-rows to form l anomalous lines while turning thevertex into Ω i − lj + l , with l ∈ N , until i − l = 0 ,(b) add k self-contractions while turning the vertexinto Ω ij +2 k , with k ∈ N , until i + j + 2 k = deg_max ,3. retain only topologically distinct diagrams.While the method is straightforward, it is indeed im-portant to discard topologically equivalent diagrams gen-erated through this brute-force procedure. Anticipating it,one can actually reduce the need to check for all of them,which is particularly beneficial given that the correspond-ing test scales factorially with the number of vertices in thediagrams. The first step in our algorithm being to generateoff-diagonal diagrams from a given straight BMBPT one,one avoids producing topologically equivalent diagrams inthe subset of children in the following way: • As normal propagators are turned anomalous, thelist of initial equivalent lines is recorded takinginto account possible vertex exchanges. Doing so,only unique permutations of normal propagators areturned into anomalous ones. • Once this first set of diagrams is generated, addingself-contractions on them can only create topologicallyequivalent diagrams when doing so on topologicallyequivalent vertices. As such, one first looks out forsuch vertices, and makes sure to add self-contractionson unique combinations of vertices.Once order- p diagrams are generated in this way, checksfor topologically equivalent diagrams are still required.This is due to the fact that by turning normal propaga-tors anomalous, one creates new topologies that might infact arise from several straight BMBPT diagrams. Somediagrams can be excluded from this check though, start-ing from the ones containing only normal propagators, i.e.,straight BMBPT ones, which have been tested beforehand.This exclusion can be extended to diagrams in which theonly anomalous lines are self-contractions, i.e., straightBMBPT diagrams with extra anomalous self-contractions,as the structure of propagators exchanged between verticeshas not been affected. Applying the discussed methodeventually results in the number of diagrams displayed inTab. 1. Once the off-diagonal BMBPT diagrams have been pro-duced, it is important to be able to represent them graph-ically. As in Ref. [1], the BMBPT diagrams are created asobjects generated via the Python package
NetworkX [26],allowing for an efficient and easy to handle storage of thenecessary information. It is then easy to read these objects14rder 0 1 2 3 4Straight BMBPT deg_max = 4 deg_max = 6 deg_max = 4 deg_max = 6
Table 1: Number of diagonal and effective off-diagonal
BMBPT diagrams generated from operators containing at most four ( deg_max = 4)or six ( deg_max = 6) legs. and extract the information necessary to produce the draw-ing instructions in the form of
FeynMP [27] commands.The routines designed in Ref. [1] only had to be adaptedto allow the drawing of anomalous propagators and self-contractions. As an example, the output displaying thedrawing instructions of the off-diagonal BMBPT diagramdisplayed in Fig. 14 is given in Fig. 10.
6. Evaluation of off-diagonal BMBPT diagrams
Having the capacity to generate all off-diagonal BMBPTFeynman diagrams of order p , the next challenge is to sys-tematically derive their expression. Doing so on the ba-sis of Feynman’s algebraic rules is rather straightforward.However, it leaves the p -tuple time integral to performin order to obtain the time-integrated expression of in-terest. In Ref. [1], an algorithm was found to overcomethis challenge for diagonal BMBPT diagrams at play instraight BMBPT without prior knowledge of the pertur-bative order or of the topology of the diagram. This even-tually led to the identification of a novel diagrammatic rule.The present section details how the method only needs tobe slightly generalized in order to realize the same objec-tive for off-diagonal
BMBPT diagrams at play in PNP-BMBPT.
Obtaining the result of p -tuple time integrals in an au-tomatic fashion was made possible via the introduction ofthe time-structure diagram underlying any given diagonalBMBPT diagram of arbitrary order and topology. We re-fer to Ref. [1] for the general theory of TSDs and onlycomment here on the specificities encountered when deal-ing with more general off-diagonal BMBPT diagrams.The key point was already alluded to in Sec. 4.5 andrelates to the impact anomalous lines may have on theTSD attributed to a given off-diagonal BMBPT diagram.The main features are • The running time labels ( τ , . . . , τ p ) are positive suchthat each Ω vertex entertains at least an ordering re-lation with the bottom vertex ˜ O ( ϕ ) independently ofthe network of lines running through the BMBPT di-agram. Consequently, the TSD remains necessarilyconnected, independently of its topology. • Contrarily to normal lines, anomalous lines do notinduce any time ordering relation. This means that,while two Ω vertices connected by at least one normalline are time ordered, it is not the case if they aresolely connected via anomalous propagators. Conse-quently, a link connecting two Ω vertices in the TSDassociated to a diagonal BMBPT diagram will disap-pear when the two vertices become only connected viaanomalous propagators in an off-diagonal partner di-agram . Whenever an Ω vertex ends up entertainingno time relation with any other due to the replace-ment of normal lines by anomalous ones, it becomesdirectly linked to the bottom vertex ˜ O ( ϕ ) in the asso-ciated TSD. • The addition of a self contraction to any given Ω ver-tex does not change the time structure of the diagramand thus the associated TSD.In conclusion, the presence of anomalous lines may, de-pending on the situation, change the TSD associated toan off-diagonal BMBPT diagram compared to the diago-nal diagram displaying the same topology. Eventually, theTSD associated to an off-diagonal BMBPT diagram canbe obtained through the following steps1. copy the off-diagonal BMBPT diagram,2. remove all the anomalous propagators,3. replace the normal propagators by links,4. add a link between the bottom vertex at time 0 andevery other vertex if such a link does not exist,5. for each pair of vertices, consider all paths linkingthem and only retain the longest one,6. match the label a q associated to a given vertex in theTSD diagram to the sum/difference of quasi-particleenergies associated with the lines entering/leaving thecorresponding vertex in the BMBPT diagram.The only difference with the procedure followed for strictlydiagonal BMBPT diagrams [1] relates to step 2 that triv-ially stipulates to strip off anomalous propagators if any. As a TSD stores only the minimal information associated withtime-ordering relations, such a disappearance may parallel the oc-curence of one or several new links with respect to the TSD of thediagonal BMBPT diagram. Hence TSDs must always be producedstarting from a given off-diagonal BMBPT diagram and not fromthe TSD of its parent BMBPT diagram. begin{fmffile}{offdiag_ex}\begin{fmfgraph*}(80,80)\fmfcmd{style_def prop_pm expr p =draw_plain p;shrink(.7);cfill (marrow (p, .25));cfill (marrow (p, .75))endshrink;enddef;}\fmfcmd{style_def prop_mm expr p =draw_plain p;shrink(.7);cfill (marrow (p, .75));cfill (marrow (reverse p, .75))endshrink;enddef;}\fmfcmd{style_def half_prop expr p =draw_plain p;shrink(.7);cfill (marrow (p, .5))endshrink;enddef;}\fmfstraight\fmftop{d1,d2,v2,d3,d4}\fmfbottom{v0}\fmf{phantom}{v0,v1}\fmfv{d.shape=square,d.filled=full,d.size=3thick}{v0}\fmf{phantom}{v1,v2}\fmfv{d.shape=circle,d.filled=full,d.size=3thick}{v1}\fmfv{d.shape=circle,d.filled=full,d.size=3thick}{v2}\fmffreeze\fmf{prop_pm,left=0.5}{v0,v1}\fmf{prop_pm,right=0.5}{v0,v1}\fmf{prop_pm,left=0.5}{v1,v2}\fmf{prop_mm,right=0.5}{v1,v2}\fmf{half_prop,right}{d3,v2}\fmf{half_prop,left}{d3,v2}\end{fmfgraph*}\end{fmffile} Figure 10:
FeynMP instructions to draw the off-diagonal BMBPTdiagram displayed in Fig. 14.
The procedure is illustrated in Fig. 11 for a third-order diagonal BMBPT diagram and for the particularoff-diagonal diagram generated from it by turning thetwo normal lines connecting vertex Ω to one of the two Ω vertices into anomalous lines. Cleared of other in-formations, the TSDs tranparently characterize the time-ordering structure underlying the diagrams. In the firstone, the three Ω ij vertices are at higher times than O such that the two Ω vertices are at higher times than Ω without being ordered with respect to one another.From the graph theory viewpoint, the corresponding TSDis a tree, i.e., it contains no cycle, with two branches. Inthe second diagram, the fact that the two lines connecting Ω to (one of the two) Ω are anomalous relaxes the time-ordering between both vertices and, as a result, changesthe nature of the associated TSD, i.e., a is now directlylinked to the bottom vertex →→ → a a a →→ → a a a Figure 11: Production of the TSDs associated with a third-order di-agonal BMBPT diagram and with an off-diagonal diagram obtainedfrom it by turning two among the eight normal lines into anomalousones.
It is mandatory to generate the TSD from its underly-ing BMBPT diagram. Indeed, only in this case can therank deg_max of the operators at play be employed to con-strain the topology of the diagrams, eventually dictatingthe topology of allowed TSDs. Furthermore, going fromdiagonal to off-diagonal BMBPT diagrams may not onlychange the nature of the TSD associated to a particulardiagram but also increase the list of active TSDs at a givenorder.With this in mind and following the above rules,the / / / TSDs of order / / / associated to off-diagonal BMBPT diagrams generated from operators with deg_max = 4 or deg_max = 6 have been produced and The two TSDs appearing in Fig. 11 are denoted, respectively, asT3.3 and T3.2 (with a cyclic permutation of ( a , a , a )) in Fig. 12. a T1.1 a a a a T2.1 T2.2 a a a a a a a a a a a a a a a T3.1 T3.2 T3.3 T3.4 T3.5
Figure 12: Zero-, first-, second- and third-order TSDs correspondingto off-diagonal BMBPT diagrams generated from operators contain-ing four or six legs at most, i.e., with deg_max = 4 or deg_max = 6. systematically displayed in Fig. 12. Interestingly, restrict-ing one-self to diagonal BMBPT diagrams and deg_max =4 , T3.4 would have to be removed from the set of activeTSDs, i.e., going from diagonal to off-diagonal BMBPT orgoing from deg_max = 4 to deg_max = 6 adds one allowedthird-order TSD. In the end, different
BMBPT diagrams of order p canhave the same TSD, i.e., the same underlying time struc-ture. At the same time, off-diagonal BMBPT diagramsoriginating from the same diagonal diagram may have different
TSDs, i.e., different underlying time structures.Once the TSD has been extracted from the BMBPT di-agram of interest, its computation follows the algorithmdetailed in Ref. [1]. In particular, the treatment of a non-tree TSD requires to turn it first into a sum of tree TSDs.Once the expression of a tree TSD of order p hasbeen obtained, the goal is to generate the actual time-integrated expression of the BMBPT diagrams associatedto it. Rather than replacing the individual factors a q , q = 1 , . . . , p , by their expressions for each BMBPT dia-gram, one introduces the notion of subdiagram, or sub-graph, to directly obtain their combinations appearing inthe denominator of the time-integrated expression of in-terest. In Ref. [1], a subdiagram of a diagonal BMBPTdiagram was defined as a diagram composed by a subsetof vertices plus the propagators that are exchanged be-tween them. As each vertex label a q in a TSD eventuallystands for the sum/difference of quasi-particle energies as-sociated with the lines entering/leaving the vertex in theassociated BMBPT diagram, the combination of these la-bels denotes the sum/difference of quasi-particle energiesassociated with the lines entering/leaving the subdiagram grouping the corresponding vertices.In the present context of off-diagonal BMBPT diagrams,one only needs to slightly generalize the notion of subdia-grams in order to apply the algorithm stipulated in Ref. [1]in order to determine the factors entering the denominatorof the time-integrated expression of the diagram. Thus asubdiagram of an off-diagonal BMBPT diagram is now de-fined as a diagram composed by a subset of vertices plusthe normal propagators that are exchanged between them.This definition is obviously consistent with the one intro-duced in Ref. [1] for strictly diagonal BMBPT diagramsgiven that the latter solely contain normal propagators.With this definition at hand, the energy denominatorresulting from a BMBPT diagram associated with a treeTSD is obtained in the following way1. Consider a vertex but the bottom one in the BMBPTdiagram,(a) determine all its descendants using the TSD,(b) form a subdiagram using the vertex and its de-scendants,(c) sum the quasi-particle energies corresponding tothe lines entering the subdiagram,(d) add the corresponding factor to the denominatorexpression,2. Go back to 1. until all vertices have been exhausted.Given that anomalous lines are excluded from the defini-tion of a subdiagram, they systematically count as enteringthe subdiagram whenever they connect to a vertex belong-ing to it.Let us illustrate the diagrammatic rule by focusingon the two third-order BMBPT diagrams displayed inFig. 13.1. The denominator in the time-integrated expression ofthe first diagram is obtained through the followingsteps(a) The vertex at time τ in the BMBPT diagramcorresponds to vertex a in the TSD. Its de-scendants are vertices a and a corresponding The results obtained in Eq. (41) for the three second-order di-agrams displayed in Fig. 6 are straightforwardly recovered by theapplication of the diagrammatic rule. k k k k k k k τ τ τ a a a k k k k k k k k k k τ τ τ a a a Figure 13: Fully-labelled third-order diagonal and off-diagonalBMBPT diagrams displayed in Fig. 11 and their associated TSDs.The bottom vertex corresponds to ˜ O ( ϕ ). to BMBPT vertices at times τ and τ , respec-tively. The sum of quasi-particle energies associ-ated to the lines entering the subgraph groupingthe three vertices is ǫ k k k k , thus, providing thefirst factor entering the denominator.(b) The vertex at time τ in the BMBPT diagramcorresponds to vertex a in the TSD. It has no de-scendant such that the corresponding subgraphreduces to itself. The sum of quasi-particle ener-gies associated to the lines entering the subgraphis ǫ k k k k , thus providing the second factor en-tering the denominator.(c) The vertex at time τ in the BMBPT diagramcorrespond to vertex a in the TSD. It has no de-scendant such that the corresponding subgraphreduces to itself. The sum of quasi-particle ener-gies associated to the lines entering the subgraphis ǫ k k k k , thus, providing the last factor enter-ing the denominator.(d) Eventually, the complete denominator reads as ǫ k k k k ǫ k k k k ǫ k k k k , where each factor corresponds to a positive sumof quasi-particle energies.2. The denominator of the second, off-diagonal, diagramcontaining two anomalous lines and corresponding toa different TSD is obtained as(a) The vertex at time τ in the BMBPT diagramcorresponds to vertex a in the TSD. Contrarilyto the previous case, vertex a is not a descen-dant of a anymore as is visible from the TSD such that the subgraph of interest solely groups a and a . The sum of quasi-particle energiesassociated to the lines entering the subgraph inthe BMBPT diagram is ǫ k k k k , thus, provid-ing the first factor entering the denominator.(b) The situation of the vertex at time τ in the off-diagonal BMBPT diagram is strictly the same asin the previous diagonal one. Consequently, theassociated factor in the denominator is ǫ k k k k .(c) As in the diagonal BMBPT diagram, the ver-tex at time τ has no descendant in the off-diagonal BMBPT diagram. Consequently, thesubgraph corresponding to vertex a reduces toitself. However, because the two anomalous linescarry two quasi-particle labels each, the sum ofquasi-particle energies associated to the lines en-tering the subgraph has now become ǫ k k k k .(d) Eventually, the complete denominator reads as ǫ k k k k ǫ k k k k ǫ k k k k , where each factor corresponds to a positive sumof quasi-particle energies.For completeness, let us work out another example high-lighting additional features of interest, i.e., the second-order off-diagonal diagram displayed in Fig. 14 togetherwith its associated TSD . Applying the diagrammaticrule, one obtains1. The vertex at time τ in the BMBPT diagram corre-sponds to vertex a in the TSD. Because it remainsone normal line connecting to the vertex at time τ , a is indeed its descendant. The subgraph of inter-est thus groups a and a . Due to the more generaldefinition of subgraphs at play in the context of off-diagonal BMBPT, the anomalous line connecting thetwo vertices is excluded from it, together with the selfcontraction on the upper vertex. Consequently, thesum of quasi-particle energies associated to the linesentering the subgraph is ǫ k k k k k k , thus, providingthe first factor entering the denominator.2. The vertex at time τ in the BMBPT diagram cor-responds to vertex a in the TSD. It has no descen-dant such that the corresponding subgraph reducesto itself, excluding the self contraction that the ver-tex exchanges with itself. The sum of quasi-particleenergies associated to the lines entering the subgraphis ǫ k k k k , thus providing the second factor enteringthe denominator.3. Eventually, the complete denominator reads as ǫ k k k k k k ǫ k k k k . It is worth noting that the TSD of the off-diagonal diagram of in-terest is unchanged compared to the diagonal diagram it is generatedfrom. However, the companion diagram with one more anomalousline joining the first and second Ω vertices relates to a different TSD. k k k k k k τ τ a a Figure 14: Fully-labelled second-order off-diagonal BMBPT diagramand its associated TSD.
7. Output of the
ADG program
The typical output of the code associated with an off-diagonal BMBPT diagram is:
Diagram 3.3:. PO . . τ →∞ ( − X k i ˜ O k k k k ( ϕ )Ω k k k k Ω k k k k × R −− k k ( ϕ ) R −− k k ( ϕ ) × Z τ d τ d τ θ ( τ − τ ) e − τ ǫ k k k k e − τ ǫ k k k k = ( − X k i ˜ O k k k k ( ϕ )Ω k k k k Ω k k k k ǫ k k k k ǫ k k k k × R −− k k ( ϕ ) R −− k k ( ϕ ) → T1:T a a a = ǫ k k k k a = ǫ k k k k
8. Connection to time-ordered diagrammatics
In Ref. [1], time-unordered and time-ordered diagram-matics emerging, respectively, from the time-dependentand the time-independent formulations of straight, i.e., di-agonal, BMBPT were compared at length. The main out-come of the analysis related to the capacity of the time-unordered diagrammatics to resum at once large classesof time-ordered diagrams. Correspondingly, it was shownthat the new diagrammatic rule allowing for the direct ob-tention of the time-integrated results on the basis of time-unordered diagrams generalizes the resolvent rule at playfor time-ordered diagrams.As in Ref. [1], the formal and numerical developmentspresented in this paper rely on the time-dependent formu-lation of PNP-BMBPT [9]. While it is traditionally more customary to design many-body perturbation theories onthe basis of a time-independent formalism [28], this taskhas so far not been attempted for PNP-BMBPT. Whilethe end result has to be the same, the partitioning ofthe complete order- p contribution to the observable O A0 will differ in both approaches. In the absence of the time-ordered diagrammatics associated to PNP-BMBPT, wecannot proceed to the same analysis as done in Ref. [1]for straight BMBPT. Leaving this analysis for the future,it can however be anticipated that time-unordered off-diagonal BMBPT diagrams will feature the same capacityto resum large classes of time-ordered diagrams at play inthe, yet-to-be-formulated , time-independent version ofPNP-BMBPT.
9. Use of the
ADG program
ADG has been designed to work on any computer witha Python3 distribution, and successfully tested on recentGNU/Linux distributions and on MacOS. Additionally toPython, setuptools and distutils packages must already beinstalled, which is the case on most standard recent distri-butions. Having pip installed eases the process but is nottechnically required. The
NumPy , NetworkX and
SciPy libraries are automatically downloaded during the installprocess. Additionally, one needs a L A TEX distribution in-stalled with the PDFL A TEX compiler for
ADG to producethe pdf file associated to the output if desired.
The easiest way to install ADG is to obtain it from thePython Package Index by entering the following com-mand pip3 install adg Provided setuptools is already installed, pip takes care ofdownloading and installing
ADG as well as
NumPy , SciPy and
NetworkX . Once a new version of
ADG is released, onecan install it by entering the command pip3 install --upgrade adg A valid partitioning relates to splitting the complete order p in asum of terms that are individually proportional to a fraction of theform 1 / ( ǫ k i ...k j . . . ǫ k u ...k v ) with p energy factors in the denominator.Any other form does not constitute a valid partitioning in the presentcontext. The fact that anomalous propagators/contractions are not diag-onal in their quasi-particle indices should lead to a rather unconven-tional time-ordered diagrammatics that shall itself lead to an inter-esting variant of the resolvent rule. As the previous version of ADG was developped with Python2,the installation will be made next to the previous version, thoughit should not cause any problem using it. Users who would wantto uninstall the previous version first might do it by entering pip2uninstall adg . https://pypi.org/project/adg/ .1.2. From the source files Once the
ADG source files are downloaded from the CPClibrary or the GitHub repository , one must enter theproject folder and either run pip3 install . or python3 setup.py install With this method, pip also takes care of downloadingand installing NumPy , NetworkX and
SciPy . The most convenient way to use
ADG is to run it in batchmode with the appropriate flags. For example, to run theprogram and generate off-diagonal BMBPT diagrams atorder 3 for example, one can use adg -o 3 -t PBMBPT -d -c where the -o flag is for the order, -t for the type of theory, -d indicates that the diagrams must be drawn and -c that ADG must compile the L A TEX output. A complete list ofthe program’s options can be obtained via the program’sdocumentation (see Sec. 9.4) or through adg -h
Currently,
ADG can be run either in relation to HF-MBPT by using -t MBPT , to straight BMBPT by using -t BMBPT or to off-diagonal BMBPT by using -t PBMBPT .Though the algorithms described in the previous sectionscan be used regardless of the diagrams’ orders,
ADG hasbeen arbitrarily restricted to order 10 or lower to avoidmajor overloads of the system. Future users are never-theless advised to first launch calculations at low orders(2, 3 or 4 typically) as the time and memory needed forcomputations rise rapidly with the perturbative order.
As an alternative to the batch mode,
ADG can be run ona terminal by entering the command adg -i
A set of questions must be answered using the keyboardto configure and launch the calculation. The interactivemode then proceeds identically to the batch mode. https://github.com/adgproject/adg Depending on the system, it might be necessary either to use the "–user" flag to install it only for a specific user or to run the previouscommand with "sudo -H" to install it system-wide.
Let us briefly recapitulate the different steps of a typical
ADG run • Run options are set either by using the command-lineflags entered by the user or during the interactive ses-sion via keyboard input. • ADG creates a list of adjacency matrices for the appro-priate theory and perturbative order using
NumPy ,and feeds them to
NetworkX that creates
MultiDi-Graph objects. • Checks are performed on the list of graphs to removetopologically equivalent or ill-defined graphs. • The list of topologically unique graphs is used to pro-duce
Diagram objects that store the graph as well assome of its associated properties depending on the the-ory (HF status, excitation level, etc.). For off-diagonalBMBPT, this is done in two-steps, first genereatingthe straight BMBPT ones. The expression associatedto the graphs are eventually extracted. • The program prints on the terminal the number of di-agrams per category and writes the L A TEX output file,the details of which depend on the options selectedby the user, as well as a list of adjacency matricesassociated to the diagrams. Other output files maybe produced, depending on the theory and the user’sinput. • If asked by the user, the program performs thePDFL A TEX compilation. • Unnecessary temporary files are removed and the pro-grams exits.
Once the source files have been downloaded, a quickstart guide is available in the
README.md file. Once
ADG isinstalled, it is possible to read its manpages through man adg or a brief description of the program and its optionsthrough adg -h
A more detailed HTML documentation can be generateddirectly from the source files by going into the docs direc-tory and running make html
The documentation is then stored in docs/build/html ,with the main file being index.html . A list of other possi-ble types of documentation format is available by running make help .4.2. Online documentation The full HTML documentation is available on-line under https://adg.readthedocs.io/ , and helpwith eventual bugs of the program can be ob-tained by opening issues on the GitHub repository at https://github.com/adgproject/adg .
10. Conclusions
The motivations underlining our work were explained atlength in Ref. [1] where the first version of the code
ADG was described. The long-term rationale relates to the pos-sibility to automatically (i) generate and (ii) algebraicallyevaluate diagrams in various quantum many-body meth-ods of interest. In Ref. [1], the focus was put on Bogoli-ubov many-body perturbation theory (BMBPT) that hasbeen recently formulated [9] and first implemented at loworders [12] to tackle (near) degenerate Fermi systems, e.g.,open-shell nuclei displaying a superfluid character. Giventhe need to tackle three-nucleon interactions, i.e., six-legvertices, in nuclear physics and the implementation of high-order contributions authorized by the rapid progress ofcomputational power, the first version of the code
ADG makes possible to generate all valid BMBPT diagrams andto evaluate their algebraic expression to be implementedin a numerical application. This is realized at an arbitraryorder p for a Hamiltonian containing both two-body (four-legs) and three-body (six-legs) interactions (vertices). Theformal advances and the numerical methods necessary toachieve this goal can be found in Ref. [1].Bogoliubov MBPT perturbatively expands the exact so-lution of the Schrödinger equation around a so-called Bo-goliubov reference state, i.e., a general product state break-ing U (1) global-gauge symmetry associated with the con-servation of good particle number in the system. Thisresults into a (hopefully small) symmetry contaminationas soon as the expansion is truncated in actual calcula-tions. Given that the breaking of a symmetry cannot ac-tually be realized in a finite quantum system, U (1) sym-metry must eventually be restored at any truncation order,which is made possible thanks to the recent formulation ofthe particle-number projected Bogoliubov many-body per-turbation theory (PNP-BMBPT) [9] that extends straightBMBPT on the basis of a more general diagrammatic ex-pansion.Consequently, the present paper details the systematicgeneration and evaluations of diagrams at play in PNP-BMBPT operated by the second version (v2.0.0) of thecode ADG . While the automated evaluation of the diagramsonly requires a mere extension of the diagrammatic ruleunrevealed in Ref. [1], the method used to first generateall allowed diagrams is different from the one employed inRef. [1]. Taking advantage of the capacity of the code
ADG to already produce all valid BMBPT diagrams of order p , the set of rules to generate all those at play in PNP-BMBPT from those appearing in BMBPT was identifiedand implemented. Eventually, the second version of the code ADG is keptflexible enough to be expanded throughout the years totackle the diagrammatics at play in yet other many-bodyformalisms that either already exist or are yet to be for-mulated.
Acknowledgments
The authors thank M. Drissi for reporting a bug in theprevious version of the program. This publication is basedon work supported by the UK Science and TechnologyFacilities Council (STFC) through grants ST/P005314/1and ST/L005516/1, the Max Planck Society, the DeutscheForschungsgemeinschaft (DFG, German Research Founda-tion) – Project-ID 279384907 – SFB 1245, and the frame-work of the Espace de Structure et de réactions NucléairesThéorique (ESNT) at CEA.
Appendix A. Reference states
Appendix A.1. Bogoliubov vacuum
We consider the Bogoliubov reference state defined as | Φ i ≡ C Y k β k | i , (A.1)where quasi-particle operators are related to particle onesthrough a Bogoliubov transformation (cid:18) ββ † (cid:19) = W † (cid:18) cc † (cid:19) = (cid:18) U † V † V T U T (cid:19) (cid:18) cc † (cid:19) . (A.2)The product state | Φ i is a vacuum for the set of quasi-particle operators, i.e., β k | Φ i = 0 for all k . Appendix A.2. Rotated Bogoliubov vacuum
One introduces the gauge-rotated Bogoliubov vacuum | Φ( ϕ ) i ≡ R ( ϕ ) | Φ i ≡ C Y k ¯ β k | i , (A.3)where rotated quasi-particle operators are defined through (cid:18) ¯ β ¯ β † (cid:19) ( ϕ ) ≡ R ( ϕ ) (cid:18) ββ † (cid:19) R − ( ϕ ) (A.4) ≡ W ϕ † (cid:18) cc † (cid:19) , with the associated Bogoliubov transformation reading as W ϕ † ≡ (cid:18) U ϕ † V ϕ † V ϕ T U ϕ T (cid:19) (A.5) = (cid:18) U † e − iϕ V † e + iϕ V T e − iϕ U T e + iϕ (cid:19) . ppendix A.3. Thouless transformation between vacua As stipulated by Eq. (A.3), | Φ( ϕ ) i is obtained from | Φ i via the unitary transformation R ( ϕ ) whose generator is A .One can rather express | Φ( ϕ ) i via a non-unitary Thoulesstransformation of | Φ i according to [29] | Φ( ϕ ) i ≡ h Φ | Φ( ϕ ) i e Z ( ϕ ) | Φ i , (A.6)where the one-body Thouless operator Z ( ϕ ) ≡ X k k Z k k ( ϕ ) β † k β † k (A.7)only contains a pure excitation part over | Φ i . The corre-sponding Thouless matrix Z ( ϕ ) ≡ N ∗ ( ϕ ) M ∗− ( ϕ ) (A.8)is expressed in terms of the Bogoliubov transformation con-necting quasi-particle operators of | Φ( ϕ ) i to those of | Φ i (cid:18) ¯ β ¯ β † (cid:19) ( ϕ ) ≡ W † ( ϕ ) (cid:18) ββ † (cid:19) , (A.9)where W † ( ϕ ) = W ϕ † W ≡ (cid:18) M † ( ϕ ) N † ( ϕ ) N T ( ϕ ) M T ( ϕ ) (cid:19) (A.10) = (cid:18) U ϕ † U + V ϕ † V V ϕ † U ∗ + U ϕ † V ∗ V ϕ T U + U ϕ T V U ϕ T U ∗ + V ϕ T V ∗ (cid:19) . Appendix A.4. Elementary contractions
The elementary contractions of quasi-particle operatorsthat are in use when employing the off-diagonal Wick the-orem [25] are given by R ( ϕ ) ≡ h Φ | β † β | Φ( ϕ ) ih Φ | Φ( ϕ ) i h Φ | β β | Φ( ϕ ) ih Φ | Φ( ϕ ) ih Φ | β † β † | Φ( ϕ ) ih Φ | Φ( ϕ ) i h Φ | β β † | Φ( ϕ ) ih Φ | Φ( ϕ ) i ! ≡ (cid:18) R + − ( ϕ ) R −− ( ϕ ) R ++ ( ϕ ) R − + ( ϕ ) (cid:19) = (cid:18) − Z ( ϕ )0 1 (cid:19) . (A.11)Most of the above contractions are easily obtained by us-ing the fact that | Φ i is the vacuum of the quasi-particleoperators, i.e., β k | Φ i = h Φ | β † k = 0 for all k . The single non-trivial (anomalous) contraction is obtained on the basis ofstandard Wick’s theorem as R −− k k ( ϕ ) = h Φ | β k β k | Φ( ϕ ) ih Φ | Φ( ϕ ) i = h Φ | β k β k e Z ( ϕ ) | Φ i = 12 X kk ′ Z kk ′ ( ϕ ) h Φ | β k β k β † k β † k ′ | Φ i = 12 (cid:0) Z k k ( ϕ ) − Z k k ( ϕ ) (cid:1) = − Z k k ( ϕ ) , (A.12)and is zero in the diagonal case, i.e., R −− k k (0) = 0 . Appendix B. Transformed operator ˜ O ( ϕ ) The gauge-dependent similarity transformed operator of O is defined through ˜ O ( ϕ ) ≡ e − Z ( ϕ ) Oe Z ( ϕ ) . (B.1)Taking as an example one term in the normal-ordered ex-pression of O , e.g., O ij ≡ i ! 1 j ! X k ...k i + j O ijk ...k i + j β † k . . . β † k i β k i + j . . . β k i +1 , its transformed partner reads as ˜ O ( ij ) ( ϕ ) ≡ e − Z ( ϕ ) O ij e Z ( ϕ ) (B.2) = 1 i ! 1 j ! X k ...k i + j O ijk ...k i + j ˜ β † k . . . ˜ β † k i ˜ β k i + j . . . ˜ β k i +1 , where the transformed quasi-particle operators are ˜ β k ( ϕ ) ≡ e − Z ( ϕ ) β k e Z ( ϕ ) = β k − [ Z ( ϕ ) , β k ] + 12! [ Z ( ϕ ) , [ Z ( ϕ ) , β k ]] + . . . = β k + X k ′ Z kk ′ ( ϕ ) β † k ′ , (B.3a) ˜ β † k ( ϕ ) ≡ e − Z ( ϕ ) β † k e Z ( ϕ ) = β † k − [ Z ( ϕ ) , β † k ] + 12! [ Z ( ϕ ) , [ Z ( ϕ ) , β † k ]] + . . . = β † k , (B.3b)were use was made of the elementary commutators h β † k β † k ′ , β k i = β † k δ k ′ k − β † k ′ δ kk , (B.4a) h β † k β † k ′ , β † k i = 0 . (B.4b)Exploiting Eq. (B.3) and normal ordering the resultingterms with respect to | Φ i , the transformed operator inEq. (B.2) is eventually written as ˜ O ( ij ) ( ϕ ) ≡ i + j X m = i j X n =0 m + n ≤ i + j m ! 1 n ! × X k ...k m + n ˜ O mn ( ij ) k ...k m + n ( ϕ ) β † k . . . β † k m β k m + n . . . β k m +1 , (B.5) The Hermitian character of an operator O is lost by the applica-tion of similarity transformation. The notation ˜ O ( ij ) ( ϕ ) denotes the transformed operator of O ij such that the upper label ( ij ) is a sole reminder of the normal-orderednature of the original operator but does not characterize the normal-ordered nature of the transformed operator. Contrarily, ˜ O mn ( ϕ ) does denote the normal-ordered part of the transformed operator ˜ O ( ϕ ) of O containing m ( n ) quasi-particle creation (annihilation)operators. i ) asthe original operator O ij and possibly up to the total num-ber of original quasi-particle operators ( i + j ). The numberof annihilation operators ranges from 0 to the original num-ber ( j ) such that the overall number of quasi-particle oper-ators is bound to remain between i and i + j in each term.One notices that the only structural difference between theoriginal and the transformed normal-ordered operators re-lates to the fact that matrix elements of the latter dependon the gauge angle. Of course, the original operator isrecovered in the unrotated limit, i.e., ˜ O (0) = O .Applying the above procedure to the complete operator O provides the normal-ordered form of the transformedoperator ˜ O ( ϕ ) ≡ ˜ O [0] ( ϕ ) + ˜ O [2] ( ϕ ) + ˜ O [4] ( ϕ ) + ˜ O [6] ( ϕ ) , (B.6)in which the term ˜ O nm ( ϕ ) collects various contributions ˜ O nm ( ij ) ( ϕ ) . Each term ˜ O nm ( ϕ ) possesses the same op-erator structure as the corresponding term in Eq. (11),except that the original matrix elements are replaced bygauge-dependent ones, e.g., O k k k k is formally replacedby ˜ O k k k k ( ϕ ) . The expressions of the matrix elementsof each normal-ordered contribution ˜ O nm ( ϕ ) in terms ofthe matrix elements of the original normal-ordered contri-butions to an operator O with deg_max = 4 can be foundin Ref. [9]. Appendix C. Changelog
Since the previous main version of
ADG , the main mod-ifications to the software have focused on: • Adding the particle-number projected BMBPT for-malism. • Fixing an error in a symmetry factor in BMBPT dia-grams. • Porting the code to Python3 while maintaining com-patibility with Python2.7. • Changing the matrices data structures to
NumPy ar-rays. • Removing deprecated calls to
NetworkX . • Various optimisations to the BMBPT diagram gener-ation process. • Reducing the memory requirements of diagrams. • Fixing several errors in the documentation. • Fixing the installation process for additional depen-dencies.
References [1] P. Arthuis, T. Duguet, A. Tichai,R.-D. Lasseri, and J.-P. Ebran,Comp. Phys. Commun , 202 (2019).[2] K. Tsukiyama, S. K. Bogner, and A. Schwenk,Phys. Rev. Lett. , 222502 (2011).[3] V. Somà, T. Duguet, and C. Barbieri,Phys. Rev. C , 064317 (2011).[4] H. Hergert, S. K. Bogner, S. Binder, A. Calci,J. Langhammer, R. Roth, and A. Schwenk,Phys. Rev. C , 034307 (2013).[5] H. Hergert, S. K. Bogner, T. D. Morris, S. Binder,A. Calci, J. Langhammer, and R. Roth,Phys. Rev. C , 041302 (2014).[6] T. Duguet, J. Phys. G , 025107 (2015).[7] S. K. Bogner, H. Hergert, J. D. Holt, A. Schwenk,S. Binder, A. Calci, J. Langhammer, and R. Roth,Phys. Rev. Lett. , 142501 (2014).[8] G. R. Jansen, J. Engel, G. Hagen, P. Navratil, andA. Signoracci, Phys. Rev. Lett. , 142502 (2014).[9] T. Duguet and A. Signoracci,J. Phys. G , 015103 (2017).[10] E. Gebrerufael, K. Vobig, H. Hergert, and R. Roth,Phys. Rev. Lett. , 152503 (2017).[11] A. Tichai, E. Gebrerufael, K. Vobig, and R. Roth,Phys. Lett. B , 448 (2018).[12] A. Tichai, P. Arthuis, T. Duguet,H. Hergert, V. Somà, and R. Roth,Phys. Lett. B , 195 (2018).[13] A. Tichai, R. Roth, and T. Duguet,Front. Phys. , 164 (2020).[14] P. Demol, M. Frosini, A. Tichai, V. Somà, andT. Duguet, (2020), arXiv:2002.02724 [nucl-th] .[15] C. Yannouleas and U. Landman,Rept. Prog. Phys. , 2067 (2007).[16] T. Papenbrock and H. A. Weidenmueller,Phys. Rev. C , 014334 (2014).[17] P. Ring and P. Schuck, The Nuclear Many-Body Prob-lem (Springer-Verlag, New-York, 1980).[18] A. Signoracci, T. Duguet, G. Hagen, and G. R.Jansen, Phys. Rev. C , 064320 (2015).[19] J. Ripoche, A. Tichai, andT. Duguet, Eur. Phys. J. A56 , 40 (2020),arXiv:1908.00765 [nucl-th] .[20] G. Hagen, T. Papenbrock, D. J. Dean,A. Schwenk, A. Nogga, M. Włoch, and P. Piecuch,Phys. Rev. C , 034302 (2007).[21] R. Roth, S. Binder, K. Vobig, A. Calci,J. Langhammer, and P. Navratil,Phys. Rev. Lett. , 052501 (2012).[22] S. Binder, J. Langhammer, A. Calci, P. Navratil, and23. Roth, Phys. Rev. C , 021303 (2013).[23] S. Binder, P. Piecuch, A. Calci, J. Lang-hammer, P. Navrátil, and R. Roth,Phys. Rev. C , 054319 (2013).[24] J. Blaizot and G. Ripka, Quantum Theory of Fi-nite Systems (MIT Press, Cambridge, Massachusetts,1986).[25] R. Balian and E. Brézin,Nuovo Cim. B , 37 (1969).[26] A. A. Hagberg, D. A. Schult, and P. J. Swart, in Proceedings of the 7th Python in Science Conference ,edited by G. Varoquaux, T. Vaught, and J. Millman(Pasadena, CA USA, 2008) pp. 11 – 15.[27] T. Ohl, Comput. Phys. Commun. , 340 (1995).[28] I. Shavitt and R. J. Bartlett, Many-Body Methods in Chemistry and Physics: MBPT and Coupled-Cluster Theory ,Cambridge Molecular Science (Cambridge UniversityPress, 2009).[29] D. J. Thouless, Annals of Physics10