#include <iostream>
#include <vector>
#include <map>
#include <set>
#include <queue>
#include <algorithm>
#include <string>
#include <cmath>
#include <cstdio>
#include <iomanip>
#include <fstream>
#include <cassert>
#include <cstring>
#include <unordered_set>
#include <unordered_map>
#include <numeric>
#include <ctime>
#include <bitset>
#include <complex>
#include <chrono>
#include <random>
#include <functional>

using namespace std;

#define int long long

const int MOD = 200003;
const int P = MOD;
const int G = 2;

namespace Fft {
	typedef double ld;
	typedef complex<ld> cd;

	void dft(vector<cd> &a, int n) {
		vector<int> rev(1 << n);
		for (int i = 1; i < (1 << n); i++) {
			rev[i] = ((i & 1) << (n - 1));
			rev[i] |= (rev[i >> 1] >> 1);
		}
		for (int i = 0; i < (1 << n); i++) {
			if (i < rev[i]) {
				swap(a[i], a[rev[i]]);
			}
		}
		for (int step = 0; step < n; step++) {
			ld alpha = M_PI / (1 << step);
			cd w(cos(alpha), sin(alpha));
			for (int i = 0; i < (1 << n); i += (1 << (step + 1))) {
				cd curw(1, 0);
				for (int j = 0; j < (1 << step); j++) {
					cd u = a[i + j];
					cd v = curw * a[i + j + (1 << step)];
					a[i + j] = u + v;
					a[i + j + (1 << step)] = u - v;
					curw *= w;
				}
			}
		}
	}

	vector<int> mul(const vector<int> &a, const vector<int> &b) {
		int to = (int)max(a.size(), b.size());
		while (__builtin_popcount(to) > 1 || to < (int)(a.size() + b.size())) {
			to++;
		}
		int p = 0;
		while ((1 << p) < to) {
			p++;
		}
		vector<cd> f(to);
		for (int i = 0; i < to; i++) {
			if (i < (int)a.size()) {
				f[i] += cd(a[i], 0);
			} 
			if (i < (int)b.size()) {
				f[i] += cd(0, b[i]);
			}
		}
		dft(f, p);
		for (int i = 0; i < to; i++) {
			f[i] *= f[i];
		}
		dft(f, p);
		reverse(f.begin() + 1, f.end());
		vector<int> ans;
		for (auto &t : f) {
			t /= to;
			ans.push_back(round(t.imag() / 2));
		}
		return ans;
	}
} //namespace Fft

signed main() {
	ios_base::sync_with_stdio(false);
	cin.tie(0);

	int n;
	cin >> n;
	vector<int> a(n);
	for (int i = 0; i < n; i++) {
		cin >> a[i];
	}
	vector<int> wr(P, -1);
	vector<int> pw(P);
	{
		int c = 1;
		for (int i = 0; i < P; i++) {
			wr[c] = i;
			pw[i] = c;
			c = (c * G) % P;
		}
	}
	vector<int> cnt(P);
	for (int i = 0; i < n; i++) {
		if (a[i] > 0) {
			cnt[wr[a[i]]]++;
		}
	}
	vector<int> res = Fft::mul(cnt, cnt);
	int ans = 0;
	{
		for (int i = 0; i < 2 * P; i++) {
			ans += (pw[i % (P - 1)] * res[i]);
		}
	}
	for (int i = 0; i < n; i++) {
		ans -= (a[i] * a[i]) % MOD;
	}
	ans /= 2;
	cout << ans << '\n';
}

### Prompt

Please provide a CPP coded solution to the problem described below:

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.

Constraints

* 2 \leq N \leq 200\,000
* 0 \leq A_i < P = 200\,003
* All values in input are integers.

Input

Input is given from Standard Input in the following format.


N
A_1 A_2 \cdots A_N


Output

Print one integer — the sum over ((A_i \cdot A_j) \bmod P).

Examples

Input

4
2019 0 2020 200002


Output

474287


Input

5
1 1 2 2 100000


Output

600013

### Response

```cpp
#include <iostream>
#include <vector>
#include <map>
#include <set>
#include <queue>
#include <algorithm>
#include <string>
#include <cmath>
#include <cstdio>
#include <iomanip>
#include <fstream>
#include <cassert>
#include <cstring>
#include <unordered_set>
#include <unordered_map>
#include <numeric>
#include <ctime>
#include <bitset>
#include <complex>
#include <chrono>
#include <random>
#include <functional>

using namespace std;

#define int long long

const int MOD = 200003;
const int P = MOD;
const int G = 2;

namespace Fft {
	typedef double ld;
	typedef complex<ld> cd;

	void dft(vector<cd> &a, int n) {
		vector<int> rev(1 << n);
		for (int i = 1; i < (1 << n); i++) {
			rev[i] = ((i & 1) << (n - 1));
			rev[i] |= (rev[i >> 1] >> 1);
		}
		for (int i = 0; i < (1 << n); i++) {
			if (i < rev[i]) {
				swap(a[i], a[rev[i]]);
			}
		}
		for (int step = 0; step < n; step++) {
			ld alpha = M_PI / (1 << step);
			cd w(cos(alpha), sin(alpha));
			for (int i = 0; i < (1 << n); i += (1 << (step + 1))) {
				cd curw(1, 0);
				for (int j = 0; j < (1 << step); j++) {
					cd u = a[i + j];
					cd v = curw * a[i + j + (1 << step)];
					a[i + j] = u + v;
					a[i + j + (1 << step)] = u - v;
					curw *= w;
				}
			}
		}
	}

	vector<int> mul(const vector<int> &a, const vector<int> &b) {
		int to = (int)max(a.size(), b.size());
		while (__builtin_popcount(to) > 1 || to < (int)(a.size() + b.size())) {
			to++;
		}
		int p = 0;
		while ((1 << p) < to) {
			p++;
		}
		vector<cd> f(to);
		for (int i = 0; i < to; i++) {
			if (i < (int)a.size()) {
				f[i] += cd(a[i], 0);
			} 
			if (i < (int)b.size()) {
				f[i] += cd(0, b[i]);
			}
		}
		dft(f, p);
		for (int i = 0; i < to; i++) {
			f[i] *= f[i];
		}
		dft(f, p);
		reverse(f.begin() + 1, f.end());
		vector<int> ans;
		for (auto &t : f) {
			t /= to;
			ans.push_back(round(t.imag() / 2));
		}
		return ans;
	}
} //namespace Fft

signed main() {
	ios_base::sync_with_stdio(false);
	cin.tie(0);

	int n;
	cin >> n;
	vector<int> a(n);
	for (int i = 0; i < n; i++) {
		cin >> a[i];
	}
	vector<int> wr(P, -1);
	vector<int> pw(P);
	{
		int c = 1;
		for (int i = 0; i < P; i++) {
			wr[c] = i;
			pw[i] = c;
			c = (c * G) % P;
		}
	}
	vector<int> cnt(P);
	for (int i = 0; i < n; i++) {
		if (a[i] > 0) {
			cnt[wr[a[i]]]++;
		}
	}
	vector<int> res = Fft::mul(cnt, cnt);
	int ans = 0;
	{
		for (int i = 0; i < 2 * P; i++) {
			ans += (pw[i % (P - 1)] * res[i]);
		}
	}
	for (int i = 0; i < n; i++) {
		ans -= (a[i] * a[i]) % MOD;
	}
	ans /= 2;
	cout << ans << '\n';
}

```