First, we note that $f(m,n) - f(m,n-1)$ does not depend on $n$, so for a fixed $m$ we can write $f(m,n) = \alpha_m n + \beta_m$ for some $\alpha_m,\beta_m$ which depend only on $m$. Imposing the boundary conditions lets us write
$$
f(m,n) = a_m (n-1) + (2^m-1)
$$
Next, substituting that into the defining equation and simplifying gives
$$
a_m = a_{m-1}(2^{m-1}-2) + 2^{m-1}
$$
For example, it is easy to see that $a_2 = 2$ and $a_3 = 8$, with
$$
f(2,n) = 2(n-1)+3\\
f(3,2) = 8(n-1)+7
$$
The usual trick for such recursions is to let (starting at $m=3$)
$$
a_m = x_m\prod_{i=2}^{m-1}(2^i-2)
$$
which gives $x_3 = a_3/(2^2-2) = 4$ and for $m>3$
$$
x_m = x_{m-1} + \frac{2^{m-1}}{\prod_{i=2}^{m-1}(2^i-2)}
\\ x_m = x_3 + \sum_{k=4}^m \frac{2^{k-1}}{\prod_{i=2}^{k-1}(2^i-2)}
$$
$x_m$ very rapidly approaches (from above) $\prod_{s=0}^m (1+2^{-s})$.
So $a_m$ approaches (from above) $$\frac13 \prod_{i=2}^m \frac{2^i-2}{1+2^{-i}}
\rightarrow \frac13\prod_{i=2}^m (2^i-3)$$.
THis not only verifies (and slightly tightens) your relation, but also gives a handle on a lower limit, which comes from taking only one $-3$ in just one term of the product (and summing over which term uses the $-3$).
(You can actually get $x_m$ and $a_m$ in closed form, but it involves Q-Pochammer symbols and elliptic functoins so that is not at all illuminating.)
