DRAFT: Deterministic Multi-sensor Measurement-adaptive Birth using Labeled Random Finite Sets

Compiled on August 27, 2024

Git Hash:

Git Branch:

Deterministic Multi-sensor Measurement-adaptive Birth using Labeled Random Finite Sets

Jennifer Bondarchuk,
Anthony Trezza,
Donald J. Bucci Jr.,

Advanced Technology Laboratories

Lockheed Martin

Cherry Hill, NJ, USA

{jennifer.a.bondarchuk, anthony.t.trezza<EMAIL_ADDRESS>
Measurement-adaptive track initiation remains a critical design requirement of many practical multi-target tracking systems.
For labeled random finite sets multi-object filters, prior work has been established to construct a labeled multi-object birth density using measurements from multiple sensors.
A truncation procedure has also been provided that leverages a stochastic Gibbs sampler to truncate the birth density for scalability.
In this work, we introduce a deterministic herded Gibbs sampling truncation solution for efficient multi-sensor adaptive track initialization.
Removing the stochastic behavior of the track initialization procedure without impacting average tracking performance enables a more robust tracking solution more suitable for safety-critical applications.
Simulation results for linear sensing scenarios are provided to verify performance.

Random finite sets,
Target tracking,
Herded Gibbs sampling,
Measurement adaptive birth

§ INTRODUCTION

MTT is a sequential estimation problem where the goal is to jointly determine the number of objects and their trajectories from noisy, discrete-time measurements.
It is distinct from state estimation filtering in the sense that object appearance and disappearance (i.e., birth and death) must be handled alongside object detection and miss-detection in superimposed false measurements (i.e., clutter).
Approaches such as GNN techniques <cit.>, JPDA <cit.>, and MHT <cit.> traditionally handle object appearance and disappearance as an implementation detail.
Labeled RFS <cit.> has recently emerged as a family of techniques that include track birth and death in the problem formulation which is then solved via unified Bayesian framework.
We direct the reader to <cit.> for a detailed survey of the field, recent advances, and example applications.

The concept of conjugacy is of critical importance in labeled RFS for maintaining computationally tractable MTT <cit.>.
In the seminal works of <cit.>, it was shown that the GLMB RFS forms a conjugate multi-object prior with itself under the labeled RFS Bayesian filtering equations as introduced in FISST <cit.>.
Under this construction the prior RFS density for appearing targets is equivalently a labeled Poisson or a LMB RFS.
In a static birth approach the prior RFS density is assumed to be time invariant and established a priori.
In contrast, an adaptive birth approach does not assume a time invariant birth distribution and must be constructed adaptively from the measurements <cit.>.

Measurement-adaptive birth techniques in RFS began with the SMC PHD and CPHD filters in <cit.>, and later extended to their GM variants in <cit.>.
Single-sensor adaptive birth techniques for labeled RFS were introduced in <cit.> for the LMB filter and <cit.> for the $\delta$-GLMB filter.
These approaches follow the concepts presented in <cit.> by constructing a component of an LMB RFS from every single-sensor measurement observed at the previous time step.
The probabilities of birth for each of these LMB components were set proportionally to a maximum birth rate parameter.
The proportionality factor was based on how much a measurement did not associate with the persisting targets from the last time step.
In contrast, the multi-sensor adaptive birth regime carries the additional complication that there are exponentially many possible measurement tuples that could be used to generate newborn targets.
In our previous work <cit.>, we proposed a stochastic Gibbs sampling approach for truncating multi-sensor measurement tuples that would likely be pruned in the subsequent $\delta$-GLMB update step.
The proposed Gibbs sampler exhibited quadratic complexity in the number of sensors, and both an SMC approximation and a closed form Gaussian solution were provided.

The drawback of using a stochastic Gibbs sampler for multi-sensor adaptive birth is that it adds non-determinism during track initiation.
Stochasticity during track initiation can cause substantial differences in the resulting track output, which can negatively impact tracking performance.
This run-to-run variability hinders verification and validation testing in production systems, especially when safety-guaranteed operation is critical <cit.>.

In this work, we draw inspiration from <cit.> to construct a deterministic Gibbs sampler for the multi-sensor adaptive birth technique presented in <cit.>.
The proposed deterministic, multi-sensor adaptive birth procedure leverages a deterministic sampling technique known as herding <cit.> as applied to Gibbs sampling <cit.>.
As we will show, this deterministic multi-sensor adaptive birth procedure converges faster, requiring fewer samples, and maintains the average tracking performance presented in <cit.>.
It has the additional benefit that for the closed-form Gaussian variation, the proposed approach is fully deterministic.

The remainder of this paper is organized as follows.
Section <ref> provides background material on the relevant labeled RFS terminology and a brief review of the multi-sensor adaptive birth technique of <cit.>.
Section <ref> provides a brief background on deterministic sampling via herding and introduces the deterministic multi-sensor adaptive birth procedure.
Finally Section <ref> compares the proposed deterministic adaptive birth technique against the stochastic version of <cit.> within the LMB and $\delta$-GLMB filters.

§ BACKGROUND

§.§ Notation

We adopt the notation of <cit.> as follows.
Lowercase letters denote vectors (e.g., $x, \textbf{x}$), whereas uppercase letters denote finite sets (e.g., $X, \textbf{X}$).
Bold letters denote labeled states and their distributions (e.g., $\textbf{x}, \textbf{X}, \boldsymbol{\pi}$) and blackboard bold letters denote spaces (e.g., $\mathbb{X}, \mathbb{Z}, \mathbb{L}, \mathbb{R}$).
Variable sequences $X_i, X_{i+1}, \dots, X_j$ are abbreviated $X_{i:j}$.
The inner product $\int f(x) g(x) dx$ is written as $\langle f, g \rangle$.
Finally, we drop the subscript notation for current time step $k$ and use subscript '$+$' to indicate the next time step ($k+1$).

§.§ Multi-sensor Multi-target Adaptive Birth Model

Our objective is to adaptively construct an LMB RFS density $\textbf{f}_{B,+}$ given the multi-object measurement sets $Z^{(s)} \in \mathbb{Z}^{(s)}$ from sensors $s \in \{1, \dots, V\}$.
Let $m^{(s)} = |Z^{(s)}|$ and $\mathbb{J}^{(s)} \triangleq \{1, \dots, m^{(s)}\}$ be an index set into $Z^{(s)}$.
Then let $\mathbb{J}^{(s)}_0 \triangleq \mathbb{J}^{(s)} \cup \{0\}$, where the entry $0$ represents a miss-detection.
The tuple $J \triangleq (j^{(1)}, \dots, j^{(V)}) \in \mathbb{J}^{(1)}_0 \times \dots \times \mathbb{J}^{(V)}_0$ is a multi-sensor measurement index with $Z_J \triangleq (z^{(1)}_{j^{(1)}}, \dots, z^{(V)}_{j^{(V)}})$ as the corresponding tuple of multi-sensor measurements.

In the single-sensor multi-target observation model of <cit.>, a target state $\textbf{x} \in \textbf{X}$ is either detected by sensor $s$ with probability $p_D^{(s)}(\textbf{x})$ or miss-detected with probability $1 - p_D^{(s)}(\textbf{x})$.
Target measurements $z^{(s)}$ are generated via the measurement likelihood function $g^{(s)}(z^{(s)}|\textbf{x})$.
The detection dynamics of a target are modeled as a multi-Bernoulli RFS with parameter set $\{(p_D^{(s)}(\textbf{x}), g^{(s)}(z^{(s)}|\textbf{x})): \textbf{x} \in \textbf{X}\}$, with conditional independence between each multi-Bernoulli RFS given $\textbf{X}$.
Each sensor generates clutter detections according to a Poisson RFS with intensity function $\kappa^{(s)}$.
The superposition of the target and clutter detections form the multi-object observation $Z^{(s)}$.
The multi-object measurement likelihood is the convolution of the target detection and clutter detection RFS distributions.
The exact form of this likelihood is omitted for brevity here, and can be found in <cit.>.
An important quantity within the multi-object likelihood relevant to this work is the per-object, per-sensor measurement pseudolikelihood

\begin{equation}\label{eq::meas_likelihood}
\psi^{s, j^{(s)}}_{Z^{(s)}}(\textbf{x}) =
\begin{cases}
\frac{p_D^{(s)}(\textbf{x})g^{(s)}(z^{(s)}_{j^{(s)}} | \textbf{x})}{\kappa^{(s)}(z^{(s)}_{j^{(s)}})} & j^{(s)} \in \mathbb{J}^{(s)} \\
1 - p_D^{(s)}(\textbf{x})                                                                            & j^{(s)} = 0
\end{cases}.
\end{equation}

As suggested in <cit.>, the birth LMB RFS density used at time ($k+1$) is modeled using the measurements at time $k$ as

\begin{equation}
\textbf{f}_{B,+} = \left\{ \left(r_{B,+}(l_+), p_{B,+}(\cdot, l_+ | Z_J) \right) \right\}_{l_+ \in \mathbb{B}_+},
\end{equation}

where $l_+ \in \mathbb{B}_+ = \{(k+1, J) : \forall J \in \mathbb{J}_0\}$ is the newborn object's label having birth probability $r_{B,+}(l_+)$.
The quantity $p_{B,+}(\cdot, l_+ | Z_J)$ is the state distribution of label $l_+$ conditioned on the measurement tuple $Z_J$.
The state distribution is given by the Chapman-Kolmogrov equation <cit.>

\begin{equation}\label{eq::pred_birth_post}
p_{B,+}(x_+, l_+ | Z_J) = \int f_+(x_+|x, l) p_B(x, l_+ | Z_J) dx,
\end{equation}

\begin{align}
p_B(x, l_+ | Z_J)       & = \frac{p_B(x, l_+) \psi^{J}_Z(x, l_+)}{\bar{\psi}^{J}_{Z}(l_+)} \label{eq::spatial_distr} \\
\bar{\psi}^{J}_{Z}(l_+) & = \langle p_B(\cdot, l_+), \psi^{J}_{Z}(\cdot, l_+)\rangle. \label{eq::psi_bar}
\end{align}

Here $f_+(x_+|x, l)$ is the Markov transition density modeling the dynamics of the object state and $p_B(x, l_+)$ is the prior distribution on the newborn object's state.
The quantity $\psi^{J}_{Z}(x, l_+)$ is the multi-sensor measurement pseudolikelihood function <cit.>.
Assuming conditional independence of the sensor measurements,

\begin{equation}\label{eq::joint_likelihood}
\psi^J_Z(\textbf{x}) \triangleq \prod\limits^V_{s=1} \psi^{s, j^{(s)}}_{Z^{(s)}}(\textbf{x}).
\end{equation}

The birth probability of label $l_+$ is modeled in <cit.> as

\begin{equation}\label{eq::birth_prob}
r_{B,+}(l_+) = \min\left(r_{B, \text{max}}, \hat{r}_{B,+}(l_+)\lambda_{B,+} \right),
\end{equation}

where $r_{B, \text{max}} \in [0, 1]$ is the maximum existence probability of a newborn target and $\lambda_{B,+}$ is the expected number of target births at time step $k+1$.
The effective birth probability $\hat{r}_{B,+}(l_+)$ is given as

\begin{equation}\label{eq::birth_prob_hat}
\hat{r}_{B,+}(l_+) = \frac{r_U(J) \bar{\psi}^{J}_{Z}(l_+)}{\sum\limits_{J' \in \mathbb{J}_0} r_U(J')\bar{\psi}^{J'}_{Z}(l_+)},
\end{equation}

where $r_U(J)$ is the probability that the multi-sensor measurement tuple $J$ does not associate with any existing targets.
This non-association probability is approximated as

\begin{equation}\label{eq::unassoc_prob}
r_{U}(J) \propto \prod_{j \in J} \left(1 - r_{A}(j)\right),
\end{equation}

where $r_{A}(j)$ is the probability that a single-sensor's measurement, $z^{(s)}_{j^{(s)}}$, is associated with existing targets in the $\delta$-GLMB posterior density.
For more information on these quantities, we direct the reader to <cit.>.

§.§ Truncating the Multi-sensor Adaptive Birth LMB

The number of newborn labels $\mathbb{B}_+$ is driven by the cardinality of $\mathbb{J}_0$ which has $O(m^V)$ elements, where $m$ is the worst case number of measurements.
It is therefore not practical to enumerate every possible entry at each time step.
In the spirit of <cit.>, a technique was proposed in <cit.> to sample labels $\mathbb{B}'_+ \subset \mathbb{B}_+$ that have high $\hat{r}_{B,+}(l_+) $ and hence are likely to persist through future measurement updates.
The sampling distribution is modeled as

\begin{equation}\label{eq::orig_sampling_distr}
p(l_+) \propto \hat{r}_{B,+}(l_+) \propto r_U(J)\bar{\psi}^{J}_{Z}(l_+).
\end{equation}

Directly sampling from the categorical distribution $p(l_+)$ requires evaluation of an exponential number of possible birth labels.
By construction of the birth set, sampling from $p(l_+)$ is equivalent to sampling multi-sensor measurement tuples from the joint distribution, $p(J) = p(j^{(1)}, \dots, j^{(V)})$.
A Gibbs sampler was thus proposed in <cit.> to achieve efficient sampling from Equation (<ref>) using the Markov transition probabilities

\begin{equation}\label{eq::cdn_likelihood}
p(j^{(s)} | J^{-s}) \propto \left(1 - r_A(j^{(s)})\right) \bar{\psi}^{J}_{Z}(l_+),
\end{equation}

where $J^{-s} = (j^{(1)}, \dots, j^{(s-1)}, j^{(s+1)}, \dots, j^{(V)})$.

The stochastic Gibbs sampler that generates the measurement tuples (i.e., birth components) with significant birth probabilities is given in <cit.>.
Similar to <cit.>, it exhibits an exponential convergence rate and does not require a burn-in period since unique solutions can be directly used in construction of the birth set.
An SMC method was provided in <cit.> for evaluating the transition probabilities and generating the state distributions.
Closed-form solutions for these quantities in the Gaussian case can be found in <cit.>.

§ HERDED GIBBS MULTI-SENSOR ADAPTIVE BIRTH

§.§ Deterministic Sampling via Herding

As discussed in <cit.>, herding is a deterministic procedure
to generate a sequence of pseudo-samples that match the moments, $\mu = \mathbb{E}_{x \sim P}\left\{ \phi(x) \right\}$, with respect to the features, $\phi(x)$, where $P$ is a discrete probability distribution.
It works using a collection of nonlinear update equations to solve a minmax problem.
The first phase selects a state such that the objective function is minimized.
The second phase updates auxiliary weights to maximize the objective function.
The $i^{\text{th}}$ sample $\gamma^{(i)} \in \mathbb{G}$ is computed as

\begin{align}
\gamma^{(i)} =  \argmax_{\gamma \in \mathbb{G}}  \langle w^{(i - 1)}, \phi(\gamma) \rangle \label{eq::herding_update1}\\
w^{(i)} = w^{(i-1)} + \mu - \phi(\gamma^{(i)}) \label{eq::herding_update2},
\end{align}

where $w \in \mathbb{W}$ is a weight vector that includes information on all previous samples, and $\phi : \mathbb{G} \rightarrow \mathbb{W}$ is a feature map.
If $\phi$ is an indicator function, it enables sampling from a probability distribution with $\mu = P$.
From here forward we will drop the superscript notation for the $i^{\text{th}}$ iteration for brevity and to avoid confusion with our superscript indexing notation.

In <cit.>, the herding procedure was applied to the discrete conditional distribution in a Gibbs sampler.
A set of all weight vectors $W$ is maintained containing the weight vector for all Markov transitions distributions.
For space efficiency, this can be maintained in a dictionary mapping $W: \mathbb{T} \times \mathbb{N} \rightarrow \mathbb{W}$,
such that $W(\gamma, n) \equiv w_{(\gamma, n)}$ is the weight vector for the transition distribution for the $n^{\text{th}}$ state of vector $\gamma$.
The dictionary $W$ is initialized as empty and $w_{(\gamma, n)}$ is initialized proportional to the transition density.
On every iteration, samples are generated by computing the argmax of $w_{(\gamma, n)}$ per Equation (<ref>) and $w_{(\gamma, n)}$ is updated according to Equation (<ref>).

§.§ Deterministic Multi-sensor Adaptive Birth

The herded Gibbs sampler for the multi-sensor adaptive birth truncation procedure is provided in Algorithm <ref>.
In practice, $p(j^{(s)} | J^{-s})$ can be computationally complex to compute for arbitrary transition densities since it requires the evaluation of multi-variate, indefinite integrals.
To alleviate the need to compute this repeatedly, we suggest caching the value in a dictionary mapping, $M : \mathbb{T} \times \mathbb{N} \rightarrow \mathbb{R}^V$ so it can be accessed in $O(1)$ on future iterations that revisit state $J$ and are transitioning for sensor $s$.

As discussed in <cit.>, tempering is required to ensure sufficient exploration such that all significant birth components are generated.
The need for tempering is especially apparent with herding Gibbs sampling, which has been observed to exhibit improved sample efficiency over its stochastic counterparts <cit.>.
To alleviate this, sensor order is lexicographically cycled at each Gibbs iteration, and is denoted by the $\text{Perm}$ function in Algorithm <ref>.
Additionally, all values of $\bar{\psi}^{J}_{Z}(l_+)$ can be capped to avoid poorly scaled transition distributions and birth probabilities, while still allowing for components with low values to be pruned.

The asymptotic complexity of Algorithm <ref> is driven by the implementation of $p(j^{(s)} | J^{-s})$ (see <cit.>).
In the worst-case, $p(j^{(s)} | J^{-s})$ needs to be evaluated the same number of times as the stochastic variant proposed in <cit.>.
The complexity of Algorithm <ref> is equivalent to the complexity of Algorithm 1 in <cit.>.

Adaptive Birth Herded Gibbs Sampler
[1]

$Z$, $r_U$, $p_{B,+}$, $T$, $k$, $V$

$\mathbb{B}'_+ = \emptyset$
$J = \left(0, \dots, 0\right)$
$W = \{\}$, $M = \{\}$
$t = 1, \dots, T$
$s \in \text{Perm}(\{1, \dots, V\})$
$j^{(s)} \in \mathbb{J}^{(s)}_0$
\(\triangleright\) Query the prior weight vector from cache
$(J, s) \notin W$
$\mu_{(J,s)} = p(j^{(s)} | J^{-s})$ <cit.>
$w_{(J,s)} = \mu_{(J,s)}$
$M(J, s) = \mu_{(J,s)}$
$w_{(J,s)} = W(J, s)$
$\mu_{(J,s)} = M(J, s)$

\(\triangleright\) Find the next sample in the Gibbs chain
$j'^{(s)} = \argmax(w_{(J,s)})$
$J' = (j^{(1)}, \dots, j^{(s-1)}, j'^{(s)}, j^{(s+1)}, \dots, j^{(V)})$
$\mathbb{B}'_+ = \mathbb{B}'_+ \cup (k, J')$

\(\triangleright\) Performing the herding update on the weights
$w_{(J,s)} = w_{(J,s)} + \mu$
$w_{(J,s)}^{(s)} = w^{(s)}_{(J,s)} - 1$
$W(J, s) = w_{(J,s)}$

§ SIMULATIONS

§.§ Overview

In this section we provide results from two simulations to demonstrate the performance of the proposed herded Gibbs adaptive birth sampler.
In the first simulation, the herded Gibbs adaptive birth sampler was compared to the stochastic variant in both the Gaussian mixture $\delta$-GLMB and LMB filters.
This simulation used the stochastic update ranked assignment procedure from <cit.>, which is used to produce simulation results in <cit.>.
This verifies the performance and convergence rate of the proposed herded Gibbs sampler.
In the second simulation, we compare fully deterministic variants of the Gaussian mixture $\delta$-GLMB and LMB filters by using a herded Gibbs sampler in the update step.
The herded Gibbs adaptive birth sampler is used with the herded Gibbs update sampler <cit.> and the stochastic Gibbs update sampler <cit.>.
This verifies the performance of a fully deterministic filter.

The scenario tracks a time-varying number of targets in planar 2D position and velocity, $x = [p_{x}, \dot{p}_x, p_y, \dot{p}_y]^T$.
Target dynamics were simulated according to a constant velocity transition model, $f_+(x_+|x) = \mathcal{N}(x_+; Fx, Gw)$ where

\begin{equation*}
F =
\begin{bmatrix}
1 & \Delta\\
0 &     1
\end{bmatrix},\qquad
G =
\begin{bmatrix}
\frac{\Delta^2}{2}\\
\Delta
\end{bmatrix},
\end{equation*}

under the discrete-time sampling interval $\Delta$ <cit.>.
The driving noise $w$ was simulated as zero-mean Gaussian with $\textnormal{Cov}(w) = \textnormal{diag}[5, 5]^T$.
The survival probability for each target was $p_s(x,l) = 0.99$.
The target trajectories were generated according to Figure <ref>, with at most 22 targets occurring at a single time step.
The birth locations of each target were fixed, but are sparse with no two targets being born from the same location.

Eight linear position sensors were simulated, each having detection probability $p_D^{(s)}(x,l) = 0.95$.
If an object was detected, the measurement $z^{(s)} = [p_x, p_y]^T$ was observed according to the single-target measurement likelihood $g(z^{(s)}|x) = \mathcal{N}(z^{(s)}; H^{(s)}x, R^{(s)})$ with $R^{(s)} = \text{diag}(10^2, 10^2)$ and, $H^{(s)} = \begin{bsmallmatrix}1 & 0\\0 & 1\end{bsmallmatrix} \otimes \begin{bsmallmatrix}1 & 0\\0 & 0\end{bsmallmatrix}$ where $\otimes$ denotes the Kronecker product.
Clutter was Poisson distributed with intensity $\kappa^{(s)}(\mathbb{Z}^{(s)}) = \lambda_c^{(s)} \mathcal{U}(\mathbb{Z}^{(s)})$ where $\mathcal{U}(\mathbb{Z}^{(s)})$ is the uniform distribution over $\mathbb{Z}^{(s)}$, and $\lambda^{(s)}_c = 15$.

Each filter was simulated over 100 Monte Carlo iterations.
The prior was modeled as a Gaussian with an uninformative covariance in position such that, $p_B(x, l_+) = \mathcal{N}(x; \mu_0, P_0)$ with $\mu_0 = \left[0, 0, 0, 0\right]^T$ and $P_0 = \text{diag}(100000^2, 50^2, 100000^2, 50^2)$.
Both adaptive birth Gibbs samplers were configured to use 250 iterations per evaluation, which is reduced from 1000 iterations as used in <cit.>.
The remaining adaptive birth parameters were kept the same.

§.§ Results

The cardinality and state estimation accuracy were estimated using the OSPA(2) metric <cit.> as depicted in Figure <ref> and Figure <ref>.
The OSPA(2) metric was computed using a distance cutoff value of $200$, a distance order of $1.0$, a sliding window length of $5$ and an expanding window weight power of $0$.

As shown in Figure <ref>, the proposed herded Gibbs adaptive birth sampler had a lower average OSPA(2) score in the $\delta$-GLMB filter simulations.
Both Gibbs samplers performed similarly for the LMB filter.
These results show that the herded Gibbs sampler converged faster, using fewer iterations, and achieved similar performance to <cit.>.

The average OSPA(2) scores in Figure <ref> is similar for both the $\delta$-GLMB and LMB filters.
The best performance is achieved using the herded Gibbs sampler for birth, with either Gibbs sampler for ranked assignment.

Single run results using deterministic multi-sensor adaptive birth model.
Labeled target position estimates (colored circles, each color corresponds to a unique track label), target trajectories (black lines), target birth locations (black circles).

[$\delta$-GLMB filter]
[LMB filter]
OSPA(2) results comparing herded and stochastic adaptive birth Gibbs samplers. These simulations were configured to use the stochastic Gibbs sampler for ranked assignment to match <cit.>.

[$\delta$-GLMB filter]
[LMB filter]
OSPA(2) results comparing herded and stochastic ranked assignment Gibbs samplers. These simulations were configured to use the herded Gibbs sampler for birth.

§ CONCLUSIONS

This paper provided a deterministic Gibbs truncation algorithm for multi-sensor measurement-adaptive birth.
This removes stochasticity from a tracking system, reducing the run-to-run variability of a tracker's performance.
The main advantage of this approach is that it guarantees reproducible results, which makes it more suitable for use in safety critical operations.
The results of the simulations show the average tracking performance was similar to or better than the stochastic Gibbs sampler, and uses fewer iterations than<cit.>.
Future work includes analysis to find the minimum required iterations to match the stochastic ranked assignment Gibbs sampler performance.


aainstitutetext: Institute of Nuclear and Particle Physics, NCSR Demokritos,
Patr. Grigoriou E’ & 27 Neapoleos Str, 15341 Agia Paraskevi,
Greecebbinstitutetext: Physics Division, National Technical University of
Athens,
Zografou Campus, Athens 15780, Greece

# One-loop Feynman integrals for $2\to 3$ scattering involving many scales
including internal masses

Nikolaos Syrrakos<EMAIL_ADDRESS>
###### Abstract

We study several multiscale one-loop five-point families of Feynman integrals.
More specifically, we employ the Simplified Differential Equations approach to
obtain results in terms of Goncharov polylogarithms of up to transcendental
weight four for families with two and three massive external legs and massless
propagators, as well as with one massive internal line and up to two massive
external legs. This is the first time this computational approach is applied
to cases involving internal masses.

###### Keywords:

Feynman integrals, QCD, NNLO Calculations

##  1 Introduction

In recent years the field of precision calculations in collider physics has
emerged as a vibrant and fruitful line of research in our attempt to
understand Nature at its most fundamental level Heinrich:2020ybq . The basic
principle of this research endeavour is to have very precise experimental
measurements of cross sections for Standard Model scattering processes
compared against theoretical predictions of equally high precision and search
for any deviations between them. Should any such deviations be established,
their analysis and physical explanation would require New Physics, giving us
an idea of what lies beyond the Standard Model.

The ever increasing demand for highly precise theoretical predictions for
scattering processes relevant to LHC searches poses a challenge in our ability
to perform higher order calculations in perturbative Quantum Field Theory
Amoroso:2020lgh . Multiloop scattering amplitudes play a fundamental role in
such calculations, encoding within their mathematical structure key
information concerning the nature of particle interactions. One major aspect
of the calculation of multiloop scattering amplitudes is the calculation of
the relevant Feynman diagrams that are involved, which can be associated
through the corresponding Feynman rules to the so-called Feynman integrals. In
the following we will use the notion of diagrams and integrals
interchangeably.

The standard approach for the calculation of these integrals involves
obtaining a complete set of Master integrals through the use of Integration-
By-Part identities Chetyrkin:1981qh , constructing a pure basis of Master
integrals Henn:2014qga  and then deriving and solving differential equations
de1 ; de2 ; de3 ; de4  in canonical form Henn:2013pwa . This approach has
yielded numerous results Kotikov:2021tai , in part due to the fact that we
have a solid understanding of the special class of functions, known as
multiple or Goncharov polylogarithms Goncharov:1998kja ; Duhr:2011zq ;
Duhr:2012fh ; Duhr:2014woa , in terms of which many Feynman integrals can be
expressed. In more complicated cases however, this class of functions is not
enough and important steps have been made in getting a better understanding of
a more general class of functions, Elliptic integrals Remiddi:2017har ;
Broedel:2017kkb ; Broedel:2017siw ; Broedel:2018iwv ; Broedel:2018qkq ;
Broedel:2019hyg ; Duhr:2019rrs , which appear in solutions of multiloop
Feynman integrals with many scales, especially when several internal masses
are introduced.

When considering multiloop Feynman integrals involving many external
particles, the current frontier lies at two-loop five-point integrals with up
to one off-shell leg and massless internal lines. For the fully massless case,
all Master integrals are by now known up to transcendental weight four
Papadopoulos:2015jft ; Gehrmann:2015bfy ; Chicherin:2017dob ; Gehrmann:2018yef
; Chicherin:2018mue ; Abreu:2018rcw ; Chicherin:2018old  and their solutions
have been implemented in a fast C++ library known as pentagon functions
Chicherin:2020oor . When one of the external particles is considered off-
shell, the planar topologies have been recently solved using two different
computational approaches for the solution of canonical differential equations,
numerically Abreu:2020jxa  and analytically Canko:2020ylt . The numerical
calculation was performed using a generalised power-series method
Francesco:2019yqt ; Hidding:2020ytt , while the analytical solution was
achieved through the use of the Simplified Differential Equations approach
Papadopoulos:2014lla , with the results given in terms of Goncharov
polylogarithms of up to transcendental weight four. These results are relevant
to many $2\to 3$ scattering processes studied experimentally at the LHC, e.g.
$W+2$ jets production. For the computation of the relevant scattering
amplitudes, one-loop five-point Feynman integrals with one off-shell leg also
have to be known up to transcendental weight four Syrrakos:2020kba . These
results were recently used for the calculation of two-loop QCD corrections to
$Wb\bar{b}$ production Badger:2021nhg . First results for one of the non-
planar topologies have also appeared using a numerical approach
Papadopoulos:2019iam . More recently the three hexabox topologies were
calculated in Abreu:2021smk  using the same approach as in Abreu:2020jxa .

While staying at the level of five-point Feynman integrals, at some point we
will have to introduce internal masses and consider more than one massive
external particle. Judging from the level of complexity of the so far
accumulated results, these Feynman integrals are expected to be highly non-
trivial to be solved using current approaches, especially the genuine two-loop
ones. To that end, we believe that it is instructive to consider first the
relevant one-loop five-point Feynman integrals with more than one off-shell
leg and/or with internal masses. The interest of these Feynman integrals is
twofold. From a more formal point of view, it is interesting to see what kind
of functions appear as solutions of the relevant canonical differential
equations and study their structure. This will give us a glimpse of the
minimum mathematical complexities and difficulties we should expect when we
consider their two-loop counterparts. From a phenomenological standpoint,
these one-loop integrals will be required for the computation of two-loop
corrections for $2\to 3$ scattering processes involving more than one massive
external particle and/or internal massive particles.

In line with the arguments presented above, we consider in this paper the
analytical calculation of several multiscale one-loop five-point Feynman
integrals. More specifically, we present analytical results in terms of
Goncharov polylogarithms of up to transcendental weight four for families with
two and three off-shell legs and massless internal lines, as well as for
families with one massive propagator and up to two external massive particles.
Our calculation is based on the Simplified Differential Equations approach
Papadopoulos:2014lla  which introduces an external dimensionless parameter $x$
in such a way that captures the off-shellness of one massive leg. The system
of canonical differential equations is constructed by differentiating a pure
basis of Master integrals in terms of $x$, regardless of the number of scales
involved in the problem. One special feature of this approach is that by
taking the limit of $x\to 1$ Papadopoulos:2015jft  for a family with $n$
massive legs, we can obtain the result for a family with $n-1$ massive legs in
an algorithmic way Canko:2020gqp . In Figure 1 we present the families of
Feynman integrals computed through the solution of simplified differential
equations in canonical form, while in Figure 2 we present the families of
integrals computed through the $x\to 1$ limit.

|  |  |
---|---|---|---

(a) C (b) E

(c) G (d) H

Figure 1: Top-sector diagrams for families computed with the SDE approach. All
external particles are incoming. Bold external (internal) lines represent
massive particles.

|
---|---

(a) D (b) F

Figure 2: Top-sector diagrams for families computed through the $x\to 1$
limit. All external particles are incoming. Bold external (internal) lines
represent massive particles.

The rest of our paper is structured as follows: in section 2 we introduce
basic notation and the kinematic configuration for each of the studied
families of Feynman integrals, in section 3 we construct pure bases and derive
and solve simplified differential equations in canonical form for all integral
families depicted in Figure 1, we present some of the resulting alphabets in
$x$ and study their structure, and solve all integral families depicted in
Figure 2 through the $x\to 1$ limit in terms of Goncharov polylogarithms of up
to transcendental weight four. In section 4 we provide an analysis of our
results, as well as numerical checks for Euclidean points and in section 5 we
summarise our findings and discuss their key features. To the best of our
knowledge these families have never before been considered in the literature,
thus their solution constitutes an original contribution. This is also the
first time that a calculation with the Simplified Differential Equations
approach involving Master integrals with internal masses is reported. Along
with this paper we provide all of our results in ancillary files. Explicit
weight-three expressions for the top-sector basis elements of families $C$ and
$H$ are given in appendix A.

##  2 Notation and kinematics

The integral families are defined through the following parametrization,

$G_{a_{1}a_{2}a_{3}a_{4}a_{5}}=\int~{}\frac{\mathrm{d}^{d}k_{1}}{i\pi^{(d/2)}}~{}\frac{\mathrm{e}^{\epsilon\gamma_{E}}}{\mathcal{D}_{1}^{a_{1}}\mathcal{D}_{2}^{a_{2}}\mathcal{D}_{3}^{a_{3}}\mathcal{D}_{4}^{a_{4}}\mathcal{D}_{5}^{a_{5}}},\quad
d=4-2\epsilon$ (1)

with

$\displaystyle\mathcal{D}_{1}=(k_{1})^{2}-n_{1}~{}m^{2},~{}\mathcal{D}_{2}=(k_{1}+q_{1})^{2},~{}\mathcal{D}_{3}=(k_{1}+q_{1}+q_{2})^{2}$
$\displaystyle\mathcal{D}_{4}=(k_{1}+q_{1}+q_{2}+q_{3})^{2}-n_{4}~{}m^{2},~{}\mathcal{D}_{5}=(k_{1}+q_{1}+q_{2}+q_{3}+q_{4})^{2}$
(2)

For the families $C,D$ and $E$ we have $n_{1}=n_{4}=0$, for the families $F$
and $G$ $n_{1}=0,n_{4}=1$ and finally for the family $H$ $n_{1}=1,n_{4}=0$.
The kinematics for the families depicted in Figure 1 is as follows,

  * •

$C\&H$: $\sum_{i=1}^{5}q_{i}=0$, $q_{i}^{2}=0$, $i=1,2,4$,
$q_{3}^{2}=m_{3}^{2},q_{5}^{2}=m_{5}^{2}$

  * •

$E$: $\sum_{i=1}^{5}q_{i}=0$, $q_{i}^{2}=0$, $i=2,4$,
$q_{1}^{2}=\bar{m}_{1}^{2},q_{3}^{2}=m_{3}^{2},q_{5}^{2}=m_{5}^{2}$

  * •

$G$: $\sum_{i=1}^{5}q_{i}=0$, $q_{i}^{2}=0$, $i=1,2,4,5$,
$q_{3}^{2}=m_{3}^{2}$

We introduce the following $x$-parametrization111We use the abbreviations
$p_{ij}=p_{i}+p_{j}$ and $p_{ijk}=p_{i}+p_{j}+p_{k}$ and similarly for $q$
later.

$q_{1}=xp_{1},~{}q_{2}=xp_{2},~{}q_{3}=p_{123}-xp_{12},~{}q_{4}=p_{4},~{}q_{5}=-p_{1234}$
(3)

The kinematics in this underline momentum parametrization is

  * •

$C\&H$: $\sum_{i=1}^{5}p_{i}=0$, $p_{i}^{2}=0$, $i=1,2,3,4$,
$p_{5}^{2}=m_{5}^{2}$

  * •

$E$: $\sum_{i=1}^{5}p_{i}=0$, $p_{i}^{2}=0$, $i=2,3,4$,
$p_{1}^{2}=m_{1}^{2},p_{5}^{2}=m_{5}^{2}$

  * •

$G$: $\sum_{i=1}^{5}p_{i}=0$, $p_{i}^{2}=0$, $i=1,2,3,4,5$

Introducing (3) results in a mapping between the kinematic invariants in the
original momentum parametrization, $q_{i}$, and the underline momentum
parametrization $\\{x,p_{i}\\}$ for each of the families $C,E,G,H$.222We use
the abbreviations $s_{ij}=q_{ij}^{2},~{}S_{ij}=p_{ij}^{2}$.

$\displaystyle C\&H:$
$\displaystyle~{}m_{3}^{2}=(x-1)\left(S_{12}x-S_{45}\right),s_{12}=S_{12}x^{2},s_{23}=S_{23}x-S_{45}x+S_{45}$
$\displaystyle
s_{34}=m_{5}^{2}(-x)+m_{5}^{2}+x\left(S_{12}(x-1)+S_{34}\right),s_{45}=S_{45},s_{15}=m_{5}^{2}(-x)+m_{5}^{2}+S_{15}x$
$\displaystyle E:$
$\displaystyle~{}s_{12}=S_{12}x^{2},s_{23}=x\left(m_{1}^{2}(x-1)+S_{23}\right)-S_{45}x+S_{45},$
$\displaystyle
s_{34}=m_{5}^{2}(-x)+m_{5}^{2}+x\left(S_{12}(x-1)+S_{34}\right),s_{45}=S_{45},$
$\displaystyle
s_{15}=x\left(m_{1}^{2}(x-1)+S_{15}\right)+m_{5}^{2}(-x)+m_{5}^{2},\bar{m}_{1}^{2}=m_{1}^{2}x^{2},m_{3}^{2}=(x-1)\left(S_{12}x-S_{45}\right)$
$\displaystyle G:$
$\displaystyle~{}s_{12}=S_{12}x^{2},s_{23}=S_{23}x-S_{45}x+S_{45},s_{34}=x\left(S_{12}(x-1)+S_{34}\right),s_{45}=S_{45},$
$\displaystyle s_{15}=S_{15}x,m_{3}^{2}=(x-1)\left(S_{12}x-S_{45}\right)$ (4)

For the families depicted in Figure 2 their definition through (1)&(2) is
obtained by taking (3) and setting $x=1$, and their kinematic configuration is
effectively the one produced by the underline momentum parametrization of the
families through which we will calculate them with the $x\to 1$ limit,
therefore we have

  * •

$D~{}(x\to 1~{}\text{of}~{}E)$: $\sum_{i=1}^{5}p_{i}=0$, $p_{i}^{2}=0$,
$i=2,3,4$, $p_{1}^{2}=m_{1}^{2},p_{5}^{2}=m_{5}^{2}$

  * •

$F~{}(x\to 1~{}\text{of}~{}G)$: $\sum_{i=1}^{5}p_{i}=0$, $p_{i}^{2}=0$,
$i=1,2,3,4,5$

##  3 Differential equations and pure solutions

In this section we will describe the analytical calculation of the integral
families considered in this paper. We will construct pure bases for the
families $C,E,G,H$ and use the Simplified Differential Equations approach
Papadopoulos:2014lla  to compute them in terms of Goncharov polylogarithms. We
will show how to obtain boundary terms Canko:2020ylt ; Canko:2020gqp  for
these canonical differential equations and present explicit results up to
transcendental weight four. A discussion on the structure of the alphabets in
$x$ for families $C$ and $H$ is also provided, along with their explicit
expressions.

We will also show how by taking the $x\to 1$ limit of the analytic solution of
families $E$ and $G$, one can obtain in an algorithmic way analytic results in
terms of Goncharov polylogarithms for the families $D$ and $F$ respectively.
Additionally, this procedure will allow us to obtain pure bases for the
families $D$ and $F$ in a straightforward manner Canko:2020gqp . The results
presented here for these last two families are up to transcendental weight
four as well.

It is important to note that introducing one internal mass does not appear to
have any effect on the efficiency of the methods that have been developed for
the determination of boundary terms, as well as taking the $x\to 1$ limit
within the Simplified Differential Equations approach.

###  3.1 Families $C,E,G,H$

Constructing pure bases for the families $C,E,G,H$ is by now a trivial
exercise. Following the reasoning of Abreu:2018rcw ; Abreu:2020jxa , the top
sector basis element at the integrand level is of the form

$\epsilon^{2}\frac{\mathcal{P}_{11111}}{\sqrt{\Delta_{5}}}\tilde{G}_{11111}$
(5)

where $\mathcal{P}_{11111}$ is the Baikov polynomial corresponding to the top
sector integral $G_{11111}$ for each family, $\tilde{G}_{11111}$ is the top
sector integrand of each family and $\Delta_{5}=\det[q_{i}\cdot q_{j}]$ is the
Gram determinant of the external momenta. The remaining pure basis elements
can be constructed through the study of the leading singularities of their
corresponding diagrams Henn:2014qga . Using Azurite Georgoudis:2016wff  and
Kira2 Klappert:2020nbg  we identify 15, 18, 16 and 18 Master Integrals for the
families $C,E,G,H$ respectively.

When considering five-point scattering, a number of square roots of the
kinematic invariants enter the differential equations of the corresponding
pure bases. These square roots originate from leading singularities of
triangles with three massive legs333For all integral families considered in
this paper the most complicated triangle Feynman integrals are the ones with
fully massive legs and one massive propagator. This one internal mass however
has no effect in the calculation of the leading singularity of the
corresponding integral. which are represented by square roots of the Källen
function $\lambda(x,y,z)=x^{2}-2xy-2xz+y^{2}-2yz+z^{2}$ and from square roots
of the Gram determinants of the five-point external momenta. The existence of
these square roots poses a challenge if one wishes to solve the differential
equations analytically in terms of Goncharov polylogarithms.

For the families considered in this subsection the following square roots
appear:

$\displaystyle r_{1}$ $\displaystyle=\sqrt{\lambda(s_{12},m_{3}^{2},s_{45})}$
(6) $\displaystyle r_{2}$
$\displaystyle=\sqrt{\lambda(s_{12},s_{34},m_{5}^{2})}$ (7) $\displaystyle
r_{3}$ $\displaystyle=\sqrt{\lambda(\bar{m}_{1}^{2},s_{23},s_{45})}$ (8)
$\displaystyle r_{4}$
$\displaystyle=\sqrt{\lambda(\bar{m}_{1}^{2},m_{5}^{2},s_{15})}$ (9)
$\displaystyle r_{5}$
$\displaystyle=\sqrt{\Delta_{5}^{C}}=\sqrt{\Delta_{5}^{H}}$ (10)
$\displaystyle r_{6}$ $\displaystyle=\sqrt{\Delta_{5}^{E}}$ (11)
$\displaystyle r_{7}$ $\displaystyle=\sqrt{\Delta_{5}^{G}}$ (12)

We should note at this point that not all square roots appear in every family
at the same time. More specifically, in families $C$ and $H$ we encounter
$r_{1},r_{2},r_{5}$, in family $E$ $r_{1},r_{2},r_{3},r_{4},r_{6}$ and in
family $G$ $r_{1},r_{7}$.

If one tries to compute these families using the standard differential
equations approach, i.e. by differentiating with respect to all kinematic
invariants, then the algebraic structure of the alphabet (i.e. the square
roots appearing as letters) of the canonical differential equation for each
family prohibits a straightforward solution in terms of Goncharov
polylogarithms. In order to achieve a result in such a form we need to find a
way to deal with these square roots. Several ideas have been put forward
recently that are able to circumvent this problem and provide solutions in
terms of Goncharov polylogarithms Heller:2019gkq ; Besier:2018jen ;
Besier:2019kco ; Bonetti:2020hqh  for specific cases, however a universal
method to treat the problem of square roots appearing in the alphabet of
canonical differential equations for multiscale families of Feynman integrals
is still missing.

It turns out that for the families $C,E,G,H$ the $x$-parametrization
introduced in (3) and the resulting mapping of the kinematic invariants (2)
rationalises all square roots with respect to $x$. This allows us to obtain a
canonical differential equation in $x$ for each of the four families
considered in this subsection,

$\partial_{x}\textbf{g}=\epsilon\left(\sum_{i=1}^{l_{max}}\frac{\textbf{M}_{i}}{x-l_{i}}\right)\textbf{g}$
(13)

where g is the pure basis for each family, $\textbf{M}_{i}$ are the residue
matrices corresponding to each letter $l_{i}$ and $l_{max}$ is the length of
the alphabet444Note here that we are following the notation of Canko:2020ylt ;
Syrrakos:2020kba ; Canko:2020gqp  when talking about the letters of the
alphabet.. The kinematic dependence is entirely contained within the letters
$l_{i}$, leaving the residue matrices $\textbf{M}_{i}$ to be solely
constructed by rational numbers. The length of the alphabet for each of the
four families considered in this subsection is
$l_{max}^{C}=14,~{}l_{max}^{E}=19,~{}l_{max}^{G}=22,~{}l_{max}^{H}=30$. The
explicit form of the alphabet for each of the four families is provided in the
ancillary files that accompany this paper. In the next subsection we will
study more closely some of these alphabets.

In order to solve (13) we need to provide boundary terms. We will follow
closely the computational framework developed in Canko:2020ylt ; Canko:2020gqp
for the determination of the relevant boundary terms. We start with the
residue matrix corresponding to the letter $\\{0\\}$, $\textbf{M}_{1}$ and
through its Jordan Decomposition we rewrite it as follows,

$\textbf{M}_{1}=\textbf{S}\textbf{D}\textbf{S}^{-1}$ (14)

Then we define the resummation matrix R as follows

$\textbf{R}=\textbf{S}e^{\epsilon\textbf{D}\log(x)}\textbf{S}^{-1}$ (15)

The next step is to use IBP identities to write the pure basis g in the
following form

$\textbf{g}=\textbf{T}\textbf{G}$ (16)

The list of Feynman integrals G is also provided in electronic form.
Furthermore, using the expansion-by-regions method Jantzen:2012mw  implemented
in the asy code which is shipped along with FIESTA4 Smirnov:2015mct , we can
obtain information for the asymptotic behaviour of the Feynman integrals in
terms of which we express the pure basis of Master integrals (16) in the limit
$x\to 0$,

${G_{i}}\mathop{=}\limits_{x\to
0}\sum\limits_{j}x^{b_{j}+a_{j}\epsilon}G^{(b_{j}+a_{j}\epsilon)}_{i}$ (17)

where $a_{j}$ and $b_{j}$ are integers and $G_{i}$ are the individual members
of the basis G of Feynman integrals in (16). As explained in Canko:2020ylt ,
we can construct the relation

$\mathbf{R}\mathbf{b}=\left.\lim_{x\rightarrow
0}\mathbf{T}\mathbf{G}\right|_{\mathcal{O}\left(x^{0+a_{j}\epsilon}\right)}$
(18)

where $\textbf{b}=\sum_{i=0}^{n}\epsilon^{i}\textbf{b}_{0}^{(i)}$ are the
boundary terms that we need to compute. The right-hand-side of (18) implies
that, apart from the terms $x^{a_{i}\epsilon}$ coming from (17), we expand
around $x=0$, keeping only terms of order $x^{0}$.

Equation (18) allows us to fix all the necessary boundary terms without the
need of any further computation for the families $C,E,H$ while for family $G$
a few regions had to be computed. Similarly to Syrrakos:2020kba , the
resulting boundary terms for all of the four families considered in this
subsection are in closed form, including some ${}_{2}F_{1}$ Hypergeometric
functions which can be easily expanded to arbitrary powers of the dimensional
regulator using HypExp Huber:2007dx . Therefore we are able to trivially
obtain solutions of (13) for the families $C,E,G,H$ in terms of Goncharov
polylogarithms of arbitrary weight.

In this paper we present explicit results for the families $C,E,G,H$ in terms
of Goncharov polylogarithms of up to transcendental weight four, which can be
written in the following compact form,

g
$\displaystyle=\epsilon^{0}\textbf{b}^{(0)}_{0}+\epsilon\bigg{(}\sum\mathcal{G}_{a}\textbf{M}_{a}\textbf{b}^{(0)}_{0}+\textbf{b}^{(1)}_{0}\bigg{)}$
$\displaystyle+\epsilon^{2}\bigg{(}\sum\mathcal{G}_{ab}\textbf{M}_{a}\textbf{M}_{b}\textbf{b}^{(0)}_{0}+\sum\mathcal{G}_{a}\textbf{M}_{a}\textbf{b}^{(1)}_{0}+\textbf{b}^{(2)}_{0}\bigg{)}$
$\displaystyle+\epsilon^{3}\bigg{(}\sum\mathcal{G}_{abc}\textbf{M}_{a}\textbf{M}_{b}\textbf{M}_{c}\textbf{b}^{(0)}_{0}+\sum\mathcal{G}_{ab}\textbf{M}_{a}\textbf{M}_{b}\textbf{b}^{(1)}_{0}+\sum\mathcal{G}_{a}\textbf{M}_{a}\textbf{b}^{(2)}_{0}+\textbf{b}^{(3)}_{0}\bigg{)}$
$\displaystyle+\epsilon^{4}\bigg{(}\sum\mathcal{G}_{abcd}\textbf{M}_{a}\textbf{M}_{b}\textbf{M}_{c}\textbf{M}_{d}\textbf{b}^{(0)}_{0}+\sum\mathcal{G}_{abc}\textbf{M}_{a}\textbf{M}_{b}\textbf{M}_{c}\textbf{b}^{(1)}_{0}$
$\displaystyle+\sum\mathcal{G}_{ab}\textbf{M}_{a}\textbf{M}_{b}\textbf{b}^{(2)}_{0}+\sum\mathcal{G}_{a}\textbf{M}_{a}\textbf{b}^{(3)}_{0}+\textbf{b}^{(4)}_{0}\bigg{)}$
(19)

were $\mathcal{G}_{ab\ldots}:=\mathcal{G}(l_{a},l_{b},\ldots;x)$ represent the
Goncharov polylogarithms. The $b_{0}^{(i)}$ terms, with $i$ indicating the
corresponding weight, consist of Zeta functions $\zeta(i)$, logarithms and
Goncharov polylogarithms of weight $i$ which have as arguments rational
functions of the underline kinematic variables $S_{ij}$.

Our results are presented in such a way that each coefficient of
$\epsilon^{k}$ has transcendental weight $k$. If we assign weight $-1$ to
$\epsilon$, then (3.1) has uniform weight zero. The closed-form expressions of
the boundary terms trivialise the extension of (3.1) to higher transcendental
weights (or higher orders in $\epsilon$).

###  3.2 The alphabet in $x$

It is instructive to have a closer look at the alphabets for some of these
families and see what lessons can be learned. We will study the alphabets of
families $C$ and $H$. We choose these families because they have the same
external kinematics but differ on the fact that family $H$ has one internal
mass. Thus it is interesting to see how the introduction of an internal mass
effects the alphabet. The conclusions drawn from the study of these families’
alphabets are similar to what can be learned from the alphabets of the rest of
the families considered in this paper.

#### Family $C$

The alphabet for this family is

$\displaystyle l_{1}\to 0,l_{2}\to
1,l_{3}\to\frac{m^{2}_{5}}{m^{2}_{5}-S_{15}},l_{4}\to\frac{S_{12}+S_{23}}{S_{12}},l_{5}\to\frac{S_{12}+S_{15}-S_{34}}{S_{12}},$
$\displaystyle
l_{6}\to\frac{-\sqrt{\Delta_{1}}+m^{2}_{5}+S_{12}-S_{34}}{2S_{12}},l_{7}\to\frac{\sqrt{\Delta_{1}}+m^{2}_{5}+S_{12}-S_{34}}{2S_{12}},l_{8}\to\frac{m^{2}_{5}-S_{45}}{m^{2}_{5}-S_{15}+S_{23}-S_{45}},$
$\displaystyle
l_{9}\to\frac{m^{2}_{5}-S_{45}}{m^{2}_{5}-S_{34}-S_{45}},l_{10}\to\frac{S_{45}}{S_{12}},l_{11}\to\frac{S_{45}}{S_{45}-S_{23}},l_{12}\to\frac{m^{2}_{5}S_{12}-S_{45}S_{12}+S_{34}S_{45}}{m^{2}_{5}S_{12}-S_{12}S_{45}},$
$\displaystyle
l_{13}\to-\frac{\sqrt{\Delta_{2}}-2m^{2}_{5}S_{12}-m^{2}_{5}S_{23}+S_{15}S_{12}-S_{23}S_{12}+2S_{45}S_{12}+S_{23}S_{34}+S_{15}S_{45}-S_{34}S_{45}}{2S_{12}\left(m^{2}_{5}-S_{15}+S_{23}-S_{45}\right)},$
$\displaystyle
l_{14}\to\frac{\sqrt{\Delta_{2}}+2m^{2}_{5}S_{12}+m^{2}_{5}S_{23}-S_{15}S_{12}+S_{23}S_{12}-2S_{45}S_{12}-S_{23}S_{34}-S_{15}S_{45}+S_{34}S_{45}}{2S_{12}\left(m^{2}_{5}-S_{15}+S_{23}-S_{45}\right)}$
(20)

We notice that two square roots appear,
$\sqrt{\Delta_{1}},\,\sqrt{\Delta_{2}}$, where $\Delta_{1},\,\Delta_{2}$ are
given by the following expressions

$\displaystyle\Delta_{1}=$
$\displaystyle-2m^{2}_{5}\left(S_{12}+S_{34}\right)+\left(m^{2}_{5}\right)^{2}+\left(S_{12}-S_{34}\right){}^{2},$
(21) $\displaystyle\Delta_{2}=$
$\displaystyle\left(m^{2}_{5}S_{23}-S_{34}S_{23}+\left(S_{34}-S_{15}\right)S_{45}\right){}^{2}+S_{12}^{2}\left(S_{15}-S_{23}\right){}^{2}$
$\displaystyle+2S_{12}\big{(}m^{2}_{5}S_{23}\left(S_{15}-S_{23}-2S_{34}\right)-S_{45}S_{15}^{2}+S_{23}S_{34}S_{15}+\left(S_{23}+S_{34}\right)S_{45}S_{15}$
$\displaystyle+S_{23}S_{34}\left(S_{45}-S_{23}\right)\big{)}$ (22)

These square roots are directly associated with $r_{2}$ and $r_{5}$
respectively. We would expect a third square root to appear, namely $r_{1}$,
as mentioned in the previous subsection. However, after introducing the
$x$-parametrization (3), $r_{1}$ is rational in all variables. More
specifically, if we express these square roots using the first set of
relations from (2), we get

$\displaystyle r_{1}=$ $\displaystyle\left(S_{12}-S_{45}\right)x$ (23)
$\displaystyle r_{2}=$ $\displaystyle\sqrt{\Delta_{1}}\,x$ (24) $\displaystyle
r_{5}=$ $\displaystyle\sqrt{\Delta_{2}}\,x^{2}$ (25)

The structure of this alphabet is similar in terms of its complexity with the
alphabet for the one-mass pentagon studied in Syrrakos:2020kba . The
difference of course is the presence of an additional square root of the
underline kinematic variables $S_{ij}$ in the present case.

#### Family $H$

The alphabet for this family is

$\displaystyle l_{1}\to 0,l_{2}\to
1,l_{3}\to\frac{m^{2}}{S_{12}},l_{4}\to-\frac{\sqrt{m^{2}}}{\sqrt{S_{12}}},l_{5}\to\frac{\sqrt{m^{2}}}{\sqrt{S_{12}}},l_{6}\to\frac{S_{12}+S_{23}}{S_{12}},$
$\displaystyle
l_{7}\to\frac{m^{2}_{5}+S_{12}-S_{34}-\sqrt{\hat{\Delta}_{1}}}{2S_{12}},l_{8}\to\frac{m^{2}_{5}+S_{12}-S_{34}+\sqrt{\hat{\Delta}_{1}}}{2S_{12}},l_{9}\to\frac{m^{2}\left(m^{2}_{5}+S_{12}-S_{34}-\sqrt{\hat{\Delta}_{1}}\right)}{2m^{2}_{5}S_{12}},$
$\displaystyle
l_{10}\to\frac{m^{2}\left(m^{2}_{5}+S_{12}-S_{34}+\sqrt{\hat{\Delta}_{1}}\right)}{2m^{2}_{5}S_{12}},l_{11}\to-\frac{m^{2}}{S_{23}-S_{45}},l_{12}\to\frac{m^{2}_{5}-S_{45}}{m^{2}_{5}-S_{34}-S_{45}},l_{13}\to\frac{m^{2}}{S_{45}},$
$\displaystyle
l_{14}\to\frac{S_{45}}{S_{12}},l_{15}\to-\frac{S_{45}}{S_{23}-S_{45}},l_{16}\to\frac{-S_{12}m^{2}-S_{12}S_{45}+\sqrt{\hat{\Delta}_{2}}}{2S_{12}\left(S_{23}-S_{45}\right)},l_{17}\to-\frac{S_{12}m^{2}+S_{12}S_{45}+\sqrt{\hat{\Delta}_{2}}}{2S_{12}\left(S_{23}-S_{45}\right)},$
$\displaystyle
l_{18}\to-\frac{-m^{2}_{5}m^{2}+S_{34}m^{2}+S_{45}m^{2}-m^{2}_{5}S_{12}+S_{12}S_{45}-S_{34}S_{45}+\sqrt{\hat{\Delta}_{3}}}{2S_{12}\left(m^{2}_{5}-S_{45}\right)},$
$\displaystyle
l_{19}\to\frac{m^{2}_{5}m^{2}-S_{34}m^{2}-S_{45}m^{2}+m^{2}_{5}S_{12}-S_{12}S_{45}+S_{34}S_{45}+\sqrt{\hat{\Delta}_{3}}}{2S_{12}\left(m^{2}_{5}-S_{45}\right)},$
$\displaystyle
l_{20}\to\frac{m^{2}}{m^{2}_{5}-S_{15}},l_{21}\to\frac{m^{2}_{5}}{m^{2}_{5}-S_{15}},l_{22}\to\frac{m^{2}_{5}-S_{45}}{m^{2}_{5}-S_{15}+S_{23}-S_{45}},l_{23}\to\frac{S_{12}+S_{15}-S_{34}}{S_{12}},$
$\displaystyle
l_{24}\to\frac{m^{2}\left(m^{2}_{5}-S_{45}\right)}{\left(m^{2}_{5}-S_{15}+S_{23}-S_{45}\right)m^{2}-m^{2}_{5}S_{23}+S_{15}S_{45}},l_{25}\to\frac{S_{12}m^{2}+m^{2}_{5}S_{12}-\sqrt{\hat{\Delta}_{4}}}{2S_{12}\left(m^{2}_{5}-S_{15}\right)},$
$\displaystyle
l_{26}\to\frac{S_{12}m^{2}+m^{2}_{5}S_{12}+\sqrt{\hat{\Delta}_{4}}}{2S_{12}\left(m^{2}_{5}-S_{15}\right)},$
$\displaystyle
l_{27}\to\frac{-S_{12}S_{15}m^{2}+m^{2}_{5}S_{23}m^{2}+S_{12}S_{23}m^{2}-S_{23}S_{34}m^{2}-S_{15}S_{45}m^{2}+S_{34}S_{45}m^{2}+\sqrt{\hat{\Delta}_{5}}}{2S_{12}\left(m^{2}_{5}S_{23}-S_{15}S_{45}\right)},$
$\displaystyle
l_{28}\to-\frac{S_{12}S_{15}m^{2}-m^{2}_{5}S_{23}m^{2}-S_{12}S_{23}m^{2}+S_{23}S_{34}m^{2}+S_{15}S_{45}m^{2}-S_{34}S_{45}m^{2}+\sqrt{\hat{\Delta}_{5}}}{2S_{12}\left(m^{2}_{5}S_{23}-S_{15}S_{45}\right)},$
$\displaystyle
l_{29}\to\frac{2m^{2}_{5}S_{12}-S_{15}S_{12}+S_{23}S_{12}-2S_{45}S_{12}+m^{2}_{5}S_{23}-S_{23}S_{34}-S_{15}S_{45}+S_{34}S_{45}+\sqrt{\hat{\Delta}_{6}}}{2S_{12}\left(m^{2}_{5}-S_{15}+S_{23}-S_{45}\right)},$
$\displaystyle
l_{30}\to-\frac{-2m^{2}_{5}S_{12}+S_{15}S_{12}-S_{23}S_{12}+2S_{45}S_{12}-m^{2}_{5}S_{23}+S_{23}S_{34}+S_{15}S_{45}-S_{34}S_{45}+\sqrt{\hat{\Delta}_{6}}}{2S_{12}\left(m^{2}_{5}-S_{15}+S_{23}-S_{45}\right)}$
(26)

The first remark that we can make for this alphabet is that we have six square
roots in the underline kinematic variables $S_{ij}$, namely
$\sqrt{\hat{\Delta}_{i}},\,i=1,\ldots 6$. The explicit expressions for the
arguments of these square roots are as follows,

$\displaystyle\hat{\Delta}_{1}=$
$\displaystyle-2m_{5}^{2}\left(S_{12}+S_{34}\right)+m_{5}^{4}+\left(S_{12}-S_{34}\right){}^{2},$
(27) $\displaystyle\hat{\Delta}_{2}=$ $\displaystyle
S_{12}\left(m^{4}S_{12}+4m^{2}S_{23}\left(S_{12}+S_{23}\right)-2m^{2}\left(S_{12}+2S_{23}\right)S_{45}+S_{12}S_{45}^{2}\right),$
(28) $\displaystyle\hat{\Delta}_{3}=$ $\displaystyle
m^{4}\left(-m_{5}^{2}+S_{34}+S_{45}\right){}^{2}+\left(m_{5}^{2}S_{12}+\left(S_{34}-S_{12}\right)S_{45}\right){}^{2}$
$\displaystyle+2m^{2}\big{(}m_{5}^{2}\left(-S_{12}\right)\left(m_{5}^{2}+S_{34}\right)+m_{5}^{2}\left(2S_{12}+S_{34}\right)S_{45}-\left(S_{12}+S_{34}\right)S_{45}^{2}$
$\displaystyle+\left(S_{12}-S_{34}\big{)}S_{34}S_{45}\right),$ (29)
$\displaystyle\hat{\Delta}_{4}=$ $\displaystyle
S_{12}\left(m^{4}S_{12}-2m_{5}^{2}m^{2}\left(S_{12}+2S_{15}-2S_{34}\right)+4m^{2}S_{15}\left(S_{12}+S_{15}-S_{34}\right)+m_{5}^{4}S_{12}\right),$
(30) $\displaystyle\hat{\Delta}_{5}=$ $\displaystyle m^{4}\,\hat{\Delta}_{6},$
(31) $\displaystyle\hat{\Delta}_{6}=$
$\displaystyle\left(m_{5}^{2}S_{23}-S_{23}S_{34}+\left(S_{34}-S_{15}\right)S_{45}\right){}^{2}+S_{12}^{2}\left(S_{15}-S_{23}\right){}^{2}$
$\displaystyle+2S_{12}\big{(}m_{5}^{2}S_{23}\left(S_{15}-S_{23}-2S_{34}\right)+S_{15}S_{23}S_{34}-S_{15}^{2}S_{45}+S_{15}\left(S_{23}+S_{34}\right)S_{45}$
$\displaystyle+S_{23}S_{34}\left(S_{45}-S_{23}\right)\big{)}$ (32)

In comparison with (3.2), we see the same two square roots associated with the
leading singularities of the massive three-point functions, $r_{2}$ and with
the Gram determinant of the external momenta, $r_{5}$, namely

$\displaystyle r_{2}=$
$\displaystyle\sqrt{\Delta_{1}}\,x=\sqrt{\hat{\Delta}_{1}}\,x$ (33)
$\displaystyle r_{5}=$
$\displaystyle\sqrt{\Delta_{2}}\,x^{2}=\sqrt{\hat{\Delta}_{6}}\,x^{2}$ (34)

We have however four more square roots,
$\sqrt{\hat{\Delta}_{2}},\,\sqrt{\hat{\Delta}_{3}},\,\sqrt{\hat{\Delta}_{4}},\,\sqrt{\hat{\Delta}_{5}}$,
which involve the internal mass $m^{2}$, and are not directly associated with
any leading singularities of any diagram from this family. We see therefore
that the introduction of an internal mass has a major impact on the complexity
of the resulting alphabet.

#### Final remarks and comparison with other methods

We have seen that introducing an internal mass can significantly increase the
complexity of the algebraic structure of an alphabet. It would be interesting
to find an explanation for the appearance of the additional square roots in
(3.2) which are not directly associated with a leading singularity. Normally,
one would expect that the square roots that appear in the alphabet also appear
in the definition of the pure basis elements, although in the case of family
$H$ the extra square roots do not appear in the pure basis definition.

On a more general note, whenever the Simplified Differential Equations
approach has been applied in conjunction with a pure basis Canko:2020ylt ;
Syrrakos:2020kba ; Canko:2020gqp , we have observed a reduction in the number
of letters that appear when compared with the alphabets that arise through the
usual method of differential equations, i.e. when one differentiates with
respect to all kinematic invariants. It would be interesting to see the
structure of the alphabets for the families studied here when one uses the
usual method of differential equations and whether this feature of alphabets
in $x$ with fewer letters still holds.

###  3.3 Families $D,F$

For these families we will obtain analytic expressions through the results of
the families $E$ and $G$. We will follow the procedure of taking the $x\to 1$
limit of our solution for a family with $n$ massive legs to obtain a pure
basis and analytic solution of a family with $n-1$ massive legs, as described
in detail in Canko:2020gqp .

For families $E$ and $G$ we exploit the shuffle properties of Goncharov
polylogarithms to write their solution (3.1) as an expansion in terms of
$\log(1-x)$ as follows

$\textbf{g}=\sum_{n\geq
0}\epsilon^{n}\sum_{i=0}^{n}\frac{1}{i!}\textbf{c}^{(n)}_{i}\log^{i}(1-x)$
(35)

with all $\textbf{c}^{(n)}_{i}$ being finite in the limit $x\to 1$. The next
step is to define the regular part of (35) at $x=1$

$\textbf{g}_{reg}=\sum_{n\geq 0}\epsilon^{n}\textbf{c}^{(n)}_{0}$ (36)

and after setting $x=1$ explicitly in (36) we may define the truncated part of
(35),

$\textbf{g}_{trunc}=\textbf{g}_{reg}(x=1)$ (37)

Having done that, we utilise the residue matrix that corresponds to the letter
$\\{1\\}$, $\textbf{M}_{2}$, and define the resummation matrix
$\tilde{\textbf{R}}$ as follows

$\tilde{\textbf{R}}=\tilde{\textbf{S}}e^{\epsilon\tilde{\textbf{D}}\log(1-x)}\tilde{\textbf{S}}^{-1}$
(38)

were $\tilde{\textbf{S}},\tilde{\textbf{D}}$ are constructed through the
Jordan decomposition of $\textbf{M}_{2}$, i.e.
$\textbf{M}_{2}=\tilde{\textbf{S}}\tilde{\textbf{D}}\tilde{\textbf{S}}^{-1}$.
The resummation matrix $\tilde{\textbf{R}}$ has terms of
$(1-x)^{a_{i}\epsilon}$, with $a_{i}$ being the eigenvalues of
$\textbf{M}_{2}$. After setting all terms $(1-x)^{a_{i}\epsilon}$ equal to
zero, we define the purely numerical matrix $\tilde{\textbf{R}}_{0}$.
Obtaining the $x\to 1$ limit of (3.1) amounts to acting with
$\tilde{\textbf{R}}_{0}$ on (37)

$\textbf{g}_{x\to 1}=\tilde{\textbf{R}}_{0}\textbf{g}_{trunc}$ (39)

Up to now we have calculated the $x\to 1$ limit for families $E$ and $G$. This
operation not only yields the result for a given family of Master integrals
with $n$ massive legs at a special limit, but also allows us to obtain results
for an independent family of Master integrals with $n-1$ massive legs. In the
case of family $E$ for example, taking the $x\to 1$ limit makes $q_{3}$ to
become massless through (3), thus yielding the kinematics for family $D$.
However, having the $x\to 1$ limit of $E$ means that we have explicit
solutions for the 18 pure basis elements of that family, whereas family $D$
has 16 basis elements. This means that out of the 18 basis elements of family
$E$ at this limit, we need to find the 16 of them that form the pure basis and
result for family $D$, i.e. at $x\to 1$ only 16 of the 18 basis elements of
family $E$ should remain linearly independent. The same reasoning holds for
family $G$ which has 16 basis elements while family $F$ has 15.

In order to find the pure bases for families $D$ and $F$ we can either use
Integration-By-Part identities to do the reduction or follow the approach
described in Canko:2020gqp . We shall use the latter method in the following.

In all cases that we have considered so far, $\tilde{\textbf{R}}_{0}$ is
always an idempotent matrix which means that, among others, it has the
following very useful property

$\tilde{\textbf{R}}_{0}^{2}=\tilde{\textbf{R}}_{0}$ (40)

Acting with $\tilde{\textbf{R}}_{0}$ on (39) and using (40) yields the
following relation

$\displaystyle\tilde{\textbf{R}}_{0}\textbf{g}_{x\to 1}$
$\displaystyle=\tilde{\textbf{R}}_{0}^{2}\textbf{g}_{trunc}$
$\displaystyle=\tilde{\textbf{R}}_{0}\textbf{g}_{trunc}$
$\displaystyle=\textbf{g}_{x\to 1}$ (41)

This relation, solved as an equation for each row, produces relations that
allow us to determine the linearly independent basis elements for families $D$
and $F$. Therefore, applying (3.3) to the pure basis of $E$ and $G$ yields in
an algorithmic way the pure bases for families $D$ and $F$.

###  3.4 On the choice of integral families

The basic rule for choosing which five-point family to consider is to have the
one-mass result Syrrakos:2020kba  as a starting point and add masses, with the
condition that their simplified differential equations in canonical form have
alphabets rational in $x$, when one uses parametrization (3). If the resulting
alphabet is not rational in $x$, then calculating the integral family through
the $x\to 1$ limit of another is considered. The exception to the above rule
is the five-point family with one internal mass and massless external legs,
which in the framework of the Simplified Differential Equations can only be
calculated through the $x\to 1$ limit of its corresponding one-mass family.

More specifically, if one tries to calculate family $D$ using (3) and deriving
differential equations in $x$, then the resulting alphabet is not rational in
$x$. Nevertheless, family $D$ can be expressed in terms of Goncharov
polylogarithms through the $x\to 1$ limit of family $E$. Apart from family $E$
there is another family with massless propagators and three massive legs, the
one which all three masses are adjacent. However, the alphabet of this family
is not rational in $x$, if one uses (3) to parametrize it.

Introducing an internal mass allows for many more families to be considered.
Apart from the ones presented in this paper, a family with one internal mass
and three massive legs (i.e. take family $H$ and regard $q_{1}$ as massive)
was considered, however its alphabet in $x$ is not rational using (3).

If the families with non-rational alphabets in $x$ can be rationalised using a
parametrization other than (3) remains an open question.

##  4 Validation

For all families computed in this paper we have made heavy use of the
Mathematica package PolyLogTools Duhr:2019tlz  for the manipulation of the
resulting Goncharov polylogarithms. As shown in (3.1), we provide explicit
results up to order $\mathcal{O}(\epsilon^{4})$. In Table 1 we provide an
analysis of our results for each family, regarding the number of Goncharov
polylogarithms that appear in each transcendental weight, where the weight is
counted as the number of $l_{i}$ indices of
$\mathcal{G}(l_{a},l_{b},\ldots;x)$. These numbers are obtained by gathering
all Goncharov polylogarithms that appear up to order
$\mathcal{O}(\epsilon^{4})$ in each integral family, and distinguishing them
according to their corresponding weight. For comparison, we perform the same
task for the top-sector basis element of each family.

A common feature of our results is that our solutions are dominated by the
number of weight-four polylogarithmic functions. Due to the universally
transcendental feature of our solutions, the $\mathcal{O}(\epsilon^{4})$ part
is expected to be the most cumbersome to calculate numerically, since weight-
four Goncharov polylogarithms take longer to calculate than lower-weight ones.
To avoid any misconceptions, it should be noted that each top-sector basis
element starts from $\mathcal{O}(\epsilon^{3})$, despite containing lower-
weight polylogarithms.

In appendix A we present explicit formulas for the weight-three part of the
pure top-sector basis elements of families $C$ and $H$ in order to give an
idea of the structure and length of the resulting expressions, as well as the
way the letters of the alphabets studied in subsection 3.2 are introduced in
the relevant solutions.

Family | W=1 | W=2 | W=3 | W=4 | Total
---|---|---|---|---|---
$C$ | 9 (0) | 54 (16) | 204 (106) | 605 (272) | 872 (394)
$E$ | 13 (0) | 87 (24) | 349 (172) | 1033 (432) | 1482 (628)
$G$ | 21 (4) | 163 (50) | 878 (329) | 2160 (884) | 3222 (1267)
$H$ | 19 (0) | 195 (42) | 1527 (616) | 5914 (2732) | 7655 (3390)
$D$ | 11 (0) | 83 (24) | 393 (192) | 1445 (656) | 1932 (872)
$F$ | 19 (4) | 151 (50) | 872 (349) | 2356 (1042) | 3398 (1445)
Table 1: Number of Goncharov polylogarithms entering the solution. Results for
the respective top-sector basis elements are in parenthesis.

Regarding the validation of our results, we have performed numerical checks of
our solution for each family against pySecDec Borowka:2017idc  for Euclidean
points. All Goncharov polylogarithms have been computed numerically using the
Ginsh command of PolyLogTools Duhr:2019tlz , as well as handyG Naterop:2019xaf
, which is a Fortran implementation of the algorithms developed in
Vollinga:2004sn . For all checks that we have performed we have found perfect
agreement. In Table 2 we also provide numerical results and timing for the
top-sector basis element for each corresponding family. We include timings
using handyG since we found that it is in general faster, although it is
restrictive in its precision compared to PolyLogTools.

Top-Sector | Time (sec) | Result
---|---|---
$C$ | 0.146897 | $-0.314547\epsilon^{4}-0.120811\epsilon^{3}$
$E$ | 0.248436 | $-0.0332408\epsilon^{4}-0.0215131\epsilon^{3}$
$G$ | 0.475048 | $-0.439003\epsilon^{4}-0.130267\epsilon^{3}$
$H$ | 1.89365 | $-0.0165223\epsilon^{4}-0.0192393\epsilon^{3}$
$D$ | 2.15734 | $-0.127286\epsilon^{4}-0.162439\epsilon^{3}$
$F$ | 0.730996 | $-0.528266\epsilon^{4}-0.33331\epsilon^{3}$
Table 2: Numerical computation of Goncharov polylogarithms using handyG with
double precision. The computations were performed on a 1,6 GHz Intel Core i5
laptop using a single CPU core.

Along with this paper we provide all of our results in the ancillary files.
More specifically, for the families $C,E,G$ and $H$, we provide the pure
basis, boundary terms in closed form, the alphabet, the Feynman integrals in
terms of which we have expressed the pure basis, the residue matrices for the
canonical differential equation, as well as the explicit result for each of
the aforementioned families in terms of Goncharov polylogarithms up to
transcendental weight four. For families $D$ and $F$ we provide the pure
basis, the alphabet and explicit results in terms of Goncharov polylogarithms
up to transcendental weight four. All ancillary files are in the form
familyname.m. We also provide for all families the Euclidean points that were
used to obtain the numerical results of Table 2 in the file EPoints.m.

##  5 Conclusions

The current frontier in the calculation of multiscale multiloop Feynman
integrals for $2\to 3$ scattering processes relevant to LHC searches lies at
two-loop five-point Feynman integrals with one off-shell leg and massless
internal lines. As of this writing, results for all planar two-loop five-point
Master integrals have been obtained using a numerical Abreu:2020jxa , as well
as an analytical approach Canko:2020ylt . Regarding the analytical results,
all planar Master integrals were expressed in terms of Goncharov
polylogarithms of up to transcendental weight four. These results, along with
analytic results for the relevant one-loop five-point Master integrals with
one off-shell leg Syrrakos:2020kba , were recently used to perform the first
fully analytic calculation of a two-loop scattering amplitude for $Wb\bar{b}$
production Badger:2021nhg . Furthermore, using a new method for calculating
Master integrals, the authors of Papadopoulos:2019iam  have computed
numerically one of the non-planar two-loop five-point families. More recently,
some of the authors of Abreu:2020jxa  presented the three hexabox topologies
in Abreu:2021smk  using the same techniques as in Abreu:2020jxa .

Looking ahead, at some point we will have to consider more complicated Feynman
integrals, involving more massive external particles and/or massive
propagators. One of the expected challenges when considering such integrals is
the introduction of many square roots in the alphabet of the resulting
canonical differential equations. It remains a non-trivial exercise to find a
universal way to handle these roots and achieve a result in terms of Goncharov
polylogarithms, however several ideas have been put forward in recent times
Heller:2019gkq ; Besier:2018jen ; Besier:2019kco ; Bonetti:2020hqh . One
should also keep in mind that even if a so-called dlog form of the
differential equations is achieved, it does not guarantee that its solution
will be in terms of Goncharov polylogarithms Brown:2020rda .

In order to get a glimpse of the complexities that lie beyond the frontier of
five-point scattering involving one off-shell leg and massless internal lines,
in this paper we studied families of one-loop five-point Feynman integrals
with two and three massive external legs and massless propagators, as well as
one-loop five-point families with one massive internal line and up to two
massive external legs.

We used the Simplified Differential Equations approach for the construction of
canonical differential equations for pure bases of the families of Figure 1
and as a special limit we obtained results for the families of Figure 2. As it
turned out, the parametrization (3) was enough to rationalise all square roots
introduced in the alphabet of families $C,E,G,H$ of Figure 1. For these
families we were also able to obtain boundary terms for the canonical
differential equations in closed form, allowing us to trivially derive
solutions for these families in terms of Goncharov polylogarithms of arbitrary
transcendental weight.

For families $D,F$ of Figure 2 we obtained analytic results through a special
limit of our solutions for families $E$ and $G$ respectively. For all families
studied in this paper we provide explicit results in terms of Goncharov
polylogarithms of up to transcendental weight four.

Regarding the structure of the resulting alphabets in $x$, we saw that when
one internal mass is introduced, square roots involving this mass arise, which
are not present in the definition of the pure basis. Further study of these
alphabets is required to pin-point the origin of these additional square
roots, which we leave for future work. It is also interesting to explore in
the future the structure of these alphabets when one employs the standard
method of differential equations, i.e. differentiating with respect to all
kinematic variables. A comparison between the two approaches might further
elucidate the effectiveness of the Simplified Differential Equations approach
in providing solutions to multiscale Feynman integrals in terms of Goncharov
polylogarithms, as well as provide an idea of whether the representation of
these integrals in terms of polylogarithmic functions is the best one for
phenomenological applications.

###### Acknowledgements.

The author acknowledges fruitful discussions with Costas G. Papadopoulos,
Dhimiter D. Canko and Lorenzo Tancredi.

##  Appendix A Explicit results at weight three

In this appendix we provide explicit results for families $C$ and $H$ for the
weight-three part of each top-sector pure basis element. Assuming that each
top-sector basis element is expressed in the following manner,

$g_{i}^{\text{family}}=\sum_{w=3}^{4}\epsilon^{w}\,\tilde{g}_{i,w}^{\text{family}}$
(42)

we give the explicit expression for the $\tilde{g}_{i,3}^{\text{family}}$
part. We introduce the following shorthand notations for brevity,

$\displaystyle\mathcal{G}_{a,b,\ldots}=\mathcal{G}(l_{a},l_{b},\ldots;x)$ (43)
$\displaystyle
L_{1}=\log\left(-S_{12}\right),L_{2}=\log\left(-S_{45}\right),L_{3}=\log\left(-m_{5}^{2}\right),$
(44) $\displaystyle
L_{4}=\log\left(m^{2}\right),L_{5}=\log\left(\frac{m^{2}-S_{45}}{m^{2}}\right),L_{6}=\log\left(\frac{m^{2}-m_{5}^{2}}{m^{2}}\right)$
(45)

###  A.1 Top sector of family C

In the following formula note that letters $l_{i},\,i=\\{6,7,13,14\\}$ contain
square roots in the underline kinematic variables $S_{ij}$, as shown
explicitly in (3.2).

$\displaystyle\tilde{g}_{15,\,3}^{C}=$
$\displaystyle\frac{1}{4}L_{1}\mathcal{G}_{13,2}-\frac{1}{4}L_{2}\mathcal{G}_{13,2}+\frac{1}{4}L_{1}\mathcal{G}_{13,3}-\frac{1}{4}L_{3}\mathcal{G}_{13,3}-\frac{1}{4}L_{1}\mathcal{G}_{13,6}+\frac{1}{4}L_{3}\mathcal{G}_{13,6}-\frac{1}{4}L_{1}\mathcal{G}_{13,7}$
$\displaystyle+\frac{1}{4}L_{3}\mathcal{G}_{13,7}-\frac{1}{4}L_{2}\mathcal{G}_{13,8}+\frac{1}{4}L_{3}\mathcal{G}_{13,8}+\frac{1}{4}L_{2}\mathcal{G}_{13,9}-\frac{1}{4}L_{3}\mathcal{G}_{13,9}+\frac{1}{4}L_{1}\mathcal{G}_{13,10}-\frac{1}{4}L_{2}\mathcal{G}_{13,10}$
$\displaystyle-\frac{1}{4}L_{1}\mathcal{G}_{13,11}+\frac{1}{4}L_{2}\mathcal{G}_{13,11}-\frac{1}{4}L_{1}\mathcal{G}_{14,2}+\frac{1}{4}L_{2}\mathcal{G}_{14,2}-\frac{1}{4}L_{1}\mathcal{G}_{14,3}+\frac{1}{4}L_{3}\mathcal{G}_{14,3}+\frac{1}{4}L_{1}\mathcal{G}_{14,6}$
$\displaystyle-\frac{1}{4}L_{3}\mathcal{G}_{14,6}+\frac{1}{4}L_{1}\mathcal{G}_{14,7}-\frac{1}{4}L_{3}\mathcal{G}_{14,7}+\frac{1}{4}L_{2}\mathcal{G}_{14,8}-\frac{1}{4}L_{3}\mathcal{G}_{14,8}-\frac{1}{4}L_{2}\mathcal{G}_{14,9}+\frac{1}{4}L_{3}\mathcal{G}_{14,9}$
$\displaystyle-\frac{1}{4}L_{1}\mathcal{G}_{14,10}+\frac{1}{4}L_{2}\mathcal{G}_{14,10}+\frac{1}{4}L_{1}\mathcal{G}_{14,11}-\frac{1}{4}L_{2}\mathcal{G}_{14,11}+\frac{1}{2}\mathcal{G}_{13,2,1}+\frac{1}{2}\mathcal{G}_{13,3,1}+\frac{1}{4}\mathcal{G}_{13,3,3}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{13,3,6}-\frac{1}{4}\mathcal{G}_{13,3,7}-\frac{1}{4}\mathcal{G}_{13,4,2}-\frac{1}{4}\mathcal{G}_{13,4,10}+\frac{1}{4}\mathcal{G}_{13,4,11}-\frac{1}{4}\mathcal{G}_{13,5,3}+\frac{1}{4}\mathcal{G}_{13,5,6}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{13,5,7}-\frac{1}{2}\mathcal{G}_{13,6,1}-\frac{1}{2}\mathcal{G}_{13,7,1}+\frac{1}{4}\mathcal{G}_{13,8,3}-\frac{1}{4}\mathcal{G}_{13,8,11}+\frac{1}{4}\mathcal{G}_{13,9,2}-\frac{1}{4}\mathcal{G}_{13,9,6}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{13,9,7}+\frac{1}{4}\mathcal{G}_{13,9,10}+\frac{1}{2}\mathcal{G}_{13,10,1}-\frac{1}{2}\mathcal{G}_{13,11,1}+\frac{1}{4}\mathcal{G}_{13,11,2}+\frac{1}{4}\mathcal{G}_{13,11,10}-\frac{1}{4}\mathcal{G}_{13,11,11}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{13,12,2}+\frac{1}{4}\mathcal{G}_{13,12,6}+\frac{1}{4}\mathcal{G}_{13,12,7}-\frac{1}{4}\mathcal{G}_{13,12,10}-\frac{1}{2}\mathcal{G}_{14,2,1}-\frac{1}{2}\mathcal{G}_{14,3,1}-\frac{1}{4}\mathcal{G}_{14,3,3}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{14,3,6}+\frac{1}{4}\mathcal{G}_{14,3,7}+\frac{1}{4}\mathcal{G}_{14,4,2}+\frac{1}{4}\mathcal{G}_{14,4,10}-\frac{1}{4}\mathcal{G}_{14,4,11}+\frac{1}{4}\mathcal{G}_{14,5,3}-\frac{1}{4}\mathcal{G}_{14,5,6}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{14,5,7}+\frac{1}{2}\mathcal{G}_{14,6,1}+\frac{1}{2}\mathcal{G}_{14,7,1}-\frac{1}{4}\mathcal{G}_{14,8,3}+\frac{1}{4}\mathcal{G}_{14,8,11}-\frac{1}{4}\mathcal{G}_{14,9,2}+\frac{1}{4}\mathcal{G}_{14,9,6}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{14,9,7}-\frac{1}{4}\mathcal{G}_{14,9,10}-\frac{1}{2}\mathcal{G}_{14,10,1}+\frac{1}{2}\mathcal{G}_{14,11,1}-\frac{1}{4}\mathcal{G}_{14,11,2}-\frac{1}{4}\mathcal{G}_{14,11,10}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{14,11,11}+\frac{1}{4}\mathcal{G}_{14,12,2}-\frac{1}{4}\mathcal{G}_{14,12,6}-\frac{1}{4}\mathcal{G}_{14,12,7}+\frac{1}{4}\mathcal{G}_{14,12,10}$
(46)

###  A.2 Top sector of family H

In the following formula note that letters
$l_{i},\,i=\\{7,8,9,10,16,17,18,19,25,26,27,28,29,30\\}$ contain square roots
in the underline kinematic variables $S_{ij}$, as shown explicitly in (3.2).

$\displaystyle g_{18,\,3}^{H}=$
$\displaystyle\frac{1}{4}L_{2}\mathcal{G}_{27,3}-\frac{1}{4}L_{4}\mathcal{G}_{27,3}-\frac{1}{4}L_{5}\mathcal{G}_{27,3}-\frac{1}{4}L_{3}\mathcal{G}_{27,9}+\frac{1}{4}L_{4}\mathcal{G}_{27,9}+\frac{1}{4}L_{6}\mathcal{G}_{27,9}-\frac{1}{4}L_{3}\mathcal{G}_{27,10}$
$\displaystyle+\frac{1}{4}L_{4}\mathcal{G}_{27,10}+\frac{1}{4}L_{6}\mathcal{G}_{27,10}-\frac{1}{4}L_{2}\mathcal{G}_{27,11}+\frac{1}{4}L_{4}\mathcal{G}_{27,11}+\frac{1}{4}L_{5}\mathcal{G}_{27,11}+\frac{1}{4}L_{2}\mathcal{G}_{27,13}-\frac{1}{4}L_{4}\mathcal{G}_{27,13}$
$\displaystyle-\frac{1}{4}L_{5}\mathcal{G}_{27,13}-\frac{1}{4}L_{2}\mathcal{G}_{27,18}+\frac{1}{4}L_{3}\mathcal{G}_{27,18}+\frac{1}{4}L_{5}\mathcal{G}_{27,18}-\frac{1}{4}L_{6}\mathcal{G}_{27,18}-\frac{1}{4}L_{2}\mathcal{G}_{27,19}+\frac{1}{4}L_{3}\mathcal{G}_{27,19}$
$\displaystyle+\frac{1}{4}L_{5}\mathcal{G}_{27,19}-\frac{1}{4}L_{6}\mathcal{G}_{27,19}+\frac{1}{4}L_{3}\mathcal{G}_{27,20}-\frac{1}{4}L_{4}\mathcal{G}_{27,20}-\frac{1}{4}L_{6}\mathcal{G}_{27,20}+\frac{1}{4}L_{2}\mathcal{G}_{27,24}-\frac{1}{4}L_{3}\mathcal{G}_{27,24}$
$\displaystyle-\frac{1}{4}L_{5}\mathcal{G}_{27,24}+\frac{1}{4}L_{6}\mathcal{G}_{27,24}-\frac{1}{4}L_{2}\mathcal{G}_{28,3}+\frac{1}{4}L_{4}\mathcal{G}_{28,3}+\frac{1}{4}L_{5}\mathcal{G}_{28,3}+\frac{1}{4}L_{3}\mathcal{G}_{28,9}-\frac{1}{4}L_{4}\mathcal{G}_{28,9}$
$\displaystyle-\frac{1}{4}L_{6}\mathcal{G}_{28,9}+\frac{1}{4}L_{3}\mathcal{G}_{28,10}-\frac{1}{4}L_{4}\mathcal{G}_{28,10}-\frac{1}{4}L_{6}\mathcal{G}_{28,10}+\frac{1}{4}L_{2}\mathcal{G}_{28,11}-\frac{1}{4}L_{4}\mathcal{G}_{28,11}-\frac{1}{4}L_{5}\mathcal{G}_{28,11}$
$\displaystyle-\frac{1}{4}L_{2}\mathcal{G}_{28,13}+\frac{1}{4}L_{4}\mathcal{G}_{28,13}+\frac{1}{4}L_{5}\mathcal{G}_{28,13}+\frac{1}{4}L_{2}\mathcal{G}_{28,18}-\frac{1}{4}L_{3}\mathcal{G}_{28,18}-\frac{1}{4}L_{5}\mathcal{G}_{28,18}+\frac{1}{4}L_{6}\mathcal{G}_{28,18}$
$\displaystyle+\frac{1}{4}L_{2}\mathcal{G}_{28,19}-\frac{1}{4}L_{3}\mathcal{G}_{28,19}-\frac{1}{4}L_{5}\mathcal{G}_{28,19}+\frac{1}{4}L_{6}\mathcal{G}_{28,19}-\frac{1}{4}L_{3}\mathcal{G}_{28,20}+\frac{1}{4}L_{4}\mathcal{G}_{28,20}+\frac{1}{4}L_{6}\mathcal{G}_{28,20}$
$\displaystyle-\frac{1}{4}L_{2}\mathcal{G}_{28,24}+\frac{1}{4}L_{3}\mathcal{G}_{28,24}+\frac{1}{4}L_{5}\mathcal{G}_{28,24}-\frac{1}{4}L_{6}\mathcal{G}_{28,24}+\frac{1}{4}L_{5}\mathcal{G}_{29,2}-\frac{1}{4}L_{6}\mathcal{G}_{29,7}-\frac{1}{4}L_{6}\mathcal{G}_{29,8}$
$\displaystyle-\frac{1}{4}L_{2}\mathcal{G}_{29,12}+\frac{1}{4}L_{3}\mathcal{G}_{29,12}+\frac{1}{4}L_{5}\mathcal{G}_{29,14}-\frac{1}{4}L_{5}\mathcal{G}_{29,15}+\frac{1}{4}L_{2}\mathcal{G}_{29,18}-\frac{1}{4}L_{3}\mathcal{G}_{29,18}-\frac{1}{4}L_{5}\mathcal{G}_{29,18}$
$\displaystyle+\frac{1}{4}L_{6}\mathcal{G}_{29,18}+\frac{1}{4}L_{2}\mathcal{G}_{29,19}-\frac{1}{4}L_{3}\mathcal{G}_{29,19}-\frac{1}{4}L_{5}\mathcal{G}_{29,19}+\frac{1}{4}L_{6}\mathcal{G}_{29,19}+\frac{1}{4}L_{6}\mathcal{G}_{29,21}+\frac{1}{4}L_{2}\mathcal{G}_{29,22}$
$\displaystyle-\frac{1}{4}L_{3}\mathcal{G}_{29,22}-\frac{1}{4}L_{2}\mathcal{G}_{29,24}+\frac{1}{4}L_{3}\mathcal{G}_{29,24}+\frac{1}{4}L_{5}\mathcal{G}_{29,24}-\frac{1}{4}L_{6}\mathcal{G}_{29,24}-\frac{1}{4}L_{5}\mathcal{G}_{30,2}+\frac{1}{4}L_{6}\mathcal{G}_{30,7}$
$\displaystyle+\frac{1}{4}L_{6}\mathcal{G}_{30,8}+\frac{1}{4}L_{2}\mathcal{G}_{30,12}-\frac{1}{4}L_{3}\mathcal{G}_{30,12}-\frac{1}{4}L_{5}\mathcal{G}_{30,14}+\frac{1}{4}L_{5}\mathcal{G}_{30,15}-\frac{1}{4}L_{2}\mathcal{G}_{30,18}+\frac{1}{4}L_{3}\mathcal{G}_{30,18}$
$\displaystyle+\frac{1}{4}L_{5}\mathcal{G}_{30,18}-\frac{1}{4}L_{6}\mathcal{G}_{30,18}-\frac{1}{4}L_{2}\mathcal{G}_{30,19}+\frac{1}{4}L_{3}\mathcal{G}_{30,19}+\frac{1}{4}L_{5}\mathcal{G}_{30,19}-\frac{1}{4}L_{6}\mathcal{G}_{30,19}-\frac{1}{4}L_{6}\mathcal{G}_{30,21}$
$\displaystyle-\frac{1}{4}L_{2}\mathcal{G}_{30,22}+\frac{1}{4}L_{3}\mathcal{G}_{30,22}+\frac{1}{4}L_{2}\mathcal{G}_{30,24}-\frac{1}{4}L_{3}\mathcal{G}_{30,24}-\frac{1}{4}L_{5}\mathcal{G}_{30,24}+\frac{1}{4}L_{6}\mathcal{G}_{30,24}+\frac{1}{4}\mathcal{G}_{27,3,2}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{27,3,4}-\frac{1}{4}\mathcal{G}_{27,3,5}+\frac{1}{4}\mathcal{G}_{27,3,14}+\frac{1}{4}\mathcal{G}_{27,9,4}+\frac{1}{4}\mathcal{G}_{27,9,5}-\frac{1}{4}\mathcal{G}_{27,9,7}-\frac{1}{4}\mathcal{G}_{27,9,8}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{27,10,4}+\frac{1}{4}\mathcal{G}_{27,10,5}-\frac{1}{4}\mathcal{G}_{27,10,7}-\frac{1}{4}\mathcal{G}_{27,10,8}-\frac{1}{4}\mathcal{G}_{27,11,15}+\frac{1}{4}\mathcal{G}_{27,13,2}-\frac{1}{4}\mathcal{G}_{27,13,4}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{27,13,5}+\frac{1}{4}\mathcal{G}_{27,13,14}-\frac{1}{4}\mathcal{G}_{27,16,2}+\frac{1}{4}\mathcal{G}_{27,16,4}+\frac{1}{4}\mathcal{G}_{27,16,5}-\frac{1}{4}\mathcal{G}_{27,16,14}+\frac{1}{4}\mathcal{G}_{27,16,15}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{27,17,2}+\frac{1}{4}\mathcal{G}_{27,17,4}+\frac{1}{4}\mathcal{G}_{27,17,5}-\frac{1}{4}\mathcal{G}_{27,17,14}+\frac{1}{4}\mathcal{G}_{27,17,15}-\frac{1}{4}\mathcal{G}_{27,18,2}+\frac{1}{4}\mathcal{G}_{27,18,7}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{27,18,8}-\frac{1}{4}\mathcal{G}_{27,18,14}-\frac{1}{4}\mathcal{G}_{27,19,2}+\frac{1}{4}\mathcal{G}_{27,19,7}+\frac{1}{4}\mathcal{G}_{27,19,8}-\frac{1}{4}\mathcal{G}_{27,19,14}+\frac{1}{4}\mathcal{G}_{27,20,21}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{27,24,15}-\frac{1}{4}\mathcal{G}_{27,24,21}-\frac{1}{4}\mathcal{G}_{27,25,4}-\frac{1}{4}\mathcal{G}_{27,25,5}+\frac{1}{4}\mathcal{G}_{27,25,7}+\frac{1}{4}\mathcal{G}_{27,25,8}-\frac{1}{4}\mathcal{G}_{27,25,21}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{27,26,4}-\frac{1}{4}\mathcal{G}_{27,26,5}+\frac{1}{4}\mathcal{G}_{27,26,7}+\frac{1}{4}\mathcal{G}_{27,26,8}-\frac{1}{4}\mathcal{G}_{27,26,21}+\frac{1}{4}\mathcal{G}_{27,29,2}-\frac{1}{4}\mathcal{G}_{27,29,7}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{27,29,8}+\frac{1}{4}\mathcal{G}_{27,29,14}-\frac{1}{4}\mathcal{G}_{27,29,15}+\frac{1}{4}\mathcal{G}_{27,29,21}+\frac{1}{4}\mathcal{G}_{27,30,2}-\frac{1}{4}\mathcal{G}_{27,30,7}-\frac{1}{4}\mathcal{G}_{27,30,8}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{27,30,14}-\frac{1}{4}\mathcal{G}_{27,30,15}+\frac{1}{4}\mathcal{G}_{27,30,21}-\frac{1}{4}\mathcal{G}_{28,3,2}+\frac{1}{4}\mathcal{G}_{28,3,4}+\frac{1}{4}\mathcal{G}_{28,3,5}-\frac{1}{4}\mathcal{G}_{28,3,14}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{28,9,4}-\frac{1}{4}\mathcal{G}_{28,9,5}+\frac{1}{4}\mathcal{G}_{28,9,7}+\frac{1}{4}\mathcal{G}_{28,9,8}-\frac{1}{4}\mathcal{G}_{28,10,4}-\frac{1}{4}\mathcal{G}_{28,10,5}+\frac{1}{4}\mathcal{G}_{28,10,7}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{28,10,8}+\frac{1}{4}\mathcal{G}_{28,11,15}-\frac{1}{4}\mathcal{G}_{28,13,2}+\frac{1}{4}\mathcal{G}_{28,13,4}+\frac{1}{4}\mathcal{G}_{28,13,5}-\frac{1}{4}\mathcal{G}_{28,13,14}+\frac{1}{4}\mathcal{G}_{28,16,2}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{28,16,4}-\frac{1}{4}\mathcal{G}_{28,16,5}+\frac{1}{4}\mathcal{G}_{28,16,14}-\frac{1}{4}\mathcal{G}_{28,16,15}+\frac{1}{4}\mathcal{G}_{28,17,2}-\frac{1}{4}\mathcal{G}_{28,17,4}-\frac{1}{4}\mathcal{G}_{28,17,5}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{28,17,14}-\frac{1}{4}\mathcal{G}_{28,17,15}+\frac{1}{4}\mathcal{G}_{28,18,2}-\frac{1}{4}\mathcal{G}_{28,18,7}-\frac{1}{4}\mathcal{G}_{28,18,8}+\frac{1}{4}\mathcal{G}_{28,18,14}+\frac{1}{4}\mathcal{G}_{28,19,2}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{28,19,7}-\frac{1}{4}\mathcal{G}_{28,19,8}+\frac{1}{4}\mathcal{G}_{28,19,14}-\frac{1}{4}\mathcal{G}_{28,20,21}-\frac{1}{4}\mathcal{G}_{28,24,15}+\frac{1}{4}\mathcal{G}_{28,24,21}+\frac{1}{4}\mathcal{G}_{28,25,4}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{28,25,5}-\frac{1}{4}\mathcal{G}_{28,25,7}-\frac{1}{4}\mathcal{G}_{28,25,8}+\frac{1}{4}\mathcal{G}_{28,25,21}+\frac{1}{4}\mathcal{G}_{28,26,4}+\frac{1}{4}\mathcal{G}_{28,26,5}-\frac{1}{4}\mathcal{G}_{28,26,7}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{28,26,8}+\frac{1}{4}\mathcal{G}_{28,26,21}-\frac{1}{4}\mathcal{G}_{28,29,2}+\frac{1}{4}\mathcal{G}_{28,29,7}+\frac{1}{4}\mathcal{G}_{28,29,8}-\frac{1}{4}\mathcal{G}_{28,29,14}+\frac{1}{4}\mathcal{G}_{28,29,15}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{28,29,21}-\frac{1}{4}\mathcal{G}_{28,30,2}+\frac{1}{4}\mathcal{G}_{28,30,7}+\frac{1}{4}\mathcal{G}_{28,30,8}-\frac{1}{4}\mathcal{G}_{28,30,14}+\frac{1}{4}\mathcal{G}_{28,30,15}-\frac{1}{4}\mathcal{G}_{28,30,21}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{29,2,4}-\frac{1}{4}\mathcal{G}_{29,2,5}+\frac{1}{4}\mathcal{G}_{29,6,2}+\frac{1}{4}\mathcal{G}_{29,6,14}-\frac{1}{4}\mathcal{G}_{29,6,15}+\frac{1}{4}\mathcal{G}_{29,7,4}+\frac{1}{4}\mathcal{G}_{29,7,5}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{29,8,4}+\frac{1}{4}\mathcal{G}_{29,8,5}-\frac{1}{4}\mathcal{G}_{29,12,2}+\frac{1}{4}\mathcal{G}_{29,12,7}+\frac{1}{4}\mathcal{G}_{29,12,8}-\frac{1}{4}\mathcal{G}_{29,12,14}-\frac{1}{4}\mathcal{G}_{29,14,4}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{29,14,5}-\frac{1}{4}\mathcal{G}_{29,16,2}+\frac{1}{4}\mathcal{G}_{29,16,4}+\frac{1}{4}\mathcal{G}_{29,16,5}-\frac{1}{4}\mathcal{G}_{29,16,14}+\frac{1}{4}\mathcal{G}_{29,16,15}-\frac{1}{4}\mathcal{G}_{29,17,2}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{29,17,4}+\frac{1}{4}\mathcal{G}_{29,17,5}-\frac{1}{4}\mathcal{G}_{29,17,14}+\frac{1}{4}\mathcal{G}_{29,17,15}+\frac{1}{4}\mathcal{G}_{29,18,2}-\frac{1}{4}\mathcal{G}_{29,18,7}-\frac{1}{4}\mathcal{G}_{29,18,8}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{29,18,14}+\frac{1}{4}\mathcal{G}_{29,19,2}-\frac{1}{4}\mathcal{G}_{29,19,7}-\frac{1}{4}\mathcal{G}_{29,19,8}+\frac{1}{4}\mathcal{G}_{29,19,14}+\frac{1}{4}\mathcal{G}_{29,22,15}-\frac{1}{4}\mathcal{G}_{29,22,21}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{29,23,7}-\frac{1}{4}\mathcal{G}_{29,23,8}+\frac{1}{4}\mathcal{G}_{29,23,21}-\frac{1}{4}\mathcal{G}_{29,24,15}+\frac{1}{4}\mathcal{G}_{29,24,21}-\frac{1}{4}\mathcal{G}_{29,25,4}-\frac{1}{4}\mathcal{G}_{29,25,5}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{29,25,7}+\frac{1}{4}\mathcal{G}_{29,25,8}-\frac{1}{4}\mathcal{G}_{29,25,21}-\frac{1}{4}\mathcal{G}_{29,26,4}-\frac{1}{4}\mathcal{G}_{29,26,5}+\frac{1}{4}\mathcal{G}_{29,26,7}+\frac{1}{4}\mathcal{G}_{29,26,8}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{29,26,21}+\frac{1}{4}\mathcal{G}_{30,2,4}+\frac{1}{4}\mathcal{G}_{30,2,5}-\frac{1}{4}\mathcal{G}_{30,6,2}-\frac{1}{4}\mathcal{G}_{30,6,14}+\frac{1}{4}\mathcal{G}_{30,6,15}-\frac{1}{4}\mathcal{G}_{30,7,4}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{30,7,5}-\frac{1}{4}\mathcal{G}_{30,8,4}-\frac{1}{4}\mathcal{G}_{30,8,5}+\frac{1}{4}\mathcal{G}_{30,12,2}-\frac{1}{4}\mathcal{G}_{30,12,7}-\frac{1}{4}\mathcal{G}_{30,12,8}+\frac{1}{4}\mathcal{G}_{30,12,14}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{30,14,4}+\frac{1}{4}\mathcal{G}_{30,14,5}+\frac{1}{4}\mathcal{G}_{30,16,2}-\frac{1}{4}\mathcal{G}_{30,16,4}-\frac{1}{4}\mathcal{G}_{30,16,5}+\frac{1}{4}\mathcal{G}_{30,16,14}-\frac{1}{4}\mathcal{G}_{30,16,15}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{30,17,2}-\frac{1}{4}\mathcal{G}_{30,17,4}-\frac{1}{4}\mathcal{G}_{30,17,5}+\frac{1}{4}\mathcal{G}_{30,17,14}-\frac{1}{4}\mathcal{G}_{30,17,15}-\frac{1}{4}\mathcal{G}_{30,18,2}+\frac{1}{4}\mathcal{G}_{30,18,7}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{30,18,8}-\frac{1}{4}\mathcal{G}_{30,18,14}-\frac{1}{4}\mathcal{G}_{30,19,2}+\frac{1}{4}\mathcal{G}_{30,19,7}+\frac{1}{4}\mathcal{G}_{30,19,8}-\frac{1}{4}\mathcal{G}_{30,19,14}-\frac{1}{4}\mathcal{G}_{30,22,15}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{30,22,21}+\frac{1}{4}\mathcal{G}_{30,23,7}+\frac{1}{4}\mathcal{G}_{30,23,8}-\frac{1}{4}\mathcal{G}_{30,23,21}+\frac{1}{4}\mathcal{G}_{30,24,15}-\frac{1}{4}\mathcal{G}_{30,24,21}+\frac{1}{4}\mathcal{G}_{30,25,4}$
$\displaystyle+\frac{1}{4}\mathcal{G}_{30,25,5}-\frac{1}{4}\mathcal{G}_{30,25,7}-\frac{1}{4}\mathcal{G}_{30,25,8}+\frac{1}{4}\mathcal{G}_{30,25,21}+\frac{1}{4}\mathcal{G}_{30,26,4}+\frac{1}{4}\mathcal{G}_{30,26,5}-\frac{1}{4}\mathcal{G}_{30,26,7}$
$\displaystyle-\frac{1}{4}\mathcal{G}_{30,26,8}+\frac{1}{4}\mathcal{G}_{30,26,21}$
(47)

## References

  * (1) G. Heinrich, “Collider Physics at the Precision Frontier,” Phys. Rept. 922 (2021), 1-69 doi:10.1016/j.physrep.2021.03.006 [arXiv:2009.00516 [hep-ph]].
  * (2) S. Amoroso, P. Azzurri, J. Bendavid, E. Bothmann, D. Britzger, H. Brooks, A. Buckley, M. Calvetti, X. Chen and M. Chiesa, et al. “Les Houches 2019: Physics at TeV Colliders: Standard Model Working Group Report,” [arXiv:2003.01700 [hep-ph]].
  * (3) K. G. Chetyrkin and F. V. Tkachov, “Integration by Parts: The Algorithm to Calculate beta Functions in 4 Loops,” Nucl. Phys. B 192 (1981), 159-204 doi:10.1016/0550-3213(81)90199-1.
  * (4) J. M. Henn, “Lectures on differential equations for Feynman integrals,” J. Phys. A 48 (2015), 153001 doi:10.1088/1751-8113/48/15/153001 [arXiv:1412.2296 [hep-ph]].
  * (5) A. V. Kotikov, “Differential equations method: New technique for massive Feynman diagrams calculation,” Phys. Lett. B 254 (1991), 158-164.
  * (6) A. V. Kotikov, “Differential equations method: The Calculation of vertex type Feynman diagrams,” Phys. Lett. B 259 (1991), 314-322.
  * (7) A. V. Kotikov, “Differential equation method: The Calculation of N point Feynman diagrams,” Phys. Lett. B 267 (1991), 123-127 [Errattum: Phys. Lett. B 295 (1992), 409].
  * (8) T. Gehrmann and E. Remiddi, “Differential equations for two loop four point functions,” Nucl. Phys. B 580 (2000), 485-518 [arXiv:hep-ph/9912329 [hep-ph]].
  * (9) J. M. Henn, “Multiloop integrals in dimensional regularization made simple,” Phys. Rev. Lett. 110 (2013), 251601 doi:10.1103/PhysRevLett.110.251601 [arXiv:1304.1806 [hep-th]].
  * (10) A. V. Kotikov, “Differential equations and Feynman integrals,” [arXiv:2102.07424 [hep-ph]].
  * (11) A. B. Goncharov, “Multiple polylogarithms, cyclotomy and modular complexes,” Math. Res. Lett. 5 (1998), 497-516 doi:10.4310/MRL.1998.v5.n4.a7 [arXiv:1105.2076 [math.AG]].
  * (12) C. Duhr, H. Gangl and J. R. Rhodes, “From polygons and symbols to polylogarithmic functions,” JHEP 10 (2012), 075 doi:10.1007/JHEP10(2012)075 [arXiv:1110.0458 [math-ph]].
  * (13) C. Duhr, “Hopf algebras, coproducts and symbols: an application to Higgs boson amplitudes,” JHEP 08 (2012), 043 doi:10.1007/JHEP08(2012)043 [arXiv:1203.0454 [hep-ph]].
  * (14) C. Duhr, “Mathematical aspects of scattering amplitudes,” doi:10.1142/9789814678766_0010 [arXiv:1411.7538 [hep-ph]].
  * (15) E. Remiddi and L. Tancredi, “An Elliptic Generalization of Multiple Polylogarithms,” Nucl. Phys. B 925 (2017), 212-251 doi:10.1016/j.nuclphysb.2017.10.007 [arXiv:1709.03622 [hep-ph]].
  * (16) J. Broedel, C. Duhr, F. Dulat and L. Tancredi, “Elliptic polylogarithms and iterated integrals on elliptic curves. Part I: general formalism,” JHEP 05 (2018), 093 doi:10.1007/JHEP05(2018)093 [arXiv:1712.07089 [hep-th]].
  * (17) J. Broedel, C. Duhr, F. Dulat and L. Tancredi, “Elliptic polylogarithms and iterated integrals on elliptic curves II: an application to the sunrise integral,” Phys. Rev. D 97 (2018) no.11, 116009 doi:10.1103/PhysRevD.97.116009 [arXiv:1712.07095 [hep-ph]].
  * (18) J. Broedel, C. Duhr, F. Dulat, B. Penante and L. Tancredi, “Elliptic symbol calculus: from elliptic polylogarithms to iterated integrals of Eisenstein series,” JHEP 08 (2018), 014 doi:10.1007/JHEP08(2018)014 [arXiv:1803.10256 [hep-th]].
  * (19) J. Broedel, C. Duhr, F. Dulat, B. Penante and L. Tancredi, “Elliptic Feynman integrals and pure functions,” JHEP 01 (2019), 023 doi:10.1007/JHEP01(2019)023 [arXiv:1809.10698 [hep-th]].
  * (20) J. Broedel, C. Duhr, F. Dulat, B. Penante and L. Tancredi, “Elliptic polylogarithms and Feynman parameter integrals,” JHEP 05 (2019), 120 doi:10.1007/JHEP05(2019)120 [arXiv:1902.09971 [hep-ph]].
  * (21) C. Duhr and L. Tancredi, “Algorithms and tools for iterated Eisenstein integrals,” JHEP 02 (2020), 105 doi:10.1007/JHEP02(2020)105 [arXiv:1912.00077 [hep-th]].
  * (22) C. G. Papadopoulos, D. Tommasini and C. Wever, “The Pentabox Master Integrals with the Simplified Differential Equations approach,” JHEP 04 (2016), 078 doi:10.1007/JHEP04(2016)078 [arXiv:1511.09404 [hep-ph]].
  * (23) T. Gehrmann, J. M. Henn and N. A. Lo Presti, “Analytic form of the two-loop planar five-gluon all-plus-helicity amplitude in QCD,” Phys. Rev. Lett. 116 (2016) no.6, 062001 [erratum: Phys. Rev. Lett. 116 (2016) no.18, 189903] doi:10.1103/PhysRevLett.116.062001 [arXiv:1511.05409 [hep-ph]].
  * (24) D. Chicherin, J. Henn and V. Mitev, “Bootstrapping pentagon functions,” JHEP 05 (2018), 164 doi:10.1007/JHEP05(2018)164 [arXiv:1712.09610 [hep-th]].
  * (25) T. Gehrmann, J. M. Henn and N. A. Lo Presti, “Pentagon functions for massless planar scattering amplitudes,” JHEP 10 (2018), 103 doi:10.1007/JHEP10(2018)103 [arXiv:1807.09812 [hep-ph]].
  * (26) D. Chicherin, T. Gehrmann, J. M. Henn, N. A. Lo Presti, V. Mitev and P. Wasser, “Analytic result for the nonplanar hexa-box integrals,” JHEP 03 (2019), 042 doi:10.1007/JHEP03(2019)042 [arXiv:1809.06240 [hep-ph]].
  * (27) S. Abreu, B. Page and M. Zeng, “Differential equations from unitarity cuts: nonplanar hexa-box integrals,” JHEP 01 (2019), 006 doi:10.1007/JHEP01(2019)006 [arXiv:1807.11522 [hep-th]].
  * (28) D. Chicherin, T. Gehrmann, J. M. Henn, P. Wasser, Y. Zhang and S. Zoia, “All Master Integrals for Three-Jet Production at Next-to-Next-to-Leading Order,” Phys. Rev. Lett. 123 (2019) no.4, 041603 doi:10.1103/PhysRevLett.123.041603 [arXiv:1812.11160 [hep-ph]].
  * (29) D. Chicherin and V. Sotnikov, “Pentagon Functions for Scattering of Five Massless Particles,” JHEP 12 (2020), 167 doi:10.1007/JHEP12(2020)167 [arXiv:2009.07803 [hep-ph]].
  * (30) S. Abreu, H. Ita, F. Moriello, B. Page, W. Tschernow and M. Zeng, “Two-Loop Integrals for Planar Five-Point One-Mass Processes,” JHEP 11 (2020), 117 doi:10.1007/JHEP11(2020)117 [arXiv:2005.04195 [hep-ph]].
  * (31) D. D. Canko, C. G. Papadopoulos and N. Syrrakos, “Analytic representation of all planar two-loop five-point Master Integrals with one off-shell leg,” JHEP 01 (2021), 199 doi:10.1007/JHEP01(2021)199 [arXiv:2009.13917 [hep-ph]].
  * (32) F. Moriello, “Generalised power series expansions for the elliptic planar families of Higgs + jet production at two loops,” JHEP 01 (2020), 150 doi:10.1007/JHEP01(2020)150 [arXiv:1907.13234 [hep-ph]].
  * (33) M. Hidding, “DiffExp, a Mathematica package for computing Feynman integrals in terms of one-dimensional series expansions,” [arXiv:2006.05510 [hep-ph]].
  * (34) C. G. Papadopoulos, “Simplified differential equations approach for Master Integrals,” JHEP 07 (2014), 088 doi:10.1007/JHEP07(2014)088 [arXiv:1401.6057 [hep-ph]].
  * (35) N. Syrrakos, “Pentagon integrals to arbitrary order in the dimensional regulator,” JHEP 06 (2021), 037 doi:10.1007/JHEP06(2021)037 [arXiv:2012.10635 [hep-ph]].
  * (36) S. Badger, H. B. Hartanto and S. Zoia, “Two-Loop QCD Corrections to $Wb\bar{b}$ Production at Hadron Colliders,” Phys. Rev. Lett. 127 (2021) no.1, 012001 doi:10.1103/PhysRevLett.127.012001 [arXiv:2102.02516 [hep-ph]].
  * (37) C. G. Papadopoulos and C. Wever, “Internal Reduction method for computing Feynman Integrals,” JHEP 02 (2020), 112 doi:10.1007/JHEP02(2020)112 [arXiv:1910.06275 [hep-ph]].
  * (38) D. D. Canko and N. Syrrakos, “Resummation methods for Master Integrals,” JHEP 02 (2021), 080 doi:10.1007/JHEP02(2021)080 [arXiv:2010.06947 [hep-ph]].
  * (39) A. Georgoudis, K. J. Larsen and Y. Zhang, “Azurite: An algebraic geometry based package for finding bases of loop integrals,” Comput. Phys. Commun. 221 (2017), 203-215 doi:10.1016/j.cpc.2017.08.013 [arXiv:1612.04252 [hep-th]].
  * (40) J. Klappert, F. Lange, P. Maierhöfer and J. Usovitsch, “Integral reduction with Kira 2.0 and finite field methods,” Comput. Phys. Commun. 266 (2021), 108024 doi:10.1016/j.cpc.2021.108024 [arXiv:2008.06494 [hep-ph]].
  * (41) M. Heller, A. von Manteuffel and R. M. Schabinger, “Multiple polylogarithms with algebraic arguments and the two-loop EW-QCD Drell-Yan master integrals,” Phys. Rev. D 102 (2020) no.1, 016025 doi:10.1103/PhysRevD.102.016025 [arXiv:1907.00491 [hep-th]].
  * (42) M. Besier, D. Van Straten and S. Weinzierl, “Rationalizing roots: an algorithmic approach,” Commun. Num. Theor. Phys. 13 (2019), 253-297 doi:10.4310/CNTP.2019.v13.n2.a1 [arXiv:1809.10983 [hep-th]].
  * (43) M. Besier, P. Wasser and S. Weinzierl, “RationalizeRoots: Software Package for the Rationalization of Square Roots,” Comput. Phys. Commun. 253 (2020), 107197 doi:10.1016/j.cpc.2020.107197 [arXiv:1910.13251 [cs.MS]].
  * (44) M. Bonetti, E. Panzer, V. A. Smirnov and L. Tancredi, “Two-loop mixed QCD-EW corrections to $gg\to Hg$,” JHEP 11 (2020), 045 doi:10.1007/JHEP11(2020)045 [arXiv:2007.09813 [hep-ph]].
  * (45) B. Jantzen, A. V. Smirnov and V. A. Smirnov, “Expansion by regions: revealing potential and Glauber regions automatically,” Eur. Phys. J. C 72 (2012), 2139 doi:10.1140/epjc/s10052-012-2139-2 [arXiv:1206.0546 [hep-ph]].
  * (46) A. V. Smirnov, “FIESTA4: Optimized Feynman integral calculations with GPU support,” Comput. Phys. Commun. 204 (2016), 189-199 doi:10.1016/j.cpc.2016.03.013 [arXiv:1511.03614 [hep-ph]].
  * (47) T. Huber and D. Maitre, “HypExp 2, Expanding Hypergeometric Functions about Half-Integer Parameters,” Comput. Phys. Commun. 178 (2008), 755-776 doi:10.1016/j.cpc.2007.12.008 [arXiv:0708.2443 [hep-ph]].
  * (48) C. Duhr and F. Dulat, “PolyLogTools — polylogs for the masses,” JHEP 08 (2019), 135 doi:10.1007/JHEP08(2019)135 [arXiv:1904.07279 [hep-th]].
  * (49) S. Borowka, G. Heinrich, S. Jahn, S. P. Jones, M. Kerner, J. Schlenk and T. Zirke, “pySecDec: a toolbox for the numerical evaluation of multi-scale integrals,” Comput. Phys. Commun. 222 (2018), 313-326 doi:10.1016/j.cpc.2017.09.015 [arXiv:1703.09692 [hep-ph]].
  * (50) L. Naterop, A. Signer and Y. Ulrich, “handyG —Rapid numerical evaluation of generalised polylogarithms in Fortran,” Comput. Phys. Commun. 253 (2020), 107165 doi:10.1016/j.cpc.2020.107165 [arXiv:1909.01656 [hep-ph]].
  * (51) J. Vollinga and S. Weinzierl, “Numerical evaluation of multiple polylogarithms,” Comput. Phys. Commun. 167 (2005), 177 doi:10.1016/j.cpc.2004.12.009 [arXiv:hep-ph/0410259 [hep-ph]].
  * (52) F. Brown and C. Duhr, “A double integral of dlog forms which is not polylogarithmic,” [arXiv:2006.09413 [hep-th]].
  * (53) S. Abreu, H. Ita, B. Page and W. Tschernow, “Two-Loop Hexa-Box Integrals for Non-Planar Five-Point One-Mass Processes,” [arXiv:2107.14180 [hep-ph]].