Challenges
Datasets
Workspaces
Papers
Discussions
Leaderboard
Log inSign up
Challenges
Datasets
Workspaces
Papers
Discussions
Leaderboard
Blog
Job Board
Q3AS

© 2026 Aqora Quantum S.A.S.

TermsPrivacyLegal Notice

Development and application of symbolic diagrammatic tools for multi-orbital correlated systems

Rayan Farid

Published May 1, 2026• 0 citations
OpenAlex

Abstract

This thesis develops and applies numerical tools for studying strongly correlated models at finite temperature. Strongly correlated materials exhibit rich emergent behavior, including high temperature superconductivity, exotic magnetism, and metal–insulator transitions. To study these phenomena, Feynman diagrammatics based on perturbation theory provides the foundation for both perturbative and non-perturbative numerical techniques. The thesis focuses on the diagrammatic evaluation of singleand two-particle response functions, with the truncation order acting as a control parameter. Recent advances in symbolic computation have led to the Algorithmic Matsubara Integration (AMI) technique, a central computational technique utilized in this work. AMI generates a symbolic representation of arbitrary Feynman diagrams and resolves the nested Matsubara sums via the residue theorem, yielding closedform semi-analytical expressions, with the remaining momentum integrals evaluated stochastically. The first part of the thesis adopts a pedagogical approach, introducing the generic many-body Hamiltonian, motivating the Hubbard model as a prototype for correlated systems, and reviewing Green’s function formalism in many-body theory. In the first application, we employ AMI to evaluate the pairing susceptibility to measure superconducting fluctuations in the single-band Hubbard model. We identify two dominant 𝑞 = [0,0] and 𝑞 = [π,π] 𝑑-wave pairing modes and map the correlation length and amplitude across the phase space 𝑇–𝑡'– 𝜇, defined by temperature 𝑇, nearest-neighbor hopping 𝑡' and chemical potential 𝜇. Next, we present an algorithm that assigns band indices to arbitrary Feynman topologies symbolically and iteratively using the non-zero elements of interaction tensor Uabcd. Applied to all connected topologies, it generates the full set of contributing diagrams for multiorbital problems. Combined with AMI, this framework termed MB-AMI enables systematic diagrammatic expansion for any low-energy multi-band system with arbitrary longrange interactions. We benchmark the method on molecular H₂ and on the Hubbard ii model under a Zeeman field. Using MB-AMI, we study the multilayer Hubbard model and the layer-resolved features of the pairing susceptibility. We report a monotonic but sublinear increase in the pairing response with layer stacking, assisted by pair-hopping interactions, with a well defined three-dimensional limit. Finally, motivated by the spin-dependent renormalization of screened interactions 𝑊ᵤᵤ and 𝑊ud, we propose a phenomenological extended Hubbard model where nearest-neighbor interactions 𝑉 are spin-dependent and independently tunable. We explore all four quadrants of the 𝑉ᵤᵤ–𝑉ud phase space, focusing particularly on the ferromagnetic-like regime 𝑉ud > 𝑉ᵤᵤ. In this regime, we find that competing ferromagnetic and antiferromagnetic spin fluctuations give rise to an exotic spin-density wave (SDW). We study the resulting fermi surface reconstruction and the emergence of triplet 𝑝-wave pairing from the induced SDW. In conclusion, we propose that merging MB-AMI with the Discrete Lehmann Representation (DLR) into a unified framework opens a pathway to performing high order, fully self-consistent diagrammatic calculations for accurate ab initio simulations of realistic materials.

Related papers

Renormalized diagrammatic perturbation theory using algorithmic Matsubara integration

Daria Gazizova

This thesis develops and applies numerical methods for computing correlated electrons at finite temperatures within perturbation theory, formulated in terms of Feynman diagrams. First we compute the fully renormalized effective interaction W in the 2D Hubbard model beyond the conventional low-order random phase approximation (RPA) and additionally resolve W into same-spin and opposite-spin channels. This analysis shows that starting from a purely local repulsive bare interaction U we get attractive effective nearest-neighbor interaction. We identify regimes of inverse temperature t and interaction U=t, where t denotes the nearest-neighbor hopping amplitude, in which this attraction appears. Such an emergent of attractive interaction is known to be able to promote superconductivity in the repulsive Hubbard model. Next, we introduce a method for incorporating higher-order diagrams within selfconsistent perturbative schemes. We merge algorithmic Matsubara integration (AMI) with discrete pole representations, including the discrete Lehmann representation and Prony (minimal pole) representation. AMI automates nested internal Matsubara frequency summations for arbitrary diagrams topologies through symbolic contour integration, which requires explicit access to the analytic pole structure of the input Green’s functions. By compressing the Green’s function into a small set of poles and corresponding weights, we regain access to its analytic pole structure, which allows us to iterate the calculation within self-consistent schemes. Benchmarks on the Bethe lattice demonstrate accurate single-shot and self-consistent self-energies on both Matsubara and real axes. For self-consistent calculations, we use two standard schemes, where AMI combined with the discrete pole representation works as a selfenergy and impurity solver. We further extend calculations to third and fourth order of interaction in self-energy diagrams. This approach yields the ability to perform self-consistent calculations at a fraction of the cost compared to other Matsubara-axis methods without moving to the real axis at every step of the iteration. Finally we demonstrate the universality of the AMI+DLR framework, by decomposing any diagram evaluation into a computational kernel that depends only on the diagram topology and DLR parameters, and a set of coefficients determined by the physical system. This separation enables precomputation and reuse of the computational kernel for a certain diagram topology and any physical problem, both in single-shot and iterative calculations. We illustrate the universality by reusing the same computational kernel for two simple systems: the 2D Hubbard square lattice and the H₂ molecule.

•May 1, 2026•0 citations

Diagrammatic Methods for Strongly-Correlated Electrons in Realistic Systems

Jia Li

This thesis is devoted to the development of numerical methods for correlated electrons in realistic systems based on perturbation theories using Feynman diagrams. Electron correlations in realistic systems give rise to a wide range of interesting physical properties. Diagrammatic expansions and its low-order approximations are key building blocks of theoretical and numerical methods for studying correlated electrons. Recent developments of the diagrammatic quantum Monte Carlo (DiagMC) methods enable efficient stochastic evaluations of diagrams and does not suffer from numerical sign problem which worsens with increasing system size. It is therefore a promising direction to apply DiagMC techniques to realistic electron systems in the pursuit of ab initio calculations of correlated materials. This thesis starts by reviewing basic concepts and definitions of realistic electron systems, the perturbation theory based on Feynman diagrams, Monte Carlo evaluations of the series, as well as self-consistent diagrammatic methods. Low-order diagrammatic methods are applied to realistic calculations and benchmarked against large collections of state-of-the-art many-body methods. The thesis then proposes a diagrammatic Monte Carlo solver for realistic molecules and impurities based on a recursive evaluation of diagrams, which is tested on realistic molecules and impurities with varying temperature, system size, and interaction strength. We observe that this method is ideal for problems with many orbitals with moderate correlations, but may suffer from series divergence at strong correlations. Next, we propose the inchworm Monte Carlo algorithm for interaction expansion, which overcomes limitations of divergent diagrammatic series by replacing the standard perturbation series with a sequence of incremental series expansions. A prototype implementation is presented for system with on-site Hubbard interactions. Finally, we discuss technical aspects of numerical representations of the Green's function, which has significant impact on the efficiency of diagrammatic methods. We review several common numerical representations, and introduce our development of the sparse sampling method, which provides compact representations of the Green's function in both imaginary time and frequency as well as fast transformation formula between the two.

•Jan 1, 2021•0 citations

Algorithmic approach to diagrammatic expansions for real-frequency evaluation of susceptibility functions

Amir Taheridehkordi, Stephanie H. Curnoe, J. P. F. LeBlanc

We systematically generate the perturbative expansion for the two-particle spin susceptibility in the Feynman diagrammatic formalism and apply this expansion to a model system---the single-band Hubbard model on a square lattice. We make use of algorithmic Matsubara integration (AMI) [A. Taheridehkordi, S. H. Curnoe, and J. P. F. LeBlanc, Phys. Rev. B 99, 035120 (2019)] to analytically evaluate Matsubara frequency summations, allowing us to symbolically impose analytic continuation to the real-frequency axis. We minimize our computational expense by applying graph invariant transformations [A. Taheridehkordi, S. H. Curnoe, and J. P. F. LeBlanc, Phys. Rev. B 101, 125109 (2020)]. We highlight extensions of the random-phase approximation and T-matrix methods that, due to AMI, become tractable. We present results for weak interaction strength where the direct perturbative expansion is convergent, and verify our results on the Matsubara axis by comparison to other numerical methods. By examining the spin susceptibility as a function of real frequency via an order-by-order expansion, we can identify precisely what role higher-order corrections play on spin susceptibility and demonstrate the utility and limitations of our approach.

arXiv:2004.11091•Physical review. B./Physical review. B•Jul 9, 2020•28 citations

Dynamical Mean-Field Theory Studies of Correlated Electrons

Yang Yu

Correlated electron systems exhibit a variety of intriguing emergent phenomena, including Mott physics, pseudogap behavior, and unconventional superconductivity. One of the most widely used numerical approaches for investigating these systems is dynamical mean-field theory (DMFT). This dissertation addresses two related problems: how to develop better quantum impurity solvers for DMFT and how to use DMFT results to gain a better understanding of correlated electron systems. The first part of this dissertation develops numerical methods for solving quantum impurity problems, which constitute the central computational task in DMFT calculations. In Chap. 4, we develop a method that combines a tensor-train numerical integration scheme with the inchworm hybridization expansion for quantum impurity problems. Benchmarks on spinless, single-orbital, and multi-orbital impurity models demonstrate the convergence of the Green’s functions to exact or reference results and identify the regimes in which tensor-train compression of Feynman diagrams and the inchworm expansion itself are effective. In Chap. 5, we introduce the concept of counterterms within the hybridization expansion and apply it to a bare hybridization expansion solver that relies on a sum-of-exponentials integration technique. Tests on spinless, multi-orbital, and realistic embedding impurity problems show that counterterms can significantly reduce the required perturbative order and improve the efficiency of hybridization-expansion calculations. The second part of this dissertation applies the cluster extension of DMFT to the single-band Hubbard model. In Chap. 6, we introduce a fluctuation decomposition for the self-energy that eliminates the Fierz ambiguity encountered in conventional fluctuation diagnostics of the Hubbard model. The resulting analysis shows that spin fluctuations make an important contribution to the pseudogap self-energy, while multi-boson contributions are also necessary for a complete description. In Chap. 7, we analyze the spin-fermion coupling and its role in d-wave pairing. Our calculations show that, for antinodal fermions in the pseudogap regime, their coupling to antiferromagnetic spin fluctuations is strongly enhanced relative to the non-interacting limit. A modified Bethe–Salpeter equation analysis further shows that this enhanced spin–fermion coupling boosts the spin-fluctuation-induced contribution to the d-wave pairing susceptibility.

•Jan 1, 2026•0 citations

Diagrammatics in the Dual Space, or There and Back Again

E. A. Stepanov

Accurately describing many-body effects in multi-orbital systems remains a major challenge in theoretical condensed matter physics. At present, there is a significant methodological gap between the numerical tools used in ab initio computational materials science and those developed to study strong electronic correlations. The former can treat realistic, large-scale systems but typically neglect many-body effects, while the latter focus on simplified models with only a few degrees of freedom, as only such models can be solved accurately in the presence of strong interactions. The purpose of this thesis is to bridge these two approaches and establish a systematic theoretical framework for realistic correlated electronic materials. This involves a full-cycle methodology that begins with constructing ab initio interacting models from density-functional theory, solving them using dynamical mean-field theory to capture local correlations, and extending beyond to incorporate non-local collective electronic fluctuations. To this end, we introduce the "dual" approach to strong correlations, which includes the dual fermion, dual boson, and dual triply irreducible local expansion methods. The central idea of the dual theories is to shift the reference point of the conventional Feynman diagrammatic expansion from a non-interacting electronic system to an interacting but exactly solvable one. Integrating out this reference system recasts the expansion in an effective dual space, where all diagrammatic building blocks are renormalized by the corresponding impurity quantities. This procedure transforms a non-perturbative expansion in the original variables into a perturbative expansion in terms of dual fermionic and bosonic fields, exact in both the weak- and strong-coupling limits. In this thesis, we collect and systematize the major developments of dual techniques achieved to date.

arXiv:2512.12389•Dec 13, 2025•0 citations

Real-frequency diagrammatic Monte Carlo at finite temperature

J. Vučičević, Michel Ferrero

Diagrammatic expansions are a central tool for treating correlated electron systems. At thermal equilibrium, they are most naturally defined within the Matsubara formalism. However, extracting any dynamic response function from a Matsubara calculation ultimately requires the ill-defined analytical continuation from the imaginary- to the real-frequency domain. It was recently proposed [A. Taheridehkordi et al., Phys. Rev. B 99, 035120 (2019)] that the internal Matsubara summations of any interaction-expansion diagram can be performed analytically by using symbolic algebra algorithms. The result of the summations is then an analytical function of the complex frequency rather than Matsubara frequency. Here we apply this principle and develop a diagrammatic Monte Carlo technique which yields results directly on the real-frequency axis. We present results for the self-energy $\mathrm{\ensuremath{\Sigma}}(\ensuremath{\omega})$ of the doped $32\ifmmode\times\else\texttimes\fi{}32$ cyclic square-lattice Hubbard model in a nontrivial parameter regime, where signatures of the pseudogap appear close to the antinode. We discuss the behavior of the perturbation series on the real-frequency axis and in particular show that one must be very careful when using the maximum entropy method on truncated perturbation series. Our approach holds great promise for future application in cases when analytical continuation is difficult and moderate-order perturbation theory may be sufficient to converge the result.

arXiv:1908.11826•Physical review. B./Physical review. B•Feb 12, 2020•35 citations
This thesis develops and applies numerical tools for studying strongly correlated models at finite temperature. Strongly correlated materials exhibit rich emergent behavior, including high temperature superconductivity, exotic magnetism, and metal–insulator transitions. To study these phenomena, Feynman diagrammatics based on perturbation theory provides the foundation for both perturbative and non-perturbative numerical techniques. The thesis focuses on the diagrammatic evaluation of singleand two-particle response functions, with the truncation order acting as a control parameter. Recent advances in symbolic computation have led to the Algorithmic Matsubara Integration (AMI) technique, a central computational technique utilized in this work. AMI generates a symbolic representation of arbitrary Feynman diagrams and resolves the nested Matsubara sums via the residue theorem, yielding closedform semi-analytical expressions, with the remaining momentum integrals evaluated stochastically. The first part of the thesis adopts a pedagogical approach, introducing the generic many-body Hamiltonian, motivating the Hubbard model as a prototype for correlated systems, and reviewing Green’s function formalism in many-body theory. In the first application, we employ AMI to evaluate the pairing susceptibility to measure superconducting fluctuations in the single-band Hubbard model. We identify two dominant 𝑞 = [0,0] and 𝑞 = [π,π] 𝑑-wave pairing modes and map the correlation length and amplitude across the phase space 𝑇–𝑡'– 𝜇, defined by temperature 𝑇, nearest-neighbor hopping 𝑡' and chemical potential 𝜇. Next, we present an algorithm that assigns band indices to arbitrary Feynman topologies symbolically and iteratively using the non-zero elements of interaction tensor Uabcd. Applied to all connected topologies, it generates the full set of contributing diagrams for multiorbital problems. Combined with AMI, this framework termed MB-AMI enables systematic diagrammatic expansion for any low-energy multi-band system with arbitrary longrange interactions. We benchmark the method on molecular H₂ and on the Hubbard ii model under a Zeeman field. Using MB-AMI, we study the multilayer Hubbard model and the layer-resolved features of the pairing susceptibility. We report a monotonic but sublinear increase in the pairing response with layer stacking, assisted by pair-hopping interactions, with a well defined three-dimensional limit. Finally, motivated by the spin-dependent renormalization of screened interactions 𝑊ᵤᵤ and 𝑊ud, we propose a phenomenological extended Hubbard model where nearest-neighbor interactions 𝑉 are spin-dependent and independently tunable. We explore all four quadrants of the 𝑉ᵤᵤ–𝑉ud phase space, focusing particularly on the ferromagnetic-like regime 𝑉ud > 𝑉ᵤᵤ. In this regime, we find that competing ferromagnetic and antiferromagnetic spin fluctuations give rise to an exotic spin-density wave (SDW). We study the resulting fermi surface reconstruction and the emergence of triplet 𝑝-wave pairing from the induced SDW. In conclusion, we propose that merging MB-AMI with the Discrete Lehmann Representation (DLR) into a unified framework opens a pathway to performing high order, fully self-consistent diagrammatic calculations for accurate ab initio simulations of realistic materials.
This thesis develops and applies numerical methods for computing correlated electrons at finite temperatures within perturbation theory, formulated in terms of Feynman diagrams. First we compute the fully renormalized effective interaction W in the 2D Hubbard model beyond the conventional low-order random phase approximation (RPA) and additionally resolve W into same-spin and opposite-spin channels. This analysis shows that starting from a purely local repulsive bare interaction U we get attractive effective nearest-neighbor interaction. We identify regimes of inverse temperature t and interaction U=t, where t denotes the nearest-neighbor hopping amplitude, in which this attraction appears. Such an emergent of attractive interaction is known to be able to promote superconductivity in the repulsive Hubbard model. Next, we introduce a method for incorporating higher-order diagrams within selfconsistent perturbative schemes. We merge algorithmic Matsubara integration (AMI) with discrete pole representations, including the discrete Lehmann representation and Prony (minimal pole) representation. AMI automates nested internal Matsubara frequency summations for arbitrary diagrams topologies through symbolic contour integration, which requires explicit access to the analytic pole structure of the input Green’s functions. By compressing the Green’s function into a small set of poles and corresponding weights, we regain access to its analytic pole structure, which allows us to iterate the calculation within self-consistent schemes. Benchmarks on the Bethe lattice demonstrate accurate single-shot and self-consistent self-energies on both Matsubara and real axes. For self-consistent calculations, we use two standard schemes, where AMI combined with the discrete pole representation works as a selfenergy and impurity solver. We further extend calculations to third and fourth order of interaction in self-energy diagrams. This approach yields the ability to perform self-consistent calculations at a fraction of the cost compared to other Matsubara-axis methods without moving to the real axis at every step of the iteration. Finally we demonstrate the universality of the AMI+DLR framework, by decomposing any diagram evaluation into a computational kernel that depends only on the diagram topology and DLR parameters, and a set of coefficients determined by the physical system. This separation enables precomputation and reuse of the computational kernel for a certain diagram topology and any physical problem, both in single-shot and iterative calculations. We illustrate the universality by reusing the same computational kernel for two simple systems: the 2D Hubbard square lattice and the H₂ molecule.
This thesis is devoted to the development of numerical methods for correlated electrons in realistic systems based on perturbation theories using Feynman diagrams. Electron correlations in realistic systems give rise to a wide range of interesting physical properties. Diagrammatic expansions and its low-order approximations are key building blocks of theoretical and numerical methods for studying correlated electrons. Recent developments of the diagrammatic quantum Monte Carlo (DiagMC) methods enable efficient stochastic evaluations of diagrams and does not suffer from numerical sign problem which worsens with increasing system size. It is therefore a promising direction to apply DiagMC techniques to realistic electron systems in the pursuit of ab initio calculations of correlated materials. This thesis starts by reviewing basic concepts and definitions of realistic electron systems, the perturbation theory based on Feynman diagrams, Monte Carlo evaluations of the series, as well as self-consistent diagrammatic methods. Low-order diagrammatic methods are applied to realistic calculations and benchmarked against large collections of state-of-the-art many-body methods. The thesis then proposes a diagrammatic Monte Carlo solver for realistic molecules and impurities based on a recursive evaluation of diagrams, which is tested on realistic molecules and impurities with varying temperature, system size, and interaction strength. We observe that this method is ideal for problems with many orbitals with moderate correlations, but may suffer from series divergence at strong correlations. Next, we propose the inchworm Monte Carlo algorithm for interaction expansion, which overcomes limitations of divergent diagrammatic series by replacing the standard perturbation series with a sequence of incremental series expansions. A prototype implementation is presented for system with on-site Hubbard interactions. Finally, we discuss technical aspects of numerical representations of the Green's function, which has significant impact on the efficiency of diagrammatic methods. We review several common numerical representations, and introduce our development of the sparse sampling method, which provides compact representations of the Green's function in both imaginary time and frequency as well as fast transformation formula between the two.
We systematically generate the perturbative expansion for the two-particle spin susceptibility in the Feynman diagrammatic formalism and apply this expansion to a model system---the single-band Hubbard model on a square lattice. We make use of algorithmic Matsubara integration (AMI) [A. Taheridehkordi, S. H. Curnoe, and J. P. F. LeBlanc, Phys. Rev. B 99, 035120 (2019)] to analytically evaluate Matsubara frequency summations, allowing us to symbolically impose analytic continuation to the real-frequency axis. We minimize our computational expense by applying graph invariant transformations [A. Taheridehkordi, S. H. Curnoe, and J. P. F. LeBlanc, Phys. Rev. B 101, 125109 (2020)]. We highlight extensions of the random-phase approximation and T-matrix methods that, due to AMI, become tractable. We present results for weak interaction strength where the direct perturbative expansion is convergent, and verify our results on the Matsubara axis by comparison to other numerical methods. By examining the spin susceptibility as a function of real frequency via an order-by-order expansion, we can identify precisely what role higher-order corrections play on spin susceptibility and demonstrate the utility and limitations of our approach.
Correlated electron systems exhibit a variety of intriguing emergent phenomena, including Mott physics, pseudogap behavior, and unconventional superconductivity. One of the most widely used numerical approaches for investigating these systems is dynamical mean-field theory (DMFT). This dissertation addresses two related problems: how to develop better quantum impurity solvers for DMFT and how to use DMFT results to gain a better understanding of correlated electron systems. The first part of this dissertation develops numerical methods for solving quantum impurity problems, which constitute the central computational task in DMFT calculations. In Chap. 4, we develop a method that combines a tensor-train numerical integration scheme with the inchworm hybridization expansion for quantum impurity problems. Benchmarks on spinless, single-orbital, and multi-orbital impurity models demonstrate the convergence of the Green’s functions to exact or reference results and identify the regimes in which tensor-train compression of Feynman diagrams and the inchworm expansion itself are effective. In Chap. 5, we introduce the concept of counterterms within the hybridization expansion and apply it to a bare hybridization expansion solver that relies on a sum-of-exponentials integration technique. Tests on spinless, multi-orbital, and realistic embedding impurity problems show that counterterms can significantly reduce the required perturbative order and improve the efficiency of hybridization-expansion calculations. The second part of this dissertation applies the cluster extension of DMFT to the single-band Hubbard model. In Chap. 6, we introduce a fluctuation decomposition for the self-energy that eliminates the Fierz ambiguity encountered in conventional fluctuation diagnostics of the Hubbard model. The resulting analysis shows that spin fluctuations make an important contribution to the pseudogap self-energy, while multi-boson contributions are also necessary for a complete description. In Chap. 7, we analyze the spin-fermion coupling and its role in d-wave pairing. Our calculations show that, for antinodal fermions in the pseudogap regime, their coupling to antiferromagnetic spin fluctuations is strongly enhanced relative to the non-interacting limit. A modified Bethe–Salpeter equation analysis further shows that this enhanced spin–fermion coupling boosts the spin-fluctuation-induced contribution to the d-wave pairing susceptibility.
Accurately describing many-body effects in multi-orbital systems remains a major challenge in theoretical condensed matter physics. At present, there is a significant methodological gap between the numerical tools used in ab initio computational materials science and those developed to study strong electronic correlations. The former can treat realistic, large-scale systems but typically neglect many-body effects, while the latter focus on simplified models with only a few degrees of freedom, as only such models can be solved accurately in the presence of strong interactions. The purpose of this thesis is to bridge these two approaches and establish a systematic theoretical framework for realistic correlated electronic materials. This involves a full-cycle methodology that begins with constructing ab initio interacting models from density-functional theory, solving them using dynamical mean-field theory to capture local correlations, and extending beyond to incorporate non-local collective electronic fluctuations. To this end, we introduce the "dual" approach to strong correlations, which includes the dual fermion, dual boson, and dual triply irreducible local expansion methods. The central idea of the dual theories is to shift the reference point of the conventional Feynman diagrammatic expansion from a non-interacting electronic system to an interacting but exactly solvable one. Integrating out this reference system recasts the expansion in an effective dual space, where all diagrammatic building blocks are renormalized by the corresponding impurity quantities. This procedure transforms a non-perturbative expansion in the original variables into a perturbative expansion in terms of dual fermionic and bosonic fields, exact in both the weak- and strong-coupling limits. In this thesis, we collect and systematize the major developments of dual techniques achieved to date.
Diagrammatic expansions are a central tool for treating correlated electron systems. At thermal equilibrium, they are most naturally defined within the Matsubara formalism. However, extracting any dynamic response function from a Matsubara calculation ultimately requires the ill-defined analytical continuation from the imaginary- to the real-frequency domain. It was recently proposed [A. Taheridehkordi et al., Phys. Rev. B 99, 035120 (2019)] that the internal Matsubara summations of any interaction-expansion diagram can be performed analytically by using symbolic algebra algorithms. The result of the summations is then an analytical function of the complex frequency rather than Matsubara frequency. Here we apply this principle and develop a diagrammatic Monte Carlo technique which yields results directly on the real-frequency axis. We present results for the self-energy $\mathrm{\ensuremath{\Sigma}}(\ensuremath{\omega})$ of the doped $32\ifmmode\times\else\texttimes\fi{}32$ cyclic square-lattice Hubbard model in a nontrivial parameter regime, where signatures of the pseudogap appear close to the antinode. We discuss the behavior of the perturbation series on the real-frequency axis and in particular show that one must be very careful when using the maximum entropy method on truncated perturbation series. Our approach holds great promise for future application in cases when analytical continuation is difficult and moderate-order perturbation theory may be sufficient to converge the result.