output_description: Print one integer — the sum over ((A_i \cdot A_j) \bmod P).
  
  * * *
input_description: Input is given from Standard Input in the following format.
  
      
      
      N
      A_1 A_2 \cdots A_N
attempt: import numpy as np
  
  
  class FastFourierTransform:
      """FFT"""
  
      def multiplySameArr(self, A):
          """multiply by FFT(A*A)"""
          N = 1
          # 1. N > deg(g) + deg(h)となる2の冪乗を選ぶ
          while N < 2 * len(A):
              N *= 2
  
          # 2. 離散フーリエ変換の計算
          nf = np.zeros(N, int)
          nf[: len(A)] = A
          nf = np.fft.fft(nf)
          nf = nf * nf
  
          # 3. 離散フーリエ逆変換
          nf = np.fft.ifft(nf)
  
          # コピー
          ans = np.zeros(N, int)
          for n in range(N):
              # 係数を1/N倍 => np.fftの場合、すでに実行済
              # ans.append( round(nf[n].real / N) )
              if nf[n].real != 0:
                  ans[n] = int(nf[n].real + 0.5)
  
          return ans
  
  
  if __name__ == "__main__":
      N = int(input())
      A = list(map(int, input().split()))
  
      P = 200_003
      rootP = 2
  
      # order of 1 to P-1 after fft
      order = [1]
      x = rootP
      while x != 1:
          order.append(x)
          x *= rootP
          x %= P
  
      ORDER_SIZE = len(order)  # 200_002
  
      # create order index
      where = np.zeros(P, int)
      for i in range(ORDER_SIZE):
          where[order[i]] = i
  
      # count non zero items
      prevCnt = np.zeros(P, int)
      resCnt = np.zeros(P, int)
      for n in range(N):
          if A[n] != 0:
              prevCnt[where[A[n]]] += 1
              # remove Ax * Ax pattern count
              resCnt[(A[n] * A[n]) % P] -= 1
  
      # fft to cnt (calc nC2 and all Ax * Ay)
      clsFFT = FastFourierTransform()
      fftRes = clsFFT.multiplySameArr(prevCnt)
  
      for r in range(len(fftRes)):
          if fftRes[r] != 0:
              # get result(Ax*Ay==resCnt.index) Cnt by summary fftRes
              resCnt[order[r % ORDER_SIZE]] += fftRes[r]
  
      ans = 0
      for p in range(P):
          if resCnt[p] != 0:
              # remove j <= i pattern count
              # and
              # calc summary value((Ax * Ay cnt) * (Ax * Ay)))
              ans += (resCnt[p] // 2) * p
  
      print(ans)
  
problem_description: Statement
  
  Let’s take a prime P = 200\,003. You are given N integers A_1, A_2, \ldots,
  A_N. Find the sum of ((A_i \cdot A_j) \bmod P) over all N \cdot (N-1) / 2
  unordered pairs of elements (i < j).
  
  Please note that the sum isn't computed modulo P.
samples: [{"input": "4\n    2019 0 2020 200002", "output": "474287\n    \n\nThe non-zero products are:\n\n  * 2019 \\cdot 2020 \\bmod P = 78320\n  * 2019 \\cdot 200002 \\bmod P = 197984\n  * 2020 \\cdot 200002 \\bmod P = 197983\n\nSo the answer is 0 + 78320 + 197984 + 0 + 0 + 197983 = 474287.\n\n* * *"}, {"input": "5\n    1 1 2 2 100000", "output": "600013"}]