Hello, Codeforces. Recently I invented a problem. Consider a set of points (x, y) with integer coordinates, 0 <= x < a и 0 <= y < b. We want to know the number of acute triangles with vertices at given points.
Trying to interpolate
We can write function f(a, b) which finds the answer and works in (ab) ^ 3
. I suggested that it is a polynomial of two variables with degree <= 6. I tried to interpolate it using this theorem. But my code finds non-zero coefficients with monoms with degrees > 6, so my hypothesis was not confirmed. I also failed with right triangles.
code (for right triangles)int stupid(int a, int b)
{
int ans = 0;
for (int x1=0;x1<a;x1++)
for (int x2=0;x2<a;x2++)
for (int x3=0;x3<a;x3++)
for (int y1=0;y1<b;y1++)
for (int y2=0;y2<b;y2++)
for (int y3=0;y3<b;y3++)
{
int a = (x1 - x2) * (x1 - x2) + (y1 - y2) * (y1 - y2);
int b = (x1 - x3) * (x1 - x3) + (y1 - y3) * (y1 - y3);
int c = (x3 - x2) * (x3 - x2) + (y3 - y2) * (y3 - y2);
if ((a + b == c) && min(a, min(b, c))) ans++;
}
return ans / 2;
}
signed main()
{
int C1 = 5, C2 = 5;
ld res = 0;
for (int a1=4;a1<=C1+3;a1++)
{
for (int a2=4;a2<=C2+3;a2++)
{
int d = 1;
for (int ai=4;ai<=C1+3;ai++) if (ai != a1) d *= a1 - ai;
for (int ai=4;ai<=C2+3;ai++) if (ai != a2) d *= a2 - ai;
res += (ld)(stupid(a1, a2)) / d;
}
}
cout << setprecision(20) << fixed << res;
}
This code finds a coefficient with a ^ (C1 - 1) * b ^ (C2 - 1)
What I would like to know:
- can this problem be solved faster than O(stupid(a, b))
- can this problem be solved in O(1)
- maybe, someone knows problems with difficult to find formulas and this method can help to find them?
UPD: Found formula for the number of right triangles with b = 2: f(a, 2) = 2 * a ^ 2 - 4
, a > 1.
UPD2: Many thanks to bronze_coder for finding O(1)
solution for constant b: OEIS A189814.
Interpolation should use ai > b ^ 2
.
UPD3: Finally wrote a O((ab) ^ 2)
solution.
Codeint fast(int a, int b)
{
int ans = 0;
for (int a1=-a+1;a1<=a-1;a1++)
{
for (int b1=-b+1;b1<=b-1;b1++)
{
if (a1 == 0 && b1 == 0) continue;
for (int a2=-a+1;a2<=a-1;a2++)
{
for (int b2=-b+1;b2<=b-1;b2++)
{
if (a2 == 0 && b2 == 0) continue;
// first condition is that dot product = 0
// second condition is that cross product < 0 not to count the same pair twice
if (b1 * b2 + a1 * a2 == 0 && a1 * b2 - b1 * a2 < 0)
{
int cnta = a - max(max(abs(a1), abs(a2)), abs(a1 - a2));
int cntb = b - max(max(abs(b1), abs(b2)), abs(b1 - b2));
if (cnta > 0 && cntb > 0) ans += cnta * cntb;
}
}
}
}
}
return ans;
}
Now I can interpolate using bigger values of a
and b
.
But it is quite strange for me that the number of right triangles is directly proportional to (ab) ^ 2
instead of (ab) ^ 3
. Now I will try to understand why the formula doesn't work with ai <= b ^ 2
.
UPD4: Code that finds a formula for f(a, b) for a > b ^ 2
and works in O(b ^ 6)
:
Spoilerint fast(int a, int b)
{
int ans = 0;
for (int a1=-a+1;a1<=a-1;a1++)
for (int b1=-b+1;b1<=b-1;b1++)
{
if (a1 == 0 && b1 == 0) continue;
for (int a2=-a+1;a2<=a-1;a2++)
for (int b2=-b+1;b2<=b-1;b2++)
{
if (a2 == 0 && b2 == 0) continue;
if (b1 * b2 + a1 * a2 == 0 && a1 * b2 - b1 * a2 < 0)
{
int cnta = a - max(max(abs(a1), abs(a2)), abs(a1 - a2));
int cntb = b - max(max(abs(b1), abs(b2)), abs(b1 - b2));
if (cnta > 0 && cntb > 0) ans += cnta * cntb;
}
}
}
return ans;
}
struct polynom
{
int x0 = 0, x1 = 0, x2 = 0;
void print()
{
cout << x2 << "x^2 + " << x1 << "x - " << -x0 << endl;
}
ld get(ld x)
{
return x2 * x * x + x1 * x + x0;
}
};
void interpolate(int b)
{
polynom P = {0, 0, 0};
vector<int> A = {b * b + 1, b * b + 2, b * b + 3};
ld res = 0;
vector<ld> B = {fast(A[0], b), fast(A[1], b), fast(A[2], b)};
for (int i=0;i<3;i++)
{
int d = 1;
for (int j=0;j<3;j++) if (A[i] != A[j]) d *= A[i] - A[j];
res += (B[i] - P.get(A[i])) / d;
}
P.x2 = res;
res = 0;
A.pop_back();
for (int i=0;i<2;i++)
{
int d = 1;
for (int j=0;j<2;j++) if (A[i] != A[j]) d *= A[i] - A[j];
res += (B[i] - P.get(A[i])) / d;
}
P.x1 = res;
P.x0 = B[0] - P.get(A[0]);
P.print();
}
I used OEIS to find a regularity between the coefficients of P, but i didn't find anything :( except x2 = b * (b - 1)
, but it was obvious.
UPD5:
Finally found a formula and solution in O(min(a, b) ^ 6)
If a < b
, let's swap them. Now we need to solve it in O(b ^ 6)
. If a <= (b - 1) ^ 2 + 1
, we can run solution in O((ab) ^ 3)
= O(b ^ 6)
. Now we need to solve it for a > (b - 1) ^ 2 + 1
.
Definitions
Consider all right triangles and divide them into 3 types:
triangles without any vertices on the line x = a - 1
triangles with some vertices on the line x = a - 1
and with two sides which are parallel to the x-axis or y-axis
other triangles
Obviously, third type triangles have exactly one vertex on the line x = a - 1
. Let C(a)
be the number of third type triangles (some function).
The number of second type triangles is (a - 1) * b * (b - 1) / 2 * 4
:
ProofLet's choose a vertical line x = z, z < a - 1
(a - 1
variants), a vertical line x = a - 1
and two horizontal lines (b * (b - 1) / 2
variants). Intersection of these lines consists of 4 points. We can throw away any of them and the triangle formed from other 3 points has the second type. Obviously, every second type triangle can be presented as an intersection of these lines without one of points.
By definition, for every vertex (x, y) of first type triangles, 0 <= x < a - 1
and 0 <= y < b
, then the number of the triangles is f(a - 1, b)
, by definition of f
.
Well, f(a, b) = f(a - 1, b) + (a - 1) * b * (b - 1) / 2 * 4 + C(a)
Let's prove that C(a)
for every a > (b - 1) ^ 2 + 1
has same values.
ProofLemma. Maximal absolute value of difference between abscissas of two vertices of a third type triangle does not exceed (b - 1) ^ 2 + 1
Let's look on a picture
I could add it only rotated( The triangle has at least two sides which are parallel to neither x-axis nor y-axis (AB and BC). Let k be the slope of line AB, then the slope of BC is equal to -1/k
. Then, the difference between abscissas of A and C does not exceed
(b - 1) / |k| + (b - 1) / |-1/k| = (b - 1) * (|k| + 1 / |k|)
.
It is a well-known fact that if we have two positive numbers and we will bring their values together, then their sum will decrease. Then, to maximize |k| + 1 / |k|
, we need minimize or maximize |k|. According to limitation of coordinates, |k| < b - 1
and 1 / |k| < b - 1
. Then the maximum value is (b - 1) * ((b - 1) + 1 / (b - 1)) = (b - 1) ^ 2 + 1
. The Lemma is proved.
We need to prove that for a > (b - 1) ^ 2 + 1
the number of third type triangles does not depend on a. It is obvious after proving Lemma because we can compare every triangle for a
with a triangle for a + 1
by shifting it to right.
C(a)
is a constant. Let c
be this constant.
Well, f(a, b) = f(a - 1, b) + (a - 1) * b * (b - 1) / 2 * 4 + c
. After some transforms, we get a formula for f(a, b)
using f((b - 1) ^ 2 + 1, b)
.
Implementationll fast(int a, int b)
{
ll ans = 0;
for (int a1=-a+1;a1<=a-1;a1++)
for (int b1=-b+1;b1<=b-1;b1++)
{
if (a1 == 0 && b1 == 0) continue;
for (int a2=-a+1;a2<=a-1;a2++)
for (int b2=-b+1;b2<=b-1;b2++)
{
if (a2 == 0 && b2 == 0) continue;
if (b1 * b2 + a1 * a2 == 0 && a1 * b2 - b1 * a2 < 0)
{
int cnta = a - max(max(abs(a1), abs(a2)), abs(a1 - a2));
int cntb = b - max(max(abs(b1), abs(b2)), abs(b1 - b2));
if (cnta > 0 && cntb > 0) ans += cnta * cntb;
}
}
}
return ans;
}
struct polynom
{
ll x0 = 0, x1 = 0, x2 = 0;
void print()
{
cout << x2 << "x^2 + " << x1 << "x - " << -x0 << endl;
}
ld get(ld x)
{
return x2 * x * x + x1 * x + x0;
}
};
ll get(int a, int b)
{
if (a < b) swap(a, b);
if (a <= (b - 1) * (b - 1) + 1) return fast(a, b);
int x = fast((b - 1) * (b - 1) + 1, b);
int y = fast((b - 1) * (b - 1) + 2, b);
// y = x + c + (a - 1) * b * (b - 1) / 2 * 4;
int c = y - x - ((b - 1) * (b - 1) + 1) * (b - 1) * b * 2;
polynom P;
P.x2 = b * b - b;
P.x1 = c - b * b + b;
P.x0 =-b*b*b*b*b*b + 5*b*b*b*b*b - 11*b*b*b*b + 13*b*b*b + b*b*(-8 - c) + b * (2 * c + 2) - 2 * c + x;
return P.get(a);
}
Unfortunately, the value of c
is different for different values of b
and I could not find a regularity between them. Interpolation and OEIS did not help me. There are some things I should do:
- Express
c
through b
, after all c
depends on b
only - Solve the problem for
a <= (b - 1) ^ 2
faster - Solve the problem for acute triangles.
It will be hard.