複素数体 \(\mathbb{C}\) との対応を考えることで与えられた式を計算しやすい形に変形します。

\(\overline{z}\) で複素数 \(z\) の共役を表します。

\(\alpha_S = S_x + i S_y,\) \(\beta_i = x_i + iy_i\) とすると、 \(\begin{aligned} (S_x - x_i)^2 + (S_y - y_i)^2 = (\alpha_S - \beta_i)\overline{(\alpha_S - \beta_i)} = (\alpha_S - \beta_i)(\overline{\alpha_S} - \overline{\beta_i}) \end{aligned}\)

です。

これにより、 \(f(S)\) は多項式 \(g,h\) を用いて

\(\begin{aligned} f(S) = g(\alpha_S) h(\overline{\alpha_S}) \end{aligned}\)

と書けます。 ここで \(t\) を \(t ^ 2 = -1 \bmod{P}\) を満たす整数とします。 このような \(t\) は \(P = 1 \bmod{4}\) であることから存在します。

これを用いて

\(C = \{a + ib \ |\ 0 \le a,b < P, a,b : 整数\}\) と \(\mathbb{Z}/P\mathbb{Z}\) の対応 \(k\) を次のように与えます。

\(\begin{aligned} k(a + ib) = a + tb \end{aligned}\)

この対応を与えた後に \(g,h\) を計算しても値は変わりません。

ところで、\(P\) が素数であることから

\(C \ni z \longmapsto (k(z), k(\overline{z}))\)

は \(C\) と \((\mathbb{Z}/P\mathbb{Z})^2\) の間の全単射を与えます。

よって、求めるべき値は 整数 \(X,Y\) \((0 \le X,Y <P)\) の組の内、

\(g(X)h(Y) = Z\) を満たすものの数と言い換えられます。

これは各 \(0 \le i < P\) に対して \(h(Y) = i\) を満たす \(Y\) の数を求めておくことで計算できます。

計算量は \(\mathcal{O}(nP)\) です。