How does one prove $\int_0^\infty \frac{\log(x)}{1 + e^{ax}} \, dx = -\frac{\log(2)(2\log(a) + \log(2))}{2a}$ for $a > 0$? Link to WolframAlpha's assertion. Here's my attempt. Using the substitution $t = ax$, we can show the integral is equal to 
$$ \frac{1}{a} \int_0^\infty \frac{\log(t)}{1 + e^t}\, dt -\frac{\log(a)}{a} \int_0^\infty \frac{dt}{1 + e^{t} } .$$
The second integral is equal to $\log(2)$ using another substitution $v = e^t$ and partial fractions.
So I'm left with the first integral. I'll switch to complex variables for notation. I make two observations:
(I) The denominator has simple poles when $z = t = (2k-1)\cdot i \pi, \, k \in \mathbb{N}.$ 
(II) The numerator has a branch point $z = 0$.
 A: Let $I$ be the integral given by 

$$\bbox[5px,border:2px solid #C0A000]{I=\int_0^\infty\frac{\log(x)}{1+e^x}\,dx} \tag 1$$

Expanding the denominator of $(1)$ in a series of $e^{-nx}$ and interchanging the order of summation and integration reveals
$$\begin{align}
I&=\int_0^\infty\frac{e^{-x}\log(x)}{1+e^{-x}}\,dx\\\\
&=\int_0^\infty\log(x)\sum_{n=0}^\infty (-1)^ne^{-(n+1)x}\,dx\\\\
&=\sum_{n=0}^\infty (-1)^n \int_0^\infty\log(x)e^{-(n+1)x}\,dx\\\\
&=\sum_{n=0}^\infty \frac{(-1)^n}{n+1}\int_0^\infty e^{-x}(\log(x)-\log(n+1))\,dx \\\\
&=\sum_0^\infty \frac{(-1)^{n+1}(\gamma+\log(n+1))}{n+1}\\\\
&=-\gamma\log(2)+\color{blue}{\sum_{n=1}^\infty \frac{(-1)^{n}\log(n)}{n}}\\\\
&=-\gamma\log(2)+\color{blue}{\eta'(1)}\\\\
&=-\gamma\log(2)+\gamma\log(2)-\frac12\log^2(2)\\\\
&=-\frac12\log^2(2)
\end{align}$$
where $\eta'(s)$ is the derivative of the Dirichlet Eta Function with $\eta'(1)=\gamma\log(2)-\frac12\log^2(2)$ (SEE THIS ANSWER).

Therefore, we have the coveted result
$$\bbox[5px,border:2px solid #C0A000]{\int_0^\infty\frac{\log(x)}{1+e^x}\,dx=-\frac12\log^2(2)}$$

which agrees with that obtained using Wolfram Alpha! 

NOTE $1$:  DIRECT EVALUATION OF THE SUM $\displaystyle \sum_{n=1}^\infty \frac{(-1)^{n}\log(n)}{n}$

Here, we provide a direct evaluation of the sum $\displaystyle \sum_{n=1}^\infty \frac{(-1)^{n}\log(n)}{n}$ facilitated by use of the Euler-Maclaurin Summation Formula (EMSF).  To that end we proceed.

First, we note that we can write for any $N\ge 1$
$$\begin{align}
\sum_{n=1}^{2N}\frac{(-1)^{n}\log(n)}{n}&=2\sum_{n=1}^{N}\frac{\log(2n)}{2n}-\sum_{n=1}^{2N}\frac{\log(n)}{n}\\\\
&=\log(2)\sum_{n=1}^N\frac1n -\sum_{n=N+1}^{2N}\frac{\log(n)}{n} \tag 2
\end{align}$$
Using the EMSF, we expand the second sum on the right-hand side of $(2)$ to obtain
$$\begin{align}
\sum_{n=N+1}^{2N}\frac{\log(n)}{n}&=\int_N^{2N}\frac{\log(x)}{x}\,dx+O\left(\frac{\log(N)}{N}\right)\\\\
&=\frac12\log^2(2N)-\frac12\log^2(N)+O\left(\frac{\log(N)}{N}\right)\\\\
&=\log(2)\log(N)+\frac12\log^2(2)+O\left(\frac{\log(N)}{N}\right) \tag 3
\end{align}$$
Substituting $(3)$ in $(2)$ reveals
$$\begin{align}
\sum_{n=1}^{\infty}\frac{(-1)^{n}\log(n)}{n}&=\lim_{N\to \infty}\sum_{n=1}^{2N}\frac{(-1)^{n}\log(n)}{n}\\\\
&=\lim_{N\to \infty}\left(\log(2)\sum_{n=1}^n\frac1n -\log(2)\log(N)-\frac12\log^2(2)+O\left(\frac{\log(N)}{N}\right)\right)\\\\
&=\log(2)\lim_{N\to \infty}\left(\sum_{n=1}^N \frac1n -\log(N)\right)-\frac12\log^2(2)\\\\
&=\gamma \log(2)-\frac12\log^2(2)
\end{align}$$
Therefore, we find that 

$$\bbox[5px,border:2px solid #C0A000]{\sum_{n=1}^{\infty}\frac{(-1)^{n}\log(n)}{n}=\gamma \log(2)-\frac12\log^2(2)}$$

A: Here is another solution, although it is not as efficient as Dr. MV's solution.
Let
$$I = \int\limits_{0}^{\infty} \frac{\ln(x)}{1+e^{ax}} \mathrm{d}x$$ and
\begin{align}
I_{1} &= \int\limits_{0}^{\infty} \frac{x^{b}}{1+e^{ax}} \mathrm{d}x = \int\limits_{0}^{\infty} \frac{x^{b}e^{-ax}}{1+e^{-ax}} \mathrm{d}x \\
&= \sum\limits_{n=0}^{\infty} (-1)^{n} \int\limits_{0}^{\infty} x^{b} e^{-(a+na)x} \mathrm{d}x
\end{align}
We designate the last integral on the right as $I_{2}$ and make the substitution $y=(a+na)x$
\begin{align}
I_{2} &= \int\limits_{0}^{\infty} x^{b} e^{-(a+na)x} \mathrm{d}x \\
&= \frac{1}{(a+na)^{b+1}} \int\limits_{0}^{\infty} y^{b} e^{-y} \mathrm{d}y \\
&= \frac{\Gamma(b+1)}{(a+na)^{b+1}} \\
&= \frac{\Gamma(b+1)}{a^{b+1}} \frac{1}{(n+1)^{b+1}}
\end{align}
Now $I_{1}$ becomes
\begin{equation}
I_{1} = \frac{\Gamma(b+1)}{a^{b+1}} \sum\limits_{n=0}^{\infty} (-1)^{n} \frac{1}{(n+1)^{b+1}} = \frac{\Gamma(b+1)}{a^{b+1}} \eta(b+1)
\end{equation}
and we have
\begin{align}
I &= \lim_{b \to 0} \frac{\partial I_{1}}{\partial b} \\
&= \lim_{b \to 0} \int\limits_{0}^{\infty} \frac{x^{b} \ln(x)}{1+e^{ax}} \mathrm{d}x \\
&= \lim_{b \to 0}\frac{\partial}{\partial b}\frac{\Gamma(b+1)\eta(b+1)}{a^{b+1}} \\
&= \lim_{b \to 0} \frac{a^{b+1} \Big[ \Gamma(b+1)\eta^{\prime}(b+1) + \Gamma^{\prime}(b+1)\eta(b+1) \Big] - \Gamma(b+1)\eta(b+1)a^{b+1}\ln(a)}{\left(a^{b+1}\right)^{2}} \\
&= \lim_{b \to 0}  \frac{\Gamma(b+1)\eta^{\prime}(b+1) + \Gamma^{\prime}(b+1)\eta(b+1) - \Gamma(b+1)\eta(b+1)\ln(a)}{a^{b+1}} \\
\tag{a}
&= \frac{1}{a} \left(\Big[\gamma \ln(2) - \frac{1}{2} \ln^{2}(2)\Big] -\gamma \ln(2) - \ln(2)\ln(a) \right) \\
&= -\frac{\ln(2)}{a} \left(\frac{1}{2} \ln(2) + \ln(a) \right)
\end{align}
In step (a) we have
\begin{align}
\lim_{b \to 0} \Gamma^{\prime}(b+1)\eta(b+1) &= \lim_{b \to 0} \Gamma^(b+1)\psi(b+1)\eta(b+1) \\
&= \Gamma(1)\psi(1)\eta(1) \\
&= -\gamma \ln(2)
\end{align}
and
\begin{align}
\lim_{b \to 0} \Gamma(b+1)\eta^{\prime}(b+1) &= \lim_{s \to 1} \eta^{\prime}(s) \\
&= \lim_{s \to 1} \sum\limits_{n=0}^{\infty} (-1)^{n} \frac{\ln(n)}{n^{s}} \\
&= \gamma \ln(2) - \frac{1}{2} \ln^{2}(2)
\end{align}
See Dr. MV's solution for a proof of this result.
Notes:


*
    
*$\Gamma(z)$ is the Gamma function.
    
*$\eta(s)$ is the Dirichlet eta function.
    
*$\zeta(s)$ is the Riemann zeta function.
    
*$\psi(z)$ is the digamma function.
    
*$\gamma$ is the Euler-Mascheroni constant.

