Could you tell us how this problem arose? In some contexts, the fact that the integral diverges means that the transform does not exist and there is no more work to be done. Of course WolframAlpha is not using that approach! Instead they are taking the
finite part of the divergent integral, or equivalently they are using the theory of distributions (generalized functions) and treating ##\cos t / t## as a pseudofunction.
I get the same expression that WolframAlpha does when I take the finite part. Since the
problem singularity is at zero, the basic idea is that you consider
$$
F(s) = \lim_{\epsilon \rightarrow \ 0^+} \int_\epsilon^\infty e^{-s t} \frac{\cos t}{t} \, dt
$$
Integrating by parts yields
$$
F(s) = \lim_{\epsilon \rightarrow \ 0^+} \left[-e^{-s \epsilon} \cos\epsilon \ln\epsilon + s \int_\epsilon^\infty e^{-s t}\, \cos t \, \ln t\,dt + \int_\epsilon^\infty e^{-s t}\, \sin t\, \ln t\, dt\right]
$$
The two integrals above converge as ##\epsilon \rightarrow 0^+##, but the first term is divergent. The
finite part is then obtained by simply throwing out the divergent term. If we let ##G(s)## denote the finite part of ##F(s)##, then we have,
$$
G(s) = s \int_0^\infty e^{-s t}\, \cos t \, \ln t \, dt + \int_0^\infty e^{-s t}\, \sin t\, \ln t \, dt
$$
which is nicely behaved by construction. I'm sure there are many ways to tackle this, but as with many integral transform problems going the direct approach isn't always the easiest. The approach I took was to Taylor expand ##\cos t## and ##\sin t## and then integrate term by term, which works because the Taylor series converge for all ##t##. Now we need to find the Laplace transform of terms that look like ##t^k \ln t##, which I computed using a trick that I have seen with similar problems in the past. Start with the definition of the Gamma function
$$ \Gamma(z) = \int_0^\infty e^{-t}\, t^{z-1}\, dt.$$
The derivative is then,
$$ \begin{eqnarray*}
\Gamma^\prime(z) & = & \int_0^\infty e^{-t}\, t^{z-1}\, \ln t \, dt \\
& = & s \int_0^\infty e^{-s u}\, (su)^{z-1}\, \ln(su) \, du \\
& = & \Gamma(z) \ln s + s^z \int_0^\infty e^{-s t}\, t^{z-1}\, \ln t \, dt
\end{eqnarray*} $$
Where in the second line I simply did a substitution ##t = su##. The last line gives us the Laplace transforms we need. I then insert this into the series, collect like terms in powers of ##s##, and the answer falls out. Note that you will need properties of ##\Gamma^\prime(z) = \Gamma(z) \psi(z)## where ##\psi(z)## is the digamma function, and the fact that ##\psi(1) = -\gamma##. See the first handful pages of
https://dlmf.nist.gov/5 to get all of the relations you need. I needed ##\psi(n+1)-\psi(n) = 1/n##, for example. At the end I had to sum a series, but it is the kind of series that most of us had to be able to identify in basic calculus class since it is directly related to the series for ##\ln(1+x)##.
Hope that helps,
Jason
