GEOXD - Editorial




Extended Euclid, Diophanite Equation Solver


Find the number of integer 3D points of the intersection between the plane ax + by + cz = d and the infinite prism x_1 \leq x \leq x_2, y_1\leq y \leq y_2.


Here problem is wrapping around with some geometry stuff but simply problem is to finding number of integer solution of equation ax + by + cz = d, where x_1 \leq x \leq x_2 and y_1 \leq y \leq y_2 and z is free. If we think about a box and a plane intersecting the box then here we need to find number of lattice points on the plane that is also inside the box. So what will be the condition for a point (x, y, z) to be inside a box whose 2 corner points are (x_1, y_1, z_1) and (x_2, y_2, z_2)? Condition will be x_1 \leq x \leq x_2, y_1 \leq y \leq y_2 and z_1 \leq z \leq z_2. Since this is a 3D stripe then, z_1 = -\infty and z_2 = +\infty. Now for satisfying the condition of point to be onto the plane also, so we basically need to find the number of integer tuple (x, y, z) such that ax + by + cz = d, x_1 \leq x \leq x_2 and y_1 \leq y \leq y_2. Now since ax + by + cz = d has 3 variables, here problem becomes a bit difficult. Same problem for diophanite equation is popular and classic problem.


So we can rewrite the given equation as follows:

gcd(a, c) \cdot f + by = d \rightarrow (i), where gcd(a, c) \cdot f = ax+cz \to (ii).
Here f is an integer.

For equation (i) we also know that y \in [y_1 , y_2], and for such y, we can also have a range of valid values of f. Let find out this below.

Let us find the values of f_g, y_g which satisfy the following:

gcd(a, c) f + by = gcd(a, b, c)

How to compute?

We can get this by using extended euclid algorithm

Then, one of the solution of equation (i) is

f_0 = f_g \frac{d}{gcd(a, b, c)} and y_0 = y_g \frac{d}{gcd(a, b, c)}.


If gcd(a, b, c) \nmid d, then there is no integer solution.

Now f_0 and g_0 is not the only solution. Other solution can be found by this arithmetic progression

f = f_0 + k \cdot \frac{b}{gcd( a, b, c )}
y = y_0 - k \cdot \frac{gcd( a, c )}{gcd( a, b, c )}.

Therefore, we get a range of values for y ( since y \in [y_1 : y_2]) as well as a range of solutions for f which will be an arithmetic progression \in [f_1, f_2] for some values f_1 and f_2. We can find f_1 and f_2 as well by the above calculation, which is not so difficult.

Suppose the first term of that arithmetic progression is a and difference between two terms is d and progression length is n. Now we can write the i th term of arithmetic progression:

T_i = u + i \cdot v, where f_1 \leq u + i \cdot v \leq f_2.

Now we have to solve our equation (ii) for every term of arithmetic progression.

Equation corresponding to T_i :

a \cdot x + c \cdot z = gcd(a, c)\cdot(u + i \cdot v) \rightarrow (iii)
where x_1 \leq x \leq x_2 and z is free.

Now problem is finding number of integer solutions of equation (iii) for all i.

Let x_g, z_g satisfy the following:

a \cdot x + c\cdot z = gcd(a, c)

This we can again find using extended euclid algorithm. Then, the solutions for equation (iii) will be:

x_0 = x_g\cdot (u + i \cdot v)
z_0 = z_g \cdot (u + i \cdot v)

Other solution can be found by the following formula:

x = x_0 + k \cdot \frac{c}{gcd(a, c)}
\implies x = x_g(u + iv) + k \cdot \frac{c}{gcd(a, c)}(replacing the value of x_0 = x_g(u + iv)).
\implies x = x_gu + x_gvi + k \cdot \frac{c}{gcd(a, c)}
\implies k = \frac{x - x_gu - x_gvi} {c / gcd(a, c )}

For i^{th} term, here lower bound is x_1. So minimum solution,

  • k_{min_i} = \lceil{\frac{x_1 - x_gu - x_gvi}{c / gcd(a, c)}}\rceil \rightarrow (iv)

For i^{th} term, here upper bound is x_2. So maximum solution,

  • k_{max_i} = \lfloor{\frac{x_2 - x_gu - x_gvi}{c / gcd(a, c)}}\rfloor \rightarrow (v)

So, number of solution of i^{th} term is (k_{max_i} - k_{min_i} + 1)

Now we need to sum this:

\sum{(k_{min_i} - k_{max_i} + 1)} = \sum{k_{min_i}} - \sum{k_{max_i}} + \sum 1

Now we will calculate \sum{k_{min_i}} and \sum{k_{max_i}} separately.

  • \sum{k_{min_i}} and \sum{k_{max_i}} can be calculated in similar way.
  • Here first observation is in equation (iv) numerator (x_1 - x_gu - x_gvi) is an arithmetic progression, where first term, p = x_1 - x_gu and difference between 2 term, q = -x_gv.
  • Now problem coverts to calculating the floor sum of an arithmetic progression.

There can be at most 5 \cdot 10^5 terms in subtask1 and 10^9 terms in subtask2. So looping through all terms will be too slow. Let p is the first term of the arithmetic progression, q is the difference between 2 terms and n is the length of the progression.

Now we need to calculate floor sum:

\sum_{i = 0}^{n - 1} \lfloor \frac{p + iq}{r} \rfloor
\implies \sum_{i = 0}^{n - 1} \frac{(p + iq) - (p + iq)\ mod\ r}{r}
\implies \frac{1}{r} \sum_{i = 0}^{n - 1} (p + iq) - \frac{1}{r} \sum_{i = 0}^{n - 1} (p + iq)\ mod\ r
\implies \frac{\frac{n}{2}(2p + (n - 1)q)}{r} - \frac{1}{r} \sum_{i = 0}^{n - 1} (p + iq)\ mod\ r

Now, here problem is to finding, S = \sum_{i = 0}^{n - 1} (p + iq)\ mod\ r

Let’s try to solve this for subtask 1.

We will break the AP (aritmetic progression) into sq = \lfloor \sqrt n \rfloor number blocks of sq length and try to sum up each blocks. Let’s take a view on i^{th}(0-based) block of AP:

  • (p + sq \cdot i \cdot q), (p + (sq \cdot i + 1)q), (p + (sq \cdot i + 2)q).....(p + (sq \cdot i + sq - 1)q)

Here j^{th}(0-based) term of i^{th} block is

  • p + (sq \cdot i + j)q = p + jq + (sq \cdot i \cdot q).

So here basically by adding (sq \cdot i \cdot q) to the j^{th} term of 0^{th} block we can get j^{th} term of i ^{th} block. Now we will take mod\ r into consideration.

Lets sum up the i^{th} block, that means:

  • \sum_{j = 0}^{sq - 1}(p + jq + (sq \cdot i \cdot q))\ mod\ r
    \implies \sum_{j = 0}^{sq - 1}((p + jq)\ mod\ r + (sq \cdot i \cdot q)\ mod\ r)\ mod\ r

Now we will use this observation that by adding (sq \cdot i \cdot q) to the j^{th} term of 0^{th} block we can get j^{th} term of i ^{th} block.

Let’s sort the 0^{th} block on the basis of (p + jq)\ mod\ r. Let \{s_1, s_2, s_3, ...., s_{sq - 1}\} is the sorted sequence, where s_j < r.

To sum up i^{th} block we need to calculate:

  • \sum_{j = 0}^{sq - 1}(s_j + (sq \cdot i \cdot q)\ mod\ r)\ mod\ r
    \implies \sum_{j = 0}^{sq - 1}(s_j + C)\ mod\ r (let C = (sq \cdot i \cdot q)\ mod\ r, here C < r)
    \implies \sum_{j = 0}^{l}(s_j + C)\ mod\ r + \sum_{j = l + 1}^{sq - 1}(s_j + C)\ mod\ r (here sequence is broken into two parts. In first part s_j + C < r where 0 \leq j \leq l and in second part r \leq s_j + C < 2r, where l + 1 \leq j \leq sq - 1)
    \implies \sum_{j = 0}^{l}(s_j + C) + \sum_{j = l + 1}^{sq - 1}(s_j + C - r)
    \implies \sum_{j = 0}^{l}(s_j + C) + \sum_{j = l + 1}^{sq - 1}(s_j + C) - \sum_{j = l + 1}^{sq - 1}r
    \implies \sum_{j = 0}^{sq - 1}(s_j + C) - r(sq - l - 1)

So, by doing a binary search we can get our l. So to solve each block we need to do a binary search on a sq length sorted sequence. There are \frac{n}{sq} = sq blocks. Complexity of solving each test case is O(\sqrt n log_2(\sqrt n)). So total time complexity for subtask1 is O(T \cdot \sqrt n log_2(\sqrt n))

Solving subtask 2

Now for subtask 2 we will try to find floor sum of AP in log_2(n), which is described briefly in the following blog with all kind of variants. So total time complexity of subtask2 is O(T \cdot log_2(n))(constant factor is too high)

In this solution we need to be concerned about some overflow issues because the first term of the AP might not be fit in long long int. So we may need to use 128 bit data type. And also we need to concern about floor sum of AP when terms are both negative and positive. We need to handle this negative and positive terms separately(This depends on how you implement the solution)


Even though I have only solved first subtask, looks like there can be another completely different approach for a full solution.

Based on this paper ( – pages 7-9 – it appears that all solutions to ax+by+cz=d equation are given by:
x = m1t + n1s + k1, y = m2t + n2s + k2, where m1,m2,n1,n2,k1,k2 are all integers that can be derived from a,b,c,d in O(log max(a,b,c,d)). (except some corner cases, where there are no solutions at all. This paper also gives a way to find out whether there’s at least 1 solution)

At this point, our problem reduces to the following (let’s replace t and s by X and Y, respectively:
How many integers X, Y exist, such that: x1 <= m1X + n1Y + k1 <= x2, y1 <= m2X + n2Y + k2 <= y2.
If you imagine points (X,Y) in a 2D plane, all real solutions to these 2 inequalities are inside a parallelogram, whose coordinates are easy to find in O(1). So we reduced our problem to the following: Find number of lattice points inside a parallelogram with rational coordinates.

Let’s split this parallelogram into 2 triangles by one of its diagonals, now we can split the task in 2 parts: (1) Find count of lattice points inside a triangle with rational coordinates of its vertices; (2) Find number of lattice points on a segment connecting 2 rational points;

Both of these can found in this paper - Idea is that triangles can be represented as rectangles with parallel sides to OX and OY, minus some other right triangles with 2 sides also parallel to OX and OY. The right triangle with OX, OY parallel sides can be represented as a “half” of a rectangle, where both sides are parallel to X and Y axis.

I omitted some details that these papers specify, but all of these methods work in O(log max_parameter) time. (Please let me know if I’m missing anything, or if anything I wrote is wrong)


By the way, using this link you can reduce the problem to two double sided (in)equations in two variables. Those aren’t trivial by any means, but they are certainly much easier than battling out all of the number theory.

How would you proceed from there?? I also made those in-equations in 2 variables…you get a polygon in x-y plane(with rational coordinates) and you need to find lattice points inside it…

I have seen this formula. But counting number of points from this equation also requires the same procedure that I have described above. Basically here finding a solution is not so difficult but finding number of solutions makes it challenging.

Basically the polygon in this case is a parallelogram. You can divide it into two triangles and proceed as in the editorial. Or else, you could attempt to solve it directly using a generalization of Pick’s theorem called Ehrhart’s quasipolynomial but that is too abstract so I couldn’t do it that way (I read of it in some obscure paper of 2003 or something).

No that polygon is not a parallelogram, you’ll get 6 in-equations(2 from each variables x, y, z)…so it is hexagon or can be pentagon in some cases…You can consider the equation:
3x+5y+7z = 18, and plot the area bounded by the inequations…

Does this work in 3D co-ordinate system?

@querty1234 Nope, you can ignore the one in z since z can take any integer value i.e. we dont care about z. Thus parallelogram. @sajib_readd no, its pick’s theorem’s analogue for rational coordinates.

Can we please have the setter’s and tester’s solutions (source codes)? Thank you.

I wonder is it possible to pursue a solution based on series? As in, the answer should be equal to the coefficient of t^d in F(t) = (t^{x_1} + t^{x_1 + 1} + ... + t ^{x_2})^a.(t^{y_1} + t^{y_1 + 1} + ... + t^{y_2})^b.(1 + t + t^2 + t^3 + ...)^c however, obviously this only captures non-negative values of z. For negative ones we can move the last term in the denominator.

Has anyone explored solving the problem from that angle?

1 Like

I had tried to solve it this way but failed badly. The series become too complicated as they are not infinite. Tried taking log and derivative but all in vain.
Also, you can determine the lowest and highest possible values of z by plugging in the limits of x and y in the plane’s equation but that only complicates the things.

1 Like