This answer describes an alternate approach to "squeezing the primes" (more accurately prime-powers) based on my related Math Overflow question on the inverse Mellin transform
$$p(x)=\mathcal{M}_s^{-1}\left[-\xi(s)\,\frac{\zeta'(s)}{s\,\zeta(s)^2}\right]\left(\frac{1}{x}\right)=\frac{1}{2 \pi i} \int\limits_{\alpha-i\, \infty}^{\alpha-i\, \infty} \left(\pi^{-\frac{s}{2}}\, (1-s)\, \Gamma\left(\frac{s}{2}+1\right)\, \frac{\zeta'(s)}{\zeta(s)}\right)\, \frac{x^s}{s}\, ds$$ $$=\underset{K\to\infty}{\text{lim}}\left(\sum\limits_{k=1}^K \Lambda(k) \left(\frac{2 \pi k^2}{x^2}-1\right) e^{-\frac{\pi k^2}{x^2}}\right)\tag{1}$$
where $\xi(s)$ is the Riemann Xi function, $\Lambda(k)$ is the von Mangoldt function, and the related Mellin transform
$$P(s)=s\, \mathcal{M}_x[p(x)](-s)=s \int\limits_0^\infty p(x)\, x^{-s-1}\, ds$$ $$=\pi^{-\frac{s}{2}}\, (1-s)\, \Gamma\left(\frac{s}{2}+1\right)\,\frac{\zeta'(s)}{\zeta(s)}\,,\quad\Re(s)>1\tag{2}$$
converges for $\Re(s)>1$.
The residue formula for $p(x)$, which I believe converges for $x>0$, is
$$p_o(x)=\log(2 \pi)+p_{-2n}(x)+p_\rho(x)\tag{3}$$
where $p_{-2n}(x)$ and $p_\rho(x)$ defined in formulas (4) and (5) below are sums over the trivial and non-trivial zeta zeros respectively.
$$p_{-2n}(x)=\sum\limits_{n=1}^{\infty} \frac{x^{-2 n} (-\pi)^n \left(\zeta'(-2 n)+n (2 n+1) \left(\zeta''(-2 n)+\zeta'(-2 n) (\psi^{(0)}(n)+2 \log(x)-\log(\pi))\right)\right)}{2 n\, n!\, \zeta'(-2 n)}\tag{4}$$
$p_\rho(x)=\sum\limits_{\rho} \left(\pi^{-\frac{\rho}{2}} \, \left(1-\rho\right)\, \Gamma\left(\frac{\rho}{2}+1\right)\, \frac{x^{\rho}}{\rho}\right)\tag{5}$
Figure (1) below illustrates formula (1) above for $p(x)$ (blue curve) seems to converge to $\log(2 \pi)$ (orange horizontal dashed line) as $x\to\infty$ which is consistent with the leading term in the explicit formula for $p_o(x)$ defined in formula (3) above.
Figure (1): Illustration of formula (1) for $p(x)$
The function $p(x)$ defined in formula (1) above is an even function of $x$, but Figure (1) above originally seemed to suggest to me the function
$$f(x)=\left\{\begin{array}{cc}
0 & x<0 \\
\frac{p(x)}{\log(2 \pi)} & x\ge 0 \\
\end{array}\right.\tag{6}$$
could potentially qualify as a CDF (Cumulative Distribution Function) related to the distribution of the primes. But I now suspect the oscillation in $p_\rho(x)$ (defined in formula (5) above) illustrated in figure (2) below causes $f'(x)$ to have an infinite number of $0$ crossings which also perhaps eventually leads to $f(x)$ having an infinite number of $1$ crossings.
Figure (2) below illustrates formula (5) above for $p_\rho(x)$ evaluated over the first $100$ pairs of non-trivial zeta-zeros which is largely dominated by the contribution of the first pair of non-trivial zeta-zeros. Note the oscillation in $p_\rho(x)$ seems to grow in magnitude as $x$ increases.
Figure (2): Illustration of formula (5) for $p_\rho(x)$
I've been struggling to define an accurate and simple asymptotic for $p(x)$ defined in formula (1) above, but the only accurate asymptotic I've come up with is $\log(2 \pi)+p_{-2n}(x)$ which is based on the residue formula for $p(x)$ illustrated in formula (3) above.
I'm investigating a potential pseudo-closed-form (meaning non-elementary functions are acceptable) for the function $p_{-2n}(x)$ defined in formula (4) above (which was the subject of my related Math StackExchange question) and have made some progress. I believe $p_{-2n}(x)$ defined in formula (4) above can also be evaluated as illustrated in formula (7) below where I use the series representation in formula (10) below to evaluate $g(x)$ and the pseudo-closed-forms in formulas (8) and (9b) below to evaluate $r(x)$ and $h(x)$, but this is not yet a closed-form much less simple.
$$p_{-2n}(x)=r(x)+h(x)+g(x)\tag{7}$$
$$r(x)=-\frac{1}{2 x^4} \left(3 \pi x^2 \text{1F1}^{(1,0,0)}\left(1,2,-\frac{\pi }{x^2}\right)-\pi ^2 \text{1F1}^{(1,0,0)}\left(2,3,-\frac{\pi }{x^2}\right)+x^4 \left(E_1\left(\frac{\pi }{x^2}\right)+\log \left(\frac{1}{x^2}\right)+2 \log (x)-2\right)+e^{-\frac{\pi }{x^2}} x^2 \left(x^2 (\gamma +2+\log (\pi ))+\left(4 \pi -2 x^2\right) \log (x)-2 \pi (\gamma -1+\log (\pi ))\right)\right)\tag{8}$$
$$h(x)=\sum\limits_{n=1}^\infty \frac{(-1)^n\, (2 n+1)}{n!}\, \pi^n \left(\gamma+\log(2 \pi)-H_{2 n}\right)\, x^{-2 n}\tag{9a}$$ $$=\frac{1}{2 x^2} e^{-\frac{\pi }{x^2}} \left(x^2 \, _1F_1^{(1,0,0)}\left(0,\frac{1}{2},\frac{\pi }{x^2}\right)-2 \pi \, _1F_1^{(1,0,0)}\left(0,\frac{3}{2},\frac{\pi }{x^2}\right)+e^{\frac{\pi }{x^2}} \left(2 \pi \, _1F_1^{(1,0,0)}\left(2,2,-\frac{\pi }{x^2}\right)-x^2 \, _1F_1^{(1,0,0)}\left(1,1,-\frac{\pi }{x^2}\right)\right)-2 \left(e^{\frac{\pi }{x^2}}-1\right) x^2 (\gamma +\log (2 \pi ))+\pi (-4 \gamma +6-4 \log (2 \pi ))\right)\tag{9b}.$$
$$g(x)=\underset{N\to\infty}{\text{lim}}\left(\sum\limits_{n=1}^N \frac{(-1)^{n+1}\, (2 n+1)}{n!}\, \pi^n\, \frac{\zeta'(2 n+1)}{\zeta(2 n+1)} x^{-2 n}\right)\tag{10}$$
Figure (3) below illustrates formula (6) for $f\left(\frac{1}{x}\right)$ above using formula (1) above for $p(x)$ evaluated at $K=1000$.
Figure (3): Illustration of formula (6) for $f\left(\frac{1}{x}\right)$ using formula (1) for $p(x)$
Figure (4) below illustrates formula (6) for $f\left(\frac{1}{x}\right)$ above using formula (3) above for $p_o(x)$ where formula (4) is evaluated over the first $100$ trivial zeta-zeros and formula (5) is evaluated over the first $100$ pairs of non-trivial zeta-zeros.
Figure (4): Illustration of formula (6) for $f\left(\frac{1}{x}\right)$ using formula (3) for $p_o(x)$
Figure (5) below illustrates formula (1) above for $p(x)$ evaluated at $K=100$ in blue and $K\in\{9,11,13,16,17\}$ in orange through brown where the blue horizontal dashed line is at $\log(2 \pi)$. Note the contribution of each prime-power comes along just in time to temporarily "boost support" of formula (1) for $p(x)$ until the contribution of the next prime-power. Omitting the contribution of a single prime-power or adding the contribution of a non-prime-power will cause a significant perturbation in the evaluation of formula (1) for $p(x)$.
Figure (5): Illustration of formula (1) for $p(x)$ evaluated at $K=100$ in blue and $K\in\{9,11,13,16,17\}$ in orange through brown
I'll note that formula (1) for $p(x)$ only converges when $K\gg x$, and when holding $K$ constant formula (1) for $p(x)$ converges to $-\psi(K)$ as $x\to\infty$ where
$$\psi(x)=\sum\limits_{n\le x} \Lambda(n)\tag{11}$$
is the second Chebyshev function.
