Finding pairs of triangular numbers whose sum and difference is triangular

There are infinitely many solutions. I'll first give infinitely many, and then talk about describing the set of all solutions.

Let $(a,u)$ be a solution to Pell's equation $$a^2-24u^2 =1$$ with $a$ and $u$ odd. All the integer solutions to this equation are of the form $a+u\sqrt{24} = (5+\sqrt{24})^n$ (using the methods in the Wikipedia link) and, when $n$ is odd, so are $a$ and $u$. The first few values of $a$ and $u$ are: $(5,1)$, $(485, 99)$, $(47525, 9701)$, $(4656965, 950599)$...

Then $(a,b,c,d) = (a, 5u, 7u, u)$ obeys $$a^2-1 = 24u^2 = c^2-b^2 = d^2-c^2$$ and $(a,b,c,d)$ are all odd.


Let $p=c/b$ and $q=d/b$, so $p^2+q^2=2$. There is a standard method for describing rational points on a conic once you know one point. In this case, $(1,1)$ is on conic. Consider the line through $(1,1)$ with slope $k/\ell$. Its other intersection with the conic is at $$\left( 1-\frac{2k(k+\ell)}{k^2+\ell^2}, 1-\frac{2k(k+\ell)}{k^2+\ell^2} \right) = \left( \frac{-k^2-2k\ell+\ell^2}{k^2+\ell^2}, \frac{k^2-2k\ell-\ell^2}{k^2+\ell^2} \right).$$ The line between two rational points always has rational slope, so every rational solution to $p^2+q^2=2$ is of the above form. Clearing denominators, we have $$(b,c,d) = u \cdot (k^2+\ell^2, -k^2-2k\ell+\ell^2, k^2-2k\ell-\ell^2).$$ In order to get $(b,c,d)$ all odd, take $k$ and $\ell$ of different parities and $u$ odd. I found $(5,7,1)$ by taking $(k,\ell) = (1,2)$. (There are also some solutions where $u$ is a half integer and $k$ and $\ell$ are both odd. I didn't check this carefully, since my goal was to find infinitely many solutions, not to classify them all, but if you want to throughly understand the problem you should consider the possibility that the three terms inside parentheses have a common factor.)

In any case, my strategy was to think about fixing $(k, \ell)$ and looking for $u$'s such that the common difference $$d^2-c^2 = c^2 - b^2 = u^2 \cdot 4 k \ell (k^2 - \ell^2)$$ would be of the form $a^2-1$. This is just Pell's equation $$a^2 - Nu^2=1$$ with $N=4 k \ell (k^2 - \ell^2)$ and the auxilliary condition that $a$ and $u$ are odd.

Pell's equation will be solvable whenever $N>0$ and not a square. There are certainly plenty of examples when $N$ is not a square. In fact, it will never be a square unless it is $0$. If $N=m^2$ were square, then $c^2=b^2+(mu)^2$ and $d^2 = (mu)^2 + c^2$. Fermat showed that there are no solutions to $x^2+y^2=z^2$ and $y^2+z^2=t^2$ with $y \neq 0$, see here.

But will there be odd solutions to Pell's equation? Parity considerations show that, since $N$ is even, $a$ will always be odd. If the primitive solution $a_0 + u_0 \sqrt{N}$ to Pell's equation has $u_0$ odd, then so will the solutions of the form $(a_0+u_0 \sqrt{N})^{2m+1}$. But if $u_0$ is even, then every solution to Pell's equation will have $u$ even.

I don't know a principle which predicts which of these two cases will occur but, as the example above shows, $u_0$ is odd at least some of the time, and that's all you need for infinitely many solutions.

Just to show off, here's one more example. I'll take $(k, \ell) = (3,4)$. This leads to the discovery that $$31^2 - 25^2 = 25^2-17^2$$ We have $N = 336$. The primitive solution to $a^2 - 336 u^2 =1$ is $55^2-3^2\cdot 336=1$. So there are infinitely many solutions of the form $$(a,b,c,d) = (a, 25u, 31u, 17u)$$ the smallest of which is $$(55, 75, 93, 51).$$