Let me carry out matters using a complex-analytical approach, as Lucia suggests, and then say where the difficulty lies.

Let $0<\beta<\alpha\leq 1$. First of all, as Lucia says,
$$\sum_{m\leq x} \frac{1}{m^\alpha} \sum_{n\leq x/m} \frac{\log(x/mn)}{n^\beta} = \frac{1}{2πi} \int_{c-i\infty}^{c+i\infty} \zeta(s+\alpha) \zeta(s+\beta) \frac{x^s}{s^2} ds$$
for $c>1$. We shift the contour of integration to the left of $\Re(s)=0$, picking up the main terms $$\begin{aligned}&\frac{y^{1-\alpha}}{(1-\alpha)^2} \zeta(1-\alpha+\beta) + \frac{y^{1-\beta}}{(1-\beta)^2} \zeta(1-\beta+\alpha)\\ &+ \zeta(\alpha) \zeta(\beta) \log y + \zeta'(\alpha) \zeta(\beta) + \zeta(\alpha) \zeta'(\beta)\end{aligned}$$ plus an error term of size $O\left(y^{\frac{1}{2} - \frac{\alpha+\beta}{2}}\right)$ along the way. We are left with the task of estimating an error term $$\frac{1}{2πi} \int_R \zeta(s+\alpha) \zeta(s+\beta) \frac{x^s}{s^2} ds,$$ where the integral is over a contour $R$ of our choice going from $-r-i\infty$ to $-r+i\infty$, say, and satisfying $\Re s\leq -r$ at all points. The error will be clearly bounded by $O(K x^{-r})$, where
$$K = \frac{1}{2πi} \int_R \frac{|\zeta(s+\alpha)| |\zeta(s+\beta)|}{|s|^2} ds.$$
The problem does reduces to estimating $K$.

Now, there are rigorous-numerics packages that include integration and the possibility to compute the zeta function $\zeta(s)$. (I currently use ARB.) However, (a) computations must obviously be finite (at least assuming mortal mathematicians), and (b) computing $\zeta(s)$ is never a walk in the park, and rigorous integration only adds to the overhead. Integrating an expression such as above from $-1/2 - i T$ to $1/2 + i T$ takes 15 minutes for $T = 10000$ (says a better programmer than I), but we should not expect to go much further than $T = 100000$ programming casually on our laptops.

The problem that remains, then, is how to bound a tail $$\frac{1}{2\pi i} \left(\int_{-r-i\infty}^{-r-i T} + \int_{-r+i T}^{-r + i \infty} \frac{|\zeta(s+\alpha)| |\zeta(s+\beta)|}{|s|^2} ds\right).$$

The most obvious approach is to use Backlund's explicit bounds (1918) on $\zeta(\sigma + it)$ (see http://iml.univ-mrs.fr/~ramare/TME-EMT/Articles/Art06.html#Size). They are of the quality $$|\zeta(\sigma + i t)| = (1+o(1)) (t/2\pi)^{(1-\sigma)/2} \log t$$ for $0\leq \sigma\leq 1$ and $$|\zeta(\sigma + i t)| = (1+o(1)) (t/2\pi)^{1/2-\sigma} \log t$$ for $-1/2\leq \sigma\leq 0$. The problem here is that convergence is painfully slow. If, say, $\alpha = 1$, $\beta=1/2$ and $r =-1/4$ (reasonable values all around), the tails will be bounded by a constant times $(\log T)^2/\sqrt{T}$. For $T=10000$, $(\log T)^2/\sqrt{T} > 0.848\dotsc$ - not exactly small; for $T=100000$, the same equals $0.419\dotsc$ - barely an improvement.

Notice, however, that why Backlund's bounds are essentially tight for $\Re s<0$, that is not the case for $0<\Re(s)<1$. Of course, they are convexity bounds, so improving on them explicitly would amount to translating into explicit terms rather non-trivial material. However, as long as we are satisfied with $r>-\beta$, what we can do instead is give $L_2$ bounds for the tails, that is, bound
$$\int_{r-i \infty}^{r-i T} \frac{|\zeta(s)|^2}{|s|^2} ds + \int_{r+iT}^{r+i\infty} \frac{|\zeta(s)|^2}{|s|^2} ds.$$
(The integral $\int_{r-i \infty}^{r-i T}$ is obviously the same.) Then we use Cauchy-Schwarz to bound the tail of the integral we were discussing.

This is non-trivial, and takes us further afield, so I will make it into a separate question: $L_2$ bounds for tails of $\zeta(s)$ on a vertical line .

4more comments