Good evening,

I'm about to write my own quadratic sieve implementation in C using GMP library for large numbers. I'm facing an issue while attempting to do the last factorization step for the number: (I try to validate an example from a lecture)

$$\begin{equation} n = 9788111 (3571*2741) \end{equation}$$

My book tells me how to calculate a solution. This congruence will only lead to a primitive solution:

$$\begin{equation} 3131^2*3174^2*3481^2 \equiv (2*5^3*7*11^3*43)^2 = 100157750^2 \end{equation}$$

$$\begin{equation} gcd(9788111, 3131*3174*3481-100157750) = 9788111 \end{equation}$$

My current implementation looks like this and it's always giving correct solutions for primitive factors.

n = 9788111

size_x = 1
size_y = 1
for smooth_nr in sieved_nrs:
size_x *= smooth_nr
size_x %= n  // cut off

size_y *= smooth_nr * smooth_nr - n  // equation x^2 - n // can't reduce

size_y = sqrt(size_y)
size_y = size_y % N

gcd(n, size_x - size_y)  // will yield a solution


But now when it comes to non-primitive solutions this algorithm will fail. For example one factor is:

$$\begin{equation} 3129^2 * 3313^2 * 3449^2 * 4426^2*4651^2 \equiv (2^2*3*5*7^2*11*17*23^2*31^2*47)^2 = 13136082114540^2 \end{equation}$$

$$\begin{equation} gcd(9788111, 3129*3313*3449*4426*4651 - 13136082114540) = 2741 \end{equation}$$

The problem right now is that the implementation will not work only for both congruences that will lead to factors. I've tried it manually several times and the square root was always wrong, using python it returned a floating point value.

So my question is: How can I modify the factorization step in order to yield the correct solution? Can you see any mistakes made here?

Thank you very much!