CF736 D2

Compared with easy Removed the condition that the coordinates are even , So now S The area of should also be considered .

S = a + b/2 - 1

2S = 2a + b - 2

a Is odd ,S Is an integer

You can get what you need at this time 2S mod 4 = b mod 4

b According to the previous formula sigma gcd(abs(xi - xj), abs(yi - yj)) Ask for ,2S Then we can use the cross product to find

Consider and easy The same thing , Put all the coordinate modules 4 discretization , There's a problem , In seeking b When , When two numbers are subtracted 4 more than 1、 model 4 more than 3 Will transform each other ,b Among the feasible solutions are 1 + 3 + 2 And so on. , So the method is right b invalid , But for 2S Because as long as you ask for the surplus 2 or 4 Can still be established .

here , We enumerate more things to make sure it works , Open an array hav[i][nw][x][y], To express with i produce b by nw, Abscissa mod 4 = x And ordinate mod 4 = y Number of elements of , violence n^2logn seek

After that 1+3+2 And so on , Then enumeration 1 and 3 Adjacent edge of , be left over 2 The opposite side of is similar to easy The version can be judged .

#include <bits/stdc++.h>
#define pii pair<int,int>
#define fi first
#define sc second
#define pb push_back
#define ll long long
#define trav(v, x) for(auto v:x)
#define VI vector<int>
#define all(x) (x).begin(),(x).end()
using namespace std;
const double eps = 1e-10;
const int N = 6050;
const ll mod = 998244353;

pii pt[N];
int hav[N][4][4][4];

signed main()
{
	ios::sync_with_stdio(0);
	cin.tie(0);
	int n;
	cin >> n;
	for(int i = 1; i <= n; i++)
	{
		cin >> pt[i].fi >> pt[i].sc;
	}
	for(int i = 1; i <= n; i++)
	{
		for(int j = 1; j <= n; j++)
		{
			if(i == j)
				continue;
			int nw = __gcd(abs(pt[i].fi - pt[j].fi), abs(pt[i].sc - pt[j].sc)) % 4;
			hav[i][nw][pt[j].fi % 4][pt[j].sc % 4]++;
			//cerr << "!!" << i << ' ' << nw <<  ' ' << pt[j].fi % 4 << ' ' << pt[j].sc % 4 << '\n';
		}
	}
	ll ans = 0, res = 0;
	//case 1
	for(int i = 1; i <= n; i++)
	{
		int x1 = pt[i].fi % 4;
		int y1 = pt[i].sc % 4;
		for(int x2 = 0; x2 < 4; x2++)
		{
			for(int y2 = 0; y2 < 4; y2++)
			{
				for(int b1 = 1; b1 <= 3; b1 += 2)
				{
					ll num1 = hav[i][b1][x2][y2];
					for(int x3 = 0; x3 < 4; x3++)
					{
						for(int y3 = 0; y3 < 4; y3++)
						{
							for(int b2 = 1; b2 <= 3; b2 += 2)
							{
								ll tmp;
								int b = __gcd(abs(x2 - x3), abs(y2 - y3)) % 4;
								b = (b + b1 + b2) % 4;
								if(b != 0 && b != 2)
									continue;
								int SS = abs((x2 - x1) * (y3 - y1) - (y2 - y1) * (x3 - x1)) % 4;
								if(b != SS)
									continue;
								if(b1 == b2 && x3 == x2 && y3 == y2)
									tmp = num1 * (num1 - 1);
								else
									tmp = num1 * hav[i][b2][x3][y3];
								ans += tmp;
							}
						}
					}
				}
			}
		}
	}
	
	// case 2
	for(int i = 1; i <= n; i++)
	{
		int x1 = pt[i].fi % 4;
		int y1 = pt[i].sc % 4;
		for(int x2 = 0; x2 < 4; x2++)
		{
			for(int y2 = 0; y2 < 4; y2++)
			{
				for(int b1 = 0; b1 <= 2; b1 += 2)
				{
					ll num1 = hav[i][b1][x2][y2];
					for(int x3 = 0; x3 < 4; x3++)
					{
						for(int y3 = 0; y3 < 4; y3++)
						{
							for(int b2 = 0; b2 <= 2; b2 += 2)
							{
								ll tmp;
								int b = __gcd(abs(x2 - x3), abs(y2 - y3)) % 4;
								b = (b + b1 + b2) % 4;
								if(b != 0 && b != 2)
									continue;
								int SS = abs((x2 - x1) * (y3 - y1) - (y2 - y1) * (x3 - x1)) % 4;
								if(b != SS)
									continue;
								if(b1 == b2 && x3 == x2 && y3 == y2)
									tmp = num1 * (num1 - 1);
								else
									tmp = num1 * hav[i][b2][x3][y3];
								res += tmp;
							}
						}
					}
				}
			}
		}
	}
	//cerr << res << '\n';
	cout << (ans / 2) + (res / 6) << '\n';
}