You have an obvious (very obvious in fact) flaw in your argument. If $M_{i} - m_{i} \leq 2|f(a)|$ and $\Delta x_{i} \leq \epsilon/(2|f(a)|)$ then we clearly have $$S(P) - s(P) \leq \sum_{i = 1}^{n}2|f(a)|\cdot\frac{\epsilon}{2|f(a)|} = n\epsilon$$ and normally when $\epsilon \to 0$ the variable $n \to \infty$ and therefore we can't say that $n\epsilon$ is small. You have removed the $n$ somehow from your calculation. If each term in a sum is less than $\epsilon$, it does not mean that the whole sum is also less than $\epsilon$.
However it is really very simple to modify your argument. When you have a sum like $\sum A_{i}B_{i}$ and you want to show that it is small (formally less than some arbitrary $\epsilon$) then you need to bound either each of $A_{i}$ or each of $B_{i}$ (but not both, you made the error of bounding both). Suppose we bound $B_{i}$ so that $B_{i} < \epsilon$ and the sum is less than $\epsilon\sum A_{i}$. So the next problem is not to bound $A_{i}$, but the whole sum $\sum A_{i}$ (and this should be made less than some constant value) and you are done. If you choose to bound $A_{i}$ (by $\epsilon$) then you need to bound the whole sum $\sum B_{i}$ (by some constant value).
In our case it is simpler to bound the $\Delta x_{i}$ by $||P||$ so that the difference $S(P) - s(P)$ is less than $||P||\sum (M_{i} - m_{i})$. Since the function is of bounded variation it is possible to bound the whole sum $\sum (M_{i} - m_{i})$ (as shown in answer by John Ma).
I hope you have understood the issue in your argument also the proper way to fix it.
If you know the standard result that a function of bounded variation can be expressed as a difference of two increasing functions then you can easily see that to establish the integrability of functions of bounded variation it is sufficient to establish the integrability of increasing functions. This is really really easy.
With your notation if $f$ is increasing then $M_{i} = f(x_{i}), m_{i} = f(x_{i - 1})$ and $\Delta x_{i} < ||P||$. So we can see that $$S(P) - s(P) \leq ||P||\sum (M_{i} - m_{i}) = ||P||\sum \{f(x_{i}) - f(x_{i - 1})\}$$ and the last sum on RHS telescopes to $f(b) - f(a)$ and hence $$S(P) - s(P) \leq ||P||\{f(b) - f(a)\}$$ Choosing $||P|| < \epsilon/\{f(b) - f(a) + 1\}$ we can see that $$S(P) - s(P) < \epsilon$$ and our job is done.
The only hard (but not too hard) part in this approach is proving the standard result mentioned in bold above.