Real Analysis Methodologies to show $\gamma =2\int_0^\infty \frac{\cos(x^2)-\cos(x)}{x}\,dx$ In THIS ANSWER, I used straightforward complex analysis to show that 

$$\gamma =2\int_0^\infty \frac{\cos(x^2)-\cos(x)}{x}\,dx \tag 1$$

where $\gamma =-\int_0^\infty \log(x) e^{-x}\,dx$ is the Euler-Mascheroni Constant.
The key in the derivation of $(1)$ was to transform the cosine terms into real exponential ones.  
To date, I have been unable to use strictly real analysis, without appealing to tabulated results of special functions (e.g., use of the $\text{Cin}(x)$ and $\text{Ci}(x)$ functions), to prove $(1)$.  
I have tried introducing a parameter and using "Feynman's Trick to augment the integral into something manageable.  Or somewhat equivalently, rewriting the integral in $(1)$ as a double integral and proceeding by exploiting Fubini-Tonelli.  

QUESTION:  What are ways to prove $(1)$ without relying on complex analysis and without simply appealing to tabulated relationships of special functions.  For example, stating that the $Ci(x)$ function is defined as $\text{Ci}(x)\equiv -\int_x^\infty\frac{\cos(t)}{t}\,dt=\gamma +\log(x) +\int_0^x \frac{\cos(t)-1}{t}\,dt$ is unsatisfactory unless one proves the latter equality.

 A: 
I thought it might be instructive to post a solution that leverages the Lemma posted by @Renascence_5.  To that end, we proceed. 

The Lemma proved in the aforementioned post is expressed as 
$$\int_0^\infty \frac{e^{-x^a}-\cos(x^a)}{x}\,dx=0 \tag 1$$
for $a>0$.

We now examine a generalized version of the integral of interest and write
$$I(a,b)=\int_0^\infty \frac{\cos(x^a)-\cos(x^b)}{x}\,dx \tag 2$$
for $a>0$ and $b>0$.

Using $(1)$ reveals that $(2)$ can be written as 
$$I(a,b)=\int_0^\infty \frac{e^{-x^a}-e^{-x^b}}{x}\,dx \tag 3$$
Next, we integrate by parts the integral in $(3)$ with $u=e^{-x^a}-e^{-x^b}$ and $v=\log(x)$ to obtain
$$\begin{align}
I(a,b)&=\int_0^\infty \left(ax^{a-1}e^{-x^a}-bx^{b-1}e^{-x^b}\right)\,\log(x)\,dx\\\\
&=\int_0^\infty ax^{a-1}e^{-x^a}\,\log(x)\,dx-\int_0^\infty bx^{b-1}e^{-x^b}\,\log(x)\,dx\\\\
&=\frac1a \int_0^\infty e^{-x}\,\log(x)\,dx-\frac1b \int_0^\infty e^{-x}\,\log(x)\,dx\\\\
&=-\left(\frac1a -\frac1b\right)\,\gamma
\end{align}$$
where we used the integral relationship $\gamma =-\int_0^\infty e^{-x}\,\log(x)$.

NOTE:
We can show that $\gamma$ as given by $\gamma=-\int_0^\infty e^{-x}\,\log(x)\,dx$ is equal to $\gamma$ as expressed by the limit $\gamma=\lim_{n\to \infty}\left(-\log(n)+\sum_{k=1}^n\frac1k\right)$ as follows.
$$\begin{align}
\int_0^\infty e^{-x}\,\log(x)\,dx&=\lim_{n\to \infty}\int_0^n \left(1-\frac xn\right)^n\,\log(x)\,dx\\\\
&=\lim_{n\to \infty} n \int_0^1 x^n \log(n(1-x))\,dx\\\\
&=\lim_{n\to \infty} n \left(\log(n) \int_0^1 x^n\,dx+\int_0^1 x^n\,\log(1-x)\,dx\right)\\\\
&=\lim_{n\to \infty} \left(\frac{n}{n+1}\log(n)+\int_0^1 x^n\,\log(1-x)\,dx\right)\\\\
&=\lim_{n\to \infty} \left(\frac{n}{n+1}\log(n)-n\sum_{k=1}^\infty \frac1k \int_0^1 x^{n+k}\,dx\right)\\\\
&=\lim_{n\to \infty} \left(\frac{n}{n+1}\log(n)-n \sum_{k=1}^\infty \frac1{k(k+n+1)}\right)\\\\
&=\lim_{n\to \infty} \frac{n}{n+1}\left(\log(n)- \sum_{k=1}^\infty \left(\frac1k-\frac1{k+n+1}\right)\right)\\\\
&=\lim_{n\to \infty} \frac{n}{n+1}\left(\log(n)- \sum_{k=1}^{n+1} \frac1k\right)\\\\
&=\lim_{n\to \infty} \left(\log(n)- \sum_{k=1}^{n} \frac1k\right)\\\\
\end{align}$$
as was to be shown!
A: It turns out that we have the following observation:

Observation. For a nice function $f : [0,\infty) \to \Bbb{C}$, we have
$$ \int_{\epsilon}^{\infty} \frac{f(x)}{x} \, dx = -f(0)\log\epsilon + c(f) + o(1) \qquad \text{as } \epsilon \to 0^+ \tag{1} $$
where the constant $c(f)$ is computed by
$$ c(f) = \lim_{R\to\infty}\left( \int_{0}^{R} \mathcal{L}f(s) \, ds - f(0)\log R\right) - f(0)\gamma. \tag{2} $$

The reasoning is surprisingly simple: First, define $g(x) = (f(x) - f(0)\mathbf{1}_{(0,1)}(x))/x$ and notice that
$$ \int_{\epsilon}^{\infty} \frac{f(x)}{x} \, dx = -f(0)\log\epsilon + \int_{\epsilon}^{\infty} g(x) \, dx. $$
Assuming that the LHS of $\text{(1)}$ exists for all $\epsilon > 0$ and that $f$ behaves nice near $x = 0$, this implies $\text{(1)}$. Next, notice that $c(f) = \mathcal{L}g(0)$ and that $-(\mathcal{L}g(s))' = \mathcal{L}f(s) - f(0) (1-e^{-s})/s$. Therefore
\begin{align*}
c(f)
&= \lim_{R\to\infty} \int_{0}^{R} (-\mathcal{L}g(s))' \, ds \\
&= \lim_{R\to\infty} \left( \int_{0}^{R} \mathcal{L}f(s) \, ds - f(0) (1 - e^{-R})\log R + f(0) \int_{0}^{R} e^{-s}\log s \, ds \right) \\
&= \lim_{R\to\infty} \left( \int_{0}^{R} \mathcal{L}f(s) \, ds - f(0) \log R \right) - f(0)\gamma.
\end{align*}

At this moment this is just a heuristic computation. For a broad class of functions for which the LHS of $\text{(1)}$ exists, however, this computation can be made rigorous. This is particularly true for our function $f(x) = \cos x$. Now plugging $\mathcal{L}f(s) = \frac{s}{s^2+1}$ shows that $c(f) = -\gamma$ and thus
$$ \int_{\epsilon}^{\infty} \frac{\cos x}{x} \, dx = -\log\epsilon - \gamma + o(1). $$
Plugging this asymptotics, we have
$$
\int_{\epsilon}^{\infty} \frac{\cos(x^2) - \cos x}{x} \, dx
= \frac{1}{2}\int_{\epsilon^2}^{\infty} \frac{\cos x}{x} \, dx - \int_{\epsilon}^{\infty} \frac{\cos x}{x} \, dx
= \frac{1}{2}\gamma + o(1)
$$
and the identity follows by letting $\epsilon \to 0^+$.

Here, the constant $c(f)$ can be thought as a regularized value of the divergent integral $\int_{0}^{\infty} \frac{f(x)}{x} \, dx$. This has the following nice properties (whenever they exist)


*

*$c$ is linear: $c(\alpha f(x) + \beta g(x)) = \alpha c(f) + \beta c(g)$.

*$c(f(x^p)) = \frac{1}{p}c(f)$ for $p > 0$,

*$c(f(px)) = c(f) - f(0)\log p$ for $p > 0$,


Together with some known values, we can easily compute other types of integrals. For instance, using the fact that $c(\cos x) = -\gamma$ and $c(e^{-x}) = -\gamma$, we have
\begin{align*}
\int_{0}^{\infty} \frac{\cos (x^p) - \exp(-x^q)}{x} \, dx
&= c\left\{ \cos (x^p) - \exp(-x^q) \right\} \\
&= \frac{1}{p}c(\cos x) - \frac{1}{q}c(e^{-x})
 = \gamma\left( \frac{1}{q} - \frac{1}{p}\right)
\end{align*}
for $p, q > 0$.
A: $\newcommand{\bbx}[1]{\,\bbox[8px,border:1px groove navy]{\displaystyle{#1}}\,}
 \newcommand{\braces}[1]{\left\lbrace\,{#1}\,\right\rbrace}
 \newcommand{\bracks}[1]{\left\lbrack\,{#1}\,\right\rbrack}
 \newcommand{\dd}{\mathrm{d}}
 \newcommand{\ds}[1]{\displaystyle{#1}}
 \newcommand{\expo}[1]{\,\mathrm{e}^{#1}\,}
 \newcommand{\ic}{\mathrm{i}}
 \newcommand{\mc}[1]{\mathcal{#1}}
 \newcommand{\mrm}[1]{\mathrm{#1}}
 \newcommand{\pars}[1]{\left(\,{#1}\,\right)}
 \newcommand{\partiald}[3][]{\frac{\partial^{#1} #2}{\partial #3^{#1}}}
 \newcommand{\root}[2][]{\,\sqrt[#1]{\,{#2}\,}\,}
 \newcommand{\totald}[3][]{\frac{\mathrm{d}^{#1} #2}{\mathrm{d} #3^{#1}}}
 \newcommand{\verts}[1]{\left\vert\,{#1}\,\right\vert}$

A '$\ds{\color{#f00}{complex\ like}}$' answer can still be useful for other users. So, that is the reason I like to put forward the following answer:

\begin{align}
&2\int_{0}^{\infty}{\cos\pars{x^{2}} - \cos\pars{x} \over x}\,\dd x =
2\lim_{\Lambda \to \infty}\bracks{%
\int_{0}^{\Lambda}{1 - \cos\pars{x} \over x}\,\dd x -
\int_{0}^{\Lambda}{1 - \cos\pars{x^{2}} \over x}\,\dd x}
\\[5mm] = &\
2\lim_{\Lambda \to \infty}\bracks{%
\int_{0}^{\Lambda}{1 - \cos\pars{x} \over x}\,\dd x -
{1 \over 2}\int_{0}^{\Lambda^{2}}{1 - \cos\pars{x} \over x}\,\dd x}
\label{1}\tag{1}
\end{align}


With $\ds{R > 0}$:

\begin{align}
&\int_{0}^{R}{1 - \cos\pars{x} \over x}\,\dd x =
\Re\int_{0}^{R}{1 - \expo{\ic x} \over x}\,\dd x
\\[5mm] = &
-\Re\int_{0}^{\pi/2}\bracks{1 - \exp\pars{\ic R\expo{\ic\theta}}}\ic\,\dd\theta -
\Re\int_{R}^{0}{1 - \expo{-y} \over y}\,\dd y
\\[5mm] = &\
-\int_{0}^{\pi/2}\sin\pars{R\cos\pars{\theta}}
\expo{-R\sin\pars{\theta}}\,\dd\theta +
\ln\pars{R}\pars{1 - \expo{-R}} - \int_{0}^{R}\ln\pars{y}\expo{-y}\,\dd y
\\[5mm] \stackrel{\mrm{as}\ R\ \to\ \infty}{\sim}\,\,\, &\
\ln\pars{R} - \int_{0}^{\infty}\ln\pars{y}\expo{-y}\,\dd y =
\bbx{\ds{\ln\pars{R} + \gamma}}
\label{2}\tag{2}
\end{align}

because

$$
\left\{\begin{array}{l}
\ds{0 < \verts{\int_{0}^{\pi/2}\sin\pars{R\cos\pars{\theta}}
\expo{-R\sin\pars{\theta}}\,\dd\theta} < \int_{0}^{\pi/2}
\exp\pars{-\,{2R \over \pi}\,\theta}\,\dd\theta =
{\pi \over 2}\,{1 - \expo{-R} \over R}}
\\[5mm]
\mbox{and}\
\ds{\int_{0}^{\infty}\ln\pars{y}\expo{-y}\,\dd y =
\left.\totald{}{\mu}\int_{0}^{\infty}y^{\mu}\expo{-y}\,\dd y\,
\right\vert_{\ \mu\ =\ 0} =
\left.\totald{\Gamma\pars{\mu + 1}}{\mu}\right\vert_{\ \mu\ =\ 0} = \Psi\pars{1} = -\gamma}
\end{array}\right.
$$

With \eqref{1} and \eqref{2}:
\begin{align}
2\int_{0}^{\infty}{\cos\pars{x^{2}} - \cos\pars{x} \over x}\,\dd x & =
2\braces{\bracks{\ln\pars{\Lambda} + \gamma} -
{1 \over 2}\bracks{\ln\pars{\Lambda^{2}} + \gamma}} = \bbx{\ds{\gamma}}
\end{align}
A: For
$$\Gamma '\left ( x \right )=\int_{0}^{\infty }e^{-t}t^{x-1}\ln t\, \mathrm{d}t$$
using
$$\ln t=\int_{0}^{\infty }\frac{e^{-s}-e^{-ts}}{s}\, \mathrm{d}s$$
we have
$$\Gamma '\left ( x \right )=\int_{0}^{\infty }e^{-t}t^{x-1}\int_{0}^{\infty }\frac{e^{-s}-e^{-ts}}{s}\, \mathrm{d}s\mathrm{d}t=\Gamma \left ( x \right )\int_{0}^{\infty }\left ( e^{-s}-\frac{1}{\left ( s+1 \right )^{x}} \right )\frac{\mathrm{d}s}{s}$$
Hence, let $x=1$ we get
$$\gamma =\int_{0}^{\infty }\left ( \frac{1}{s+1 }-e^{-s} \right )\frac{\mathrm{d}s}{s}$$
let $s=t^k,~k>0$, we get
$$\gamma =\int_{0}^{\infty }\left ( \frac{1}{t^{k}+1 }-e^{-t^{k}} \right )\frac{k\, \mathrm{d}t}{t}$$
So,let $k=a,b$
$$\frac{\gamma}{a} =\int_{0}^{\infty }\left ( \frac{1}{t^{a}+1 }-e^{-t^{a}} \right )\frac{ \mathrm{d}t}{t}~~,~~\frac{\gamma}{b} =\int_{0}^{\infty }\left ( \frac{1}{t^{b}+1 }-e^{-t^{b}} \right )\frac{ \mathrm{d}t}{t}$$
hence
$$\frac{\gamma}{b}-\frac{\gamma}{a} =\int_{0}^{\infty }\left [\left ( \frac{1}{t^{b}+1 }- \frac{1}{t^{a}+1 } \right )+\left ( e^{-x^a}-e^{-x^b} \right )  \right ]\frac{ \mathrm{d}t}{t}$$
then
$$\int_{0}^{\infty }\left ( \frac{1}{t^{b}+1 }- \frac{1}{t^{a}+1 } \right )\frac{\mathrm{d}t}{t}=\int_{0}^{1}\left ( \frac{1}{t^{b}+1 }- \frac{1}{t^{a}+1 } \right )\frac{\mathrm{d}t}{t}+\int_{1}^{\infty }\left ( \frac{1}{t^{b}+1 }- \frac{1}{t^{a}+1 } \right )\frac{\mathrm{d}t}{t}$$
let $t\rightarrow \dfrac{1}{t}$,we get
$$\int_{0}^{1}\left ( \frac{1}{t^{b}+1 }- \frac{1}{t^{a}+1 } \right )\frac{\mathrm{d}t}{t}=-\int_{1}^{\infty }\left ( \frac{1}{t^{b}+1 }- \frac{1}{t^{a}+1 } \right )\frac{\mathrm{d}t}{t}$$
So
$$\int_{0}^{\infty }\left ( \frac{1}{t^{b}+1 }- \frac{1}{t^{a}+1 } \right )\frac{\mathrm{d}t}{t}=0$$
and
$$\left ( \frac{1}{b}-\frac{1}{a} \right )\gamma =\int_{0}^{\infty }\frac{e^{-t^a}-e^{-t^b}}{t}\, \mathrm{d}t\tag1$$
Lemma:

$$\int_{0}^{\infty }\frac{e^{-t^a}-\cos t^a}{t}\, \mathrm{d}t=0~,~a>0$$

Proof:
Let $$f\left ( x \right )=\int_{0}^{\infty }\frac{e^{-t}-\cos t}{t}\, e^{-xt}\, \mathrm{d}t$$
so
$$f'\left ( x \right )=\int_{0}^{\infty }\left ( \cos t-e^{-t} \right )e^{-xt}\, \mathrm{d}t=\frac{x}{1+x^2}-\frac{1}{1+x}$$
hence
$$\int_{0}^{\infty }f'\left ( x \right ) \mathrm{d}x=\ln\frac{\sqrt{1+x^2}}{1+x}\Bigg|_{0}^{\infty }=0=f\left ( \infty  \right )-f\left ( 0 \right )$$
It's easy to see that $f\left ( \infty  \right )=0$,so
$$f\left ( 0 \right )=\int_{0}^{\infty }\frac{e^{-t}-\cos t}{t}\, \mathrm{d}t=0$$
Let $x^a\to t$, we get
$$a\int_{0}^{\infty }\frac{e^{-t^{a}}-\cos t^{a}}{t}\, \mathrm{d}t=0\Rightarrow \int_{0}^{\infty }\frac{e^{-t^{a}}-\cos t^{a}}{t}\, \mathrm{d}t=0\tag2$$
Now using $(1)$ and $(2)$, we get
$$\Large\boxed{\color{Blue} {\int_{0}^{\infty }\frac{\cos x^{a}-\cos x^b}{x}\, \mathrm{d}x=\left ( \frac{1}{b}-\frac{1}{a} \right )\gamma }}$$
