#include<cmath>
#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
const double pai=acos(-1.0);
typedef long long ll;
const int mod=200003,g=2;
const int N=1600010;
struct complex{
    double real,imaginary;
    complex(double x=0,double y=0):real(x),imaginary(y){}
    complex operator=(const double b){
        real=b;
        return *this;
    }
    complex operator+(const complex b){
        return complex(real+b.real,imaginary+b.imaginary);
    }
    complex operator-(const complex b){
        return complex(real-b.real,imaginary-b.imaginary);
    }
    complex operator*(const complex b){
        return complex(real*b.real-imaginary*b.imaginary,real*b.imaginary+imaginary*b.real);
    }
    complex operator*=(const complex b){
        *this=*this*b;
        return *this;
    }
};
int n,t[N],r[N];
complex f[N];
ll ans;
void FFT(complex *a,int limit,int type){
    for(int i=0;i<limit;i++)
        if(i<r[i])
            swap(a[i],a[r[i]]);
    for(int mid=1;mid<limit;mid<<=1){
        const double x=pai/(const int)mid;
        complex wn=complex(cos(x),type*sin(x));
        for(int R=mid<<1,j=0;j<limit;j+=R){
            complex w=complex(1,0);
            for(int k=0;k<mid;k++,w*=wn){
                complex y=a[k+j],z=w*a[k+j+mid];
                a[k+j]=y+z;
                a[k+j+mid]=y-z;
            }
        }
    }
    if(type==-1)
    	for(int i=0;i<limit;i++)
	    	f[i].real/=limit;
    return;
}
int main(){
	scanf("%d",&n);
	for(int i=g,l=1;i!=1;i=1ll*i*g%mod,l++)
		t[i]=l;
	t[1]=0;
	for(int i=1;i<=n;i++){
		int x;
		scanf("%d",&x);
		ans-=1ll*x*x%mod;
		if(x)
			f[t[x]].real=f[t[x]].real+1;
	}
	int limit=1;
	while(limit<=mod*2)
		limit<<=1;
	for(int i=0;i<limit;i++)
		r[i]=(r[i>>1]>>1)|((i&1)?(limit>>1):0);
	FFT(f,limit,1);
	for(int i=0;i<limit;i++)
		f[i]=f[i]*f[i];
	FFT(f,limit,-1);
	int p=1;
	for(int i=0;i<limit;i++){
		ans+=((ll)(f[i].real+0.5))*p;
		p=2ll*p%mod;
	}
	printf("%lld\n",ans/2);
	return 0;
}
### 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<cmath>
#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
const double pai=acos(-1.0);
typedef long long ll;
const int mod=200003,g=2;
const int N=1600010;
struct complex{
    double real,imaginary;
    complex(double x=0,double y=0):real(x),imaginary(y){}
    complex operator=(const double b){
        real=b;
        return *this;
    }
    complex operator+(const complex b){
        return complex(real+b.real,imaginary+b.imaginary);
    }
    complex operator-(const complex b){
        return complex(real-b.real,imaginary-b.imaginary);
    }
    complex operator*(const complex b){
        return complex(real*b.real-imaginary*b.imaginary,real*b.imaginary+imaginary*b.real);
    }
    complex operator*=(const complex b){
        *this=*this*b;
        return *this;
    }
};
int n,t[N],r[N];
complex f[N];
ll ans;
void FFT(complex *a,int limit,int type){
    for(int i=0;i<limit;i++)
        if(i<r[i])
            swap(a[i],a[r[i]]);
    for(int mid=1;mid<limit;mid<<=1){
        const double x=pai/(const int)mid;
        complex wn=complex(cos(x),type*sin(x));
        for(int R=mid<<1,j=0;j<limit;j+=R){
            complex w=complex(1,0);
            for(int k=0;k<mid;k++,w*=wn){
                complex y=a[k+j],z=w*a[k+j+mid];
                a[k+j]=y+z;
                a[k+j+mid]=y-z;
            }
        }
    }
    if(type==-1)
    	for(int i=0;i<limit;i++)
	    	f[i].real/=limit;
    return;
}
int main(){
	scanf("%d",&n);
	for(int i=g,l=1;i!=1;i=1ll*i*g%mod,l++)
		t[i]=l;
	t[1]=0;
	for(int i=1;i<=n;i++){
		int x;
		scanf("%d",&x);
		ans-=1ll*x*x%mod;
		if(x)
			f[t[x]].real=f[t[x]].real+1;
	}
	int limit=1;
	while(limit<=mod*2)
		limit<<=1;
	for(int i=0;i<limit;i++)
		r[i]=(r[i>>1]>>1)|((i&1)?(limit>>1):0);
	FFT(f,limit,1);
	for(int i=0;i<limit;i++)
		f[i]=f[i]*f[i];
	FFT(f,limit,-1);
	int p=1;
	for(int i=0;i<limit;i++){
		ans+=((ll)(f[i].real+0.5))*p;
		p=2ll*p%mod;
	}
	printf("%lld\n",ans/2);
	return 0;
}
```