Prove $\int_0^\infty \frac{x^{k-1} + x^{-k-1}}{x^a + x^{-a}}dx = \frac{\pi}{a \cos(\frac{\pi k}{2a})}$.

Question

I need help in proving this identity

$$\int_0^\infty \frac{x^{k-1} + x^{-k-1}}{x^a + x^{-a}}dx = \frac{\pi}{a \cos(\frac{\pi k}{2a})}$$

for $0<k<a.$ It might be done using residues, but I don't know which contour to choose.

Answer 1

Let $x=e^u$; then the integral is equal to

$$\int_{-\infty}^{\infty} du\, \frac{\cosh{k u}}{\cosh{a u}} $$

This integral may be evaluated by considering the following contour integral

$$\oint_C dz \, \frac{\cosh{k z}}{\cosh{a z}} $$

where $C$ is the rectangle having vertices $-R$, $R$, $R+i \pi/a$, and $-R+i \pi/a$. The contour integral is equal to

$$\int_{-R}^R dx \, \frac{\cosh{k x}}{\cosh{a x}} + i \int_0^{\pi/a} dy \, \frac{\cosh{k (R+i y)}}{\cosh{a (R+i y)}} \\ + \int_{-R}^R dx \, \frac{\cosh{[k (x+i \pi/a)]}}{\cosh{a x}}+i \int_{\pi/a}^0 dy \, \frac{\cosh{k (-R+i y)}}{\cosh{a (-R+i y)}}$$

Note that, because $0 \lt k \lt a$, the second and fourth integrals vanish as $R \to \infty$. Thus, in this limit, the contour integral is equal to

$$[1+ \cos{(\pi k/a)} ]\int_{-\infty}^{\infty} dx \, \frac{\cosh{k x}}{\cosh{a x}} -i \sin{(\pi k/a)} \int_{-\infty}^{\infty} dx \, \frac{\sinh{k x}}{\cosh{a x}} $$

Note that the second integral is zero because it is an odd integral over a symmetric interval.

By the residue theorem, the contour integral is also equal to $i 2 \pi$ times the residue of the integrand at $z=i \pi/(2 a)$, or $(2 \pi/a) \, \cos{[\pi k/(2 a)]} $. The result is therefore

$$\int_{-\infty}^{\infty} dx \, \frac{\cosh{k x}}{\cosh{a x}} = \frac{\pi}{a} \frac{2 \cos{[\pi k/(2 a)]}}{1+ \cos{(\pi k/a)}} = \frac{\pi}{a} \sec{\left (\frac{\pi k}{2 a} \right )} $$

as was to be shown.

Answer 2

\begin{align*} I &= \int^\infty_0 \frac{x^{k - 1} + x^{-kn - 1}}{x^a + x^{-a}} \, dx\\ &= \int^\infty_0 \frac{x^{k - 1} + x^{-k - 1}}{x^{-a} (x^{2a} + 1)} \, dx\\ &= \int^\infty_0 \frac{x^{k + a -1} + x^{a - k - 1}}{x^{2a} + 1} \, dx \end{align*} Let $t = x^{2a}, x = t^{\frac{1}{2a}}, dx = \frac{1}{2a} t^{\frac{1 - 2a}{2a}} \, dt$ while the limits of integration remain unchanged. Thus \begin{align*} I &= \int^\infty_0 \frac{t^{\frac{k + a - 1}{2a}} + t^{\frac{a - k - 1}{2a}}}{1 + t} \cdot \frac{1}{2a} t^{\frac{1 - 2a}{2a}} \, dt\\ &= \frac{1}{2a} \int^\infty_0 \int^\infty_0 \frac{t^{\frac{k - a}{2a}} + t^{\frac{-a-k}{2a}}}{1 + t} \, dt\\ &= \frac{1}{2a} \int^\infty_0 \frac{t^{\frac{k + a}{2a} - 1}}{(1 + t)^{\frac{k + a}{2a} + \frac{a - k}{2a}}} \, dt + \frac{1}{2a} \int^\infty_0 \frac{t^{\frac{a - k}{2a} - 1}}{(1 + t)^{\frac{a - k}{2a} + \frac{k + a}{2a}}} \, dt\\ &= \frac{1}{2a} \mbox{B} \left (\frac{k + a}{2a}, \frac{a - k}{2a} \right ) + \frac{1}{2a} \mbox{B} \left (\frac{k - a}{2a}, \frac{a + k}{2a} \right ) \end{align*} Here $\displaystyle{\mbox{B} (m,n)}$ is the beta function. Since $\mbox{B}(m,n) = \mbox{B}(n,m)$ we have \begin{align*} I &= \frac{1}{a} \mbox{B} \left (\frac{k + a}{2a}, \frac{a - k}{2a} \right )\\ &= \frac{1}{a} \cdot \frac{\Gamma \left (\frac{k + a}{2a} \right ) \Gamma \left (\frac{a - k}{2a} \right )}{\Gamma \left (\frac{k + a}{2a} + \frac{a - k}{2a} \right )}, \quad \mbox{since} \,\, \mbox{B}(m,n) = \frac{\Gamma (m) \Gamma (n)}{\Gamma (m + n)}\\ &= \frac{1}{a} \Gamma \left (\frac{1}{2} + \frac{k}{2a} \right ) \Gamma \left (\frac{1}{2} - \frac{k}{2a} \right ), \quad \mbox{since} \,\, \Gamma (1) = 1\\ &= \frac{\pi}{a \cos \left (\frac{k\pi}{2} \right )} \end{align*} as required to show. Note in the last line the following reflection formula for the gamma function has been used $$\Gamma \left (\frac{1}{2} + x \right ) \Gamma \left (\frac{1}{2} - x \right ) = \frac{\pi}{\cos (x \pi)}.$$

Answer 3

From this post, we know that $$ \int_0^\infty\frac{u^{\alpha-1}}{1+u^\beta}dx=\frac{\pi}{\beta\sin(\frac{\pi\alpha}{\beta})}.$$ So \begin{eqnarray} \int_0^\infty \frac{x^{k-1} + x^{-k-1}}{x^a + x^{-a}}dx &=& \int_0^\infty \frac{x^{a+k-1}}{1 + x^{2a}}dx+\int_0^\infty \frac{x^{a-k-1}}{1 + x^{2a}}dx\\ &=&\frac{\pi}{2a\sin(\frac{(a+k-1)\pi}{2a})}+\frac{\pi}{2a\sin(\frac{(a-k-1)\pi}{2a})}\\ &=& \frac{\pi}{a \cos(\frac{\pi k}{2a})} \end{eqnarray}

Answer 4

Suppose we seek to evaluate

$$\int_0^\infty \frac{x^{k-1}+x^{-k-1}}{x^a+x^{-a}} dx$$

where $0\lt k \lt a$ directly using residues without additional substitutions.

We evaluate this for $a$ rational so that $a=p/q$ and extend by continuity to $a$ real. We take $p$ and $q$ in lowest terms.

We use a pizza slice contour of radius $R$ with the branch cut of the logarithm on the positive real axis and argument from $0$ to $2\pi.$ Moreover we treat the two summands $J_1$ and $J_2$ appearing in the numerator in turn.

We set the sloped component of the slice at angle $\pi q/p$ and we obtain with $x$ real $$(x\exp(\pi i q/p))^a = \exp(a\log(x\exp(\pi i q/p))) \\= \exp(a\log x + \pi i) = -x^a$$

and $$(x\exp(\pi i q/p))^{-a} = \exp(-a\log(x\exp(\pi i q/p))) \\= \exp(-a\log x - \pi i) = -x^{-a}.$$

We also have on the sloped segment for $J_1$ $$(x\exp(\pi i q/p))^{k-1} = \exp((k-1)\log(x\exp(\pi i q/p))) \\ = x^{k-1} \exp(\pi i q (k-1)/p)$$

and for $J_2$

$$(x\exp(\pi i q/p))^{-k-1} = \exp((-k-1)\log(x\exp(\pi i q/p))) \\ = x^{-k-1} \exp(\pi i q (-k-1)/p).$$

Now by ML on the curved segment of radius $R$ of $J_1$ we get the bound $$\pi q/p \times R \times \frac{R^{k-1}}{R^a-R^{-a}} \rightarrow 0$$ as $R\rightarrow\infty$ because $k \lt a.$

Similarly for $J_2$ we obtain $$\pi q/p \times R \times \frac{R^{-k-1}}{R^a-R^{-a}} \rightarrow 0$$ as $R\rightarrow\infty$ as well.

We also need to apply ML on the circular indentation segment of radius $\epsilon$ excluding the origin where we get for $J_1$ with $\epsilon$ going to zero

$$\pi q/p \times \epsilon \times \frac{\epsilon^{k-1}}{(1/\epsilon)^a-\epsilon^{a}} \rightarrow 0$$

as $\epsilon\rightarrow 0$ because the numerator vanishes and the first term in the denominator dominates, going to infinity.

Similarly for $J_2$ we obtain $$\pi q/p \times \epsilon \times \frac{\epsilon^{-k-1}}{(1/\epsilon)^a-\epsilon^{a}} \rightarrow 0$$

as $\epsilon\rightarrow 0$ because the numerator is $(1/\epsilon)^k$ and the denominator is $(1/\epsilon)^a$ and $k\lt a.$

The zeros of $x^a + x^{-a} = x^{-a} (x^{2a}+1)$ and hence the poles are at $$\rho_m = \exp((2\pi i m + \pi i) /(2p/q)) = \exp(\pi i m q/p + \pi i q/2/p).$$

The first pole at $\rho_0 = \exp(\pi i q/2/p)$ is inside the slice but the next one $\rho_1 = \exp(3\pi i q/2/p)$ is outside the slice already.

We have periodicity when $m$ is a multiple of $2p$ and $\rho_0$ re-appears so $\rho_0$ is in fact the only pole we have to consider. (Here we have used the fact that $p/q$ is in lowest terms.)

We get the following residue for $J_1$:

$$\frac{\rho_0^{k-1}}{a\rho_0^{a-1}-a\rho_0^{-a-1}} = \frac{\rho_0^{k}}{a\rho_0^{a}-a\rho_0^{-a}} = \frac{\rho_0^{k}}{a\rho_0^{a}+a\rho_0^{a}} = \frac{1}{2a} \rho_0^{k-a} = - \frac{i}{2a} \rho_0^{k}.$$

Similarly we get for $J_2$ $$\frac{\rho_0^{-k-1}}{a\rho_0^{a-1}-a\rho_0^{-a-1}} = \frac{\rho_0^{-k}}{a\rho_0^{a}-a\rho_0^{-a}} = \frac{\rho_0^{-k}}{a\rho_0^{a}+a\rho_0^{a}} = \frac{1}{2a} \rho_0^{-k-a} = - \frac{i}{2a} \rho_0^{-k}.$$

Collecting everything we have

$$J_1(1+\rho_0^{2k}) = - 2\pi i \frac{i}{2a} \rho_0^{k} \quad\text{or}\quad J_1 = \frac{\pi}{a} \frac{\rho_0^k}{1+\rho_0^{2k}} = \frac{\pi}{a} \frac{1}{\rho_0^{-k}+\rho_0^{k}} $$

and $$J_2(1+\rho_0^{-2k}) = - 2\pi i \frac{i}{2a} \rho_0^{-k} \quad\text{or}\quad J_2 = \frac{\pi}{a} \frac{\rho_0^{-k}}{1+\rho_0^{-2k}} = \frac{\pi}{a} \frac{1}{\rho_0^k+\rho_0^{-k}}$$

Now we have $$\frac{1}{\rho_0^k+\rho_0^{-k}} = \frac{1}{\exp(\pi i k/a/2) + \exp(-\pi i k/a/2)} = \frac{1}{2} \frac{1}{\cos(\pi k/a/2)}.$$

Adding these yields the end result

$$\frac{\pi}{a} \frac{1}{\cos(\pi k/a/2)}.$$

Important remark. When verifying this computation with a CAS we need to re-define the logarithm and the power function, e.g. in Maple we use

LOG :=
proc(z)
    local w;

    w := log(z);

    if evalf(Im(w)) < 0 then
        w := w + 2*Pi*I;
    fi;

    w;
end;

POW := (v, q) -> exp(q*LOG(v));


Q :=
proc(p, q, k)
    local a, cont, rho, res;

    a := p/q;

    rho := exp(Pi*I/a/2);
    res := POW(rho, k)/(a*POW(rho, a)-a*POW(rho, -a));
    # res := POW(rho, k-a)/2/a;

    2*Pi*I*res/(1 + exp(I*Pi/a*k));
    # 2*Pi*I*res/(1 + rho^(2*k));
end;

QQ := (p, q, k) -> 
int(x^(k-1)/(x^(p/q) + x^(-p/q)), x=0..infinity);

Of course all the powers of $\rho_0$ that appear in the above are calculated with this branch of the logarithm.

Observation. The angle $\pi q/p$ is taken modulo $2\pi$ when defining the contour (which does not make multiple turns) and when computing the bounds on the two arcs.

Highly significant addendum. It so happens I missed a crucial fact when I wrote the above, namely that the branch of the logarithm to use depends on the value of $a$! We must pick a branch where $$\rho_0^{2a} = -1$$ holds, otherwise the computation will not go through. E.g. suppose you want to solve the equation $\sqrt{x}+1 = 0.$ The branch to choose here must produce $-1$ as the root of $1$. We have for all branches that $\sqrt{x} = \exp(1/2\log(x)).$ Now if $\log(1) = 2\pi i$ we get the desired result. Apply this observation to obtain the following modified Maple routine:

LOG :=
proc(z, sheet)
    local w;

    w := log(z);

    if evalf(Im(w)) < 0 then
        w := w + 2*Pi*I;
    fi;

    w + sheet*2*Pi*I;
end;

POW := (v, q, sheet) -> exp(q*LOG(v, sheet));


Q :=
proc(p, q, k)
    local a, cont, rho, res, sheet, val;

    a := p/q;

    rho := exp(Pi*I/a/2);

    for sheet from 0 to 2*q-1 do
        val := 2*a*Im(LOG(rho, sheet))/Pi;

        if type(val, odd) then
            break;
        fi;
    od;

    if sheet = 2*q then return FAIL fi;

    res := 
    POW(rho, k, sheet)/
    (a*POW(rho, a, sheet)-a*POW(rho, -a, sheet));
    # res := POW(rho, k-a)/2/a;

    [sheet, 2*Pi*I*res/(1 + exp(I*Pi/a*k))];
    # 2*Pi*I*res/(1 + rho^(2*k));
end;

QQ := (p, q, k) ->
int(x^(k-1)/(x^(p/q) + x^(-p/q)), x=0..infinity);

TEST :=
proc()
    local cases, inst, p, q, k;

    cases :=
    [[3, 64, 1/41], [5, 64, 1/41], [1, 32, 1/41],
     [1, 10, 1/11], [5, 3, 1], [7, 3, 2],
     [5, 48, 1/19], [7, 48, 1/13],
     [11, 16, 1/3], [13, 16, 1/8],
     [3, 16, 1/10], [1, 16, 1/32],
     [19, 4, 11/3], [21, 7, 32/19]];

    for inst in cases do
        p := inst[1]; q := inst[2]; k := inst[3];

        print(evalf(Q(p, q, k)),
              evalf(QQ(p, q, k)));
    od;
end;

This produces the following results for the test cases:

> TEST();
                                            -8
          [5., 48.97919294 + 0.2273921483 10   I], 48.97919294

                                            -9
          [3., 22.79234062 + 0.6739119228 10   I], 22.79234062

                                            -7
          [8., 148.7070193 + 0.1069863986 10   I], 148.7070193

                                            -8
          [2., 110.3747402 - 0.8833572709 10   I], 110.3747400

                                            -9
          [0., 1.603438998 + 0.5642199691 10   I], 1.603438998

                                            -9
          [0., 3.025326263 + 0.7216433442 10   I], 3.025326264

                                            -8
          [2., 21.50435498 + 0.1066992697 10   I], 21.50435498

                                            -8
          [1., 15.93501901 - 0.5378737091 10   I], 15.93501900

                                            -8
          [0., 3.156953436 + 0.1156017429 10   I], 3.156953434

                                            -9
          [0., 1.991146897 + 0.1577177352 10   I], 1.991146897

                                            -9
          [1., 12.52009747 + 0.6173768878 10   I], 12.52009747

                  [4., 35.54306350 + 0. I], 35.54306350

                                            -9
         [0., 0.9431219183 + 0.5926932608 10   I], 0.9431219179

                                            -9
         [0., 0.8236262644 - 0.4232839786 10   I], 0.8236262640

> quit
memory used=23.4MB, alloc=44.3MB, time=0.48