Quit[]

<< MB.m

<< AMBRE.m

Fullintegral[{1}, {PR[k1 - k2, m, n1] 
            PR[k1, 0, n2] PR[k1 + p1, m, n3] PR[k2 + p1, 0, n4] PR[k2, m, 
              n5] }, {k2, k1}];

invariants = {p1^2 -> m^2};

IntPart[1]

SubLoop[integral]

IntPart[2]

repr = SubLoop[integral];

BarnesLemma[repr, 1]

fin = BarnesLemma[%, 2]

F01101OS = fin /. n1 -> 1 /. n2 -> 1 /. 
            n3 -> 1 /. n4 -> 1 /. n5 -> 1

F01101OS // InputForm

rules = MBoptimizedRules[F01101OS, eps -> 0, {}, {eps}]

integrals = MBcontinue[F01101OS, eps -> 0, rules]

ser4 = 
        MBexpand[{integrals}, 
          Exp[2 eps EulerGamma], {eps, 0, 2}];

num1 = MBintegrate[ser4, {m -> 1}]

F01101OS = (+ep*(-ln[2]*ln[2]*ln[2]*ln[2] - 
                      12*ln[2]*ln[2]*zt[2] + 12*ln[2]*zt[2] - 3*zt[3] + 
                      63/2*zt[4] - 24*Li[4, 1/2]) + 
                ep^2*(-2*ln[2]*ln[2]*ln[2]*ln[2] + 
                      6/5*ln[2]*ln[2]*ln[2]*ln[2]*ln[2] + 
                      24*ln[2]*ln[2]*ln[2]*zt[2] - 24*ln[2]*ln[2]*zt[2] + 
                      24*ln[2]*zt[2] - 93*ln[2]*zt[4] + 18*zt[2]*zt[3] - 
                      6*zt[3] + 63*zt[4] + 609/4*zt[5] - 48*Li[4, 1/2] - 
                      144*Li[5, 1/2]) + 6*ln[2]*zt[2] - 3/2*zt[3]) /. 
            Li[x_, y_] -> PolyLog[x, y] /. zt[x_] -> Zeta[x] /. 
        ln[x_] -> Log[x]

res = 
        Normal[Series[-(F01101OS) Gamma[1 + ep]^2, {ep, 0, 2}]] /. 
            EulerGamma -> 0 // Expand;

Expand[%] /. m -> 1 // N

num1[[1]]