Maximum integer not in $\{ ax+by : \gcd(a,b) = 1 \land a,b \ge 0 \}$

With a heuristic argument, I can show that, if $z$ is an integer with $z \gtrsim \zeta(2)\,xy\,\ln(xy)$, then $z$ can be written as $ax+by$ with $a,b\in\mathbb{Z}_{\geq 0}$ such that $\gcd(a,b)=1$. This would mean that $$xy-1\leq m(x,y)\lesssim \zeta(2)\,xy\,\ln(xy)\,,$$ if $m(x,y)$ is the largest integer not in $\big\{ax+by\,\big|\,a,b\in\mathbb{Z}_{\geq 0}\,\text{ and }\gcd(a,b)=1\big\}$. (Here, $\zeta$ is the Riemann zeta-function.)

Note that integer solutions $(a,b)$ to the equation $ax+by=z$ are of the form $$(a,b)=\left(a_0+ky,b_0-kx\right)$$ for some $a_0,b_0\in\mathbb{Z}$, where $k\in\mathbb{Z}$ is arbitrary. Of these, there are approximately $\frac{z}{xy}$ solutions $(a,b)$ with $a,b\geq 0$.

Consider a prime $p\nmid \gcd(xy,z)$. The probability that $p$ divides $a_0+ky$ is $\frac{1}{p}$ and the probability that $p$ divides $b_0-kx$ is $\frac{1}{p}$. Hence, approximately, the probability that $p$ divides both $a_0+ky$ and $b_0-kx$ is $\frac{1}{p^2}$. (Note that the event that $a_0+ky$ is divisible by $p$ and the event that $b_0-kx$ is divisible by $p$ are not quite independent; hence, this argument is not valid, but as I said, this is a heuristic argument.)

Consider a prime $p\mid \gcd(xy,z)$. Then, either $p\mid x$ or $p\mid y$, but not both. Without loss of generality, assume that $p\mid x$. Then, as $p\mid z$, we have $p\mid b_0$ and $p$ always divides $b_0-kx$. The probability that $p$ also divides $a_0+ky$ is $\frac{1}{p}$.

Consequently, among nonnegative integer solutions $(a,b)$ to $ax+by=z$, approximately $$\frac{z}{xy}\,\prod_{\substack{{p\text{ prime}}\\{p\nmid \gcd(xy,z)}}}\,\left(1-\frac{1}{p^2}\right)\,\prod_{\substack{{p\text{ prime}}\\{p\mid \gcd(xy,z)}}}\,\left(1-\frac{1}{p}\right)=\frac{z}{xy}\,\frac{\displaystyle\prod_{p\text{ prime}}\,\left(1-\frac{1}{p^2}\right)}{\displaystyle \prod_{\substack{{p\text{ prime}}\\{p\mid \gcd(xy,z)}}}\,\left(1+\frac{1}{p}\right)}$$ pairs have the property that $\gcd(a,b)=1$. Hence, we need $$\frac{z}{xy}\,\frac{\displaystyle\prod_{p\text{ prime}}\,\left(1-\frac{1}{p^2}\right)}{\displaystyle \prod_{\substack{{p\text{ prime}}\\{p\mid \gcd(xy,z)}}}\,\left(1+\frac{1}{p}\right)}\gtrsim 1\,.$$ However, note that $$\prod_{\substack{{p\text{ prime}}\\{p\mid \gcd(xy,z)}}}\,\left(1+\frac{1}{p}\right)\leq\sum_{j=1}^{\gcd(xy,z)}\,\frac{1}{j}\leq\sum_{j=1}^{xy}\,\frac{1}{j}\approx \ln(xy)$$ and $$\prod_{p\text{ prime}}\,\left(1-\frac{1}{p^2}\right)=\frac{1}{\zeta(2)}\,.$$ Since $z\gtrsim\zeta(2)\,xy\,\ln(xy)$, we obtain $$\frac{z}{xy}\,\frac{\displaystyle\prod_{p\text{ prime}}\,\left(1-\frac{1}{p^2}\right)}{\displaystyle \prod_{\substack{{p\text{ prime}}\\{p\mid \gcd(xy,z)}}}\,\left(1+\frac{1}{p}\right)}\gtrsim \frac{z}{\zeta(2)\,xy\,\ln(xy)}\gtrsim 1\,.$$

Here is a plot to illustrate that this upper bound of $m(x,y)$ works well for most values of $x$ and $y$. In this plot, $x,y\in\{1,2,\ldots,100\}$. The magenta curve is given by $m(x,y)=\zeta(2)\,xy\,\ln(xy)$, the red curve is given by $m(x,y)=xy\,\ln(xy)$, and the cyan line is the lower bound $m(x,y)=xy-1$. The scattered blue plot is the actual $xy$-versus-$m(x,y)$ plot. Hence, I now believe that $\displaystyle\limsup_{xy\to\infty}\,\frac{m(x,y)}{xy\,\ln(xy)}$ is very close to $\zeta(2)$.

enter image description here

For $x,y\in\{1,2,\ldots,100\}$ with $x\leq y$, there are only $18$ pairs with $m(x,y)>\zeta(2)\,xy\,\ln(xy)$. For the sake of convenience, write $f(x,y)$ for $\zeta(2)\,xy\,\ln(xy)$. These pairs $(x,y)$ along with the associated $m(x,y)$ are listed below:
(1) $m(3,7)=110$ with $f(3,7)\approx 105.17$,
(2) $m(4,15)=462$ with $f(4,15)\approx 404.10$,
(3) $m(4,37)=1386$ with $f(4,37)\approx 1216.57$,
(4) $m(4,41)=1386$ with $f(4,41)\approx 1375.79$,
(5) $m(4,63)=2310$ with $f(4,63)\approx 2292.08$,
(6) $m(8,37)=2772$ with $f(8,37)\approx 2770.64$,
(7) $m(8,73)=6270$ with $f(8,73)\approx 6119.19$,
(8) $m(11,23)=2310$ with $f(11,23)\approx 2302.82$,
(9) $m(11,58)=6930$ with $f(11,58)\approx 6777.82$,
(10) $m(14,45)=7410$ with $f(14,45)\approx 6679.75$,
(11) $m(19,20)=3990$ with $f(19,20)\approx 3713.05$,
(12) $m(21,61)=15510$ with $f(21,61)\approx 15077.6$,
(13) $m(22,25)=6090$ with $f(22,25)\approx 5708.67$,
(14) $m(28,45)=14820$ with $f(28,45)\approx 14796.1$,
(15) $m(35,86)=40590$ with $f(35,86)\approx 39658.0$,
(16) $m(41,55)=30030$ with $f(41,55)\approx 28639.4$,
(17) $m(61,89)=87780$ with $f(61,89)\approx 76796.6$, and
(18) $m(67,89)=90090$ with $f(67,89)\approx 85270.6$.
The pair $(4,15)$ is the one most deviated from the $xy$-vs-$f(x,y)$ curve, followed by $(61,89)$ and $(14,45)$, respectively. I originally thought that these extraordinary pairs $(x,y)$ would have many prime divisors, especially small primes, but apparently, the data indicate that the prime divisors of $xy$ are irrelevant.