For a perfect square $n$, can we calculate integers $0I want to implement a fast algorithm (avoiding or at least minimizing bruteforce) which for a given square number $n$ calculates a series of positive integers $y$ up to a certain limit $l$ such that $n+y^2$ is again a perfect square?
I can obtain the trivial case $y=\frac{1}{2}(n-1)$ as described in the answer and comments on this question using the equation $(x-y)(x+y)=x^2-y^2$.
Example: Let us take $n=225(=15^2)$, we obtain one solution $y=\frac{1}{2}(225-1)=112$ directly. But how we can obtain the whole list $(8, 20, 36, 112)$?
My ideas: As described here I came across the so called "Norm-Form Equation", which is a diophantine equation of the form $x^2-Dy^2=n$ and we have the special case $D=1$. To deal with this Norm-Form Equation, in his book Richard A. Mollin involves the divisor function.
Refering to our example, if we can deduce the other $y$-values $8,20,36$ from the (directly obtained) value $y=112$, it would allow me to speed up my bruteforce search drastically.
Note: If we cannot avoid bruteforce completely, I would be greateful for an approach that at least minimizes the search for these $y$-values.
 A: For the given problem, you neither need Pell-like equations nor Pythagorean triples. You just need , as mentioned , the divisors of $\ n\ $ :
You want the solutions of $$n+y^2=z^2$$ with positive integers $\ y\ $ and $\ z\ $
This means $$n=(z-y)(z+y)$$
Hence every pair of positive integers $a,b$ with $ab=n$ , where $a,b$ have the same parity gives a solution : $y=\frac{a-b}{2}$ , $z=\frac{a+b}{2}$
To ensure that $y$ is positive, $a$ must exceed $\sqrt{n}$
In PARI/GP, the following routine determines the pairs of positive integer solutions :
gp > n=225;fordiv(n,s,t=n/s;if(s^2>n,if(Mod(s,2)==Mod(t,2),[y,z]=[(s-t)/2,(s+t)/2];print(y," ",z))))
8 17
20 25
36 39
112 113
gp >

A: We begin with a formula that will generate
$\space (n,y,z):(n+y^2=z^2\space$ given two parameters which can be found given any number in the triple.
$$ n=(m^2-k^2)^2\qquad y=2mk\qquad z=m^2+k^2$$
We solve the $\space n$-function for $\space k,\space$ in terms of $\space m \space$ and $\space n.\quad$ Then we test a defined range of $\space m$-values and, any that yield an integer $\space k\space$ indicate the $\space (m,k)\space$ values needed to generate a triple.
$$n=(m^2-k^2)^2\implies k=\sqrt{m^2-\sqrt{n}}\\
\text{for}\qquad  \sqrt{\sqrt{n}+1} 
\le m \le 
\frac{\sqrt{n}+1}{2}$$
The lower limit ensures $k\in\mathbb{N}$ and the upper limit ensures $m> k$.
$$n=15^2\implies \sqrt{15+1}=4\le m \le \frac{15+1}{2} =8\\\text{and we find}\quad m\in\{4,8\}\implies k \in\{1,7\} $$
$$F(4,1)=(225,8,17)\qquad F(8,7)=(225,112,113) $$
But what about the other two  $\space y$-values? We find a triple for each of the prime factors of
$\space \sqrt{n}\space$ and then multiply the results by the cofactor(s). $\quad 225=3^2\times5^2\space$ so
$$n=3^2\implies \sqrt{3+1}=2\le m \le \frac{3+1}{2} =2\\
\text{and we find}\quad m\in\{2\}\implies k \in\{1\} $$
$$F(2,1)=(3^2,4,5)\longrightarrow 5^2(3^2,4,5)=(225,20,25) $$
$$n=5^2\implies \sqrt{5+1}=3\le m \le \frac{5+1}{2} =3\\
\text{and we find}\quad m\in\{3\}\implies k \in\{2\} $$
$$F(3,2)=(5^2,12,13)\longrightarrow 
3^2(5^2,12,13)=(225,36,39) $$
