PROBLEM LINK:
DIFFICULTY:
HARD
PREREQUISITES:
Sweep line algorithm, Computing square roots
PROBLEM:
You are given three points A, B and C and N pairs (x[i], y[i]). You can shift coordinates of A by at most K of these pairs. Each pair can be used at most once. What maximal value can get |AB| + |BC| + |CA| in this situation?
EXPLANATION:
Clearly |BC| is constant. So we should consider only |AB| + |AC| and add |BC| in the end. The key idea is to maximize some linear function f(X, Y) = U * X + V * Y instead of maximizing the value |AB| + |AC| that involves square roots. Here (X, Y) are the coordinates of the new position of A and U, V are some fixed parameters. Namely, the following statement holds.
Lemma. There exists some real numbers U and V such that the set of pairs, applied to the point A, that maximizes the function f(X, Y) = U * X + V * Y, also maximizes the sum of the distances |AB| + |AC|.
Proof. Consider the set of at most K pairs that maximizes the sum of distances |AB| + |AC|. Let this maximum be D and the corresponding point A, on which this maximum is achieved, be A’. So for each other set of at most K pairs this sum of distances is not greater than D. Note, that the set of points A for which |AB| + |AC| <= D is a set of points enclosed by an ellipse with foci at B and C that passes through the point A’. Consider the tangent line L to this ellipse at the point A’. Let its equation be U * X + V * Y + W. Since ellipse encloses the convex set and all the points A, corresponding to other sets of at most K pairs, lies inside this set, all these points A will also lie at the same half-plane with respect to the line L. Hence, we can assume, changing signs of U, V and W if necessary, that U * X + V * Y + W <= 0 for all other points A. And since for the point A’ we have equality here it exactly means that the point A’ maximizes the linear function f(X, Y) = U * X + V * Y. We can think on this from another position. Namely, the set of pairs that maximizes considered linear function also maximizes the sum of distances. So we are done.
As we see, the proof do not provide any hints on how to find this U and V. So in order to solve the problem we ideally should try all possible values of U and V, maximize the corresponding linear function and update the actual answer by the value |AB| + |AC| for the optimal point A of the current linear function. Of course, it’s unreal plan, so we should decrease the set of different pairs (U, V), that we will check, to some relatively small set of candidates.
For this consider the problem of maximization of linear function f(X, Y) = U * X + V * Y in detail. At first it is clear that we can remove all completely zero pairs (x[i], y[i]). So now we assume that every given pair is non-zero. Denote by F[i] the value U * x[i] + V * y[i]. If we apply to A the set of pairs with indexes i1, …, iL and obtain the point (X, Y), then by linearity of f we get f(X, Y) = f(Ax, Ay) + F[i<sub>1</sub>] + … + F[i<sub>L</sub>]. It is clear now, how to find the set of pairs that maximizes f(X, Y). If the number of positive values among F[1], …, F[N] is no more than K than choose all pairs with positive value, otherwise choose the K largest by F[i] pairs.
This observation shows that the only essential things that we should know for the given U and V are, which of F[i] are positive, and what their relative order. Clearly F[i] = 0 if and only if (U, V) is parallel to (y[i], -x[i]) and F[i] = F[j] if and only if (U, V) is parallel to (y[i] - y[j], x[j] - x[i]). If for some i and j we have equal pairs then F[i] = F[j] for all U and V and we shouldn’t consider direction produced by this pair. Lines parallel to these directions and passing through the origin divide the plane into several sectors. Clearly the signs of F[i] and their relative order are constant for all (U, V) in each particular sector. So we can choose just one particular pair (U, V) from each sector, find the optimal set of at most K pairs in O(N * log N) time (since we in general need to sort them) and update the actual answer, which is the sum of squares |AB| + |AC|. In this way we get the O(N3 * log N) solution since there are O(N2) sectors.
We can improve this solution to O(N2 * log N) solution. For this we need to maintain the signs and relative order of F[i] during iterations through the sectors without actual values of F[i], as well as, the optimal point.
Let’s sort all critical directions by polar angle (each critical line produce two opposite critical directions). This can be made by comparing at first half-planes of two directions and then by the sign of their cross-product. I suggest you to not use attractive atan2() function in C++ that calculates the polar angle, since in general situation it may leads to precision issues. In this problem it will be enough. But if coordinates of directions can be up to 109 then there can be two different directions that can’t be distinguished by atan2(). For example
**atan2(109, 109 + 1) = 0.78539816389744830936...**; **atan2(109 - 1, 109) = 0.78539816389744830986...**. As you see, they differ only at 19-th digit after the decimal point which is out of `double` precision.At the beginning we directly calculate relative order and signs distribution of F[i] for some direction (U, V) in the last sector (it is in fact the previous to the first one). For example one can take (U, V) = (2 * maxX + 1, -1), where maxX = 106 is the maximal possible absolute value of x[i] and y[i].
Next we consider critical directions one-by-one. If for the current direction F[i] = 0 then we change the corresponding sign. If for this direction F[i] = F[j] for some i and j we swap their positions. During both these transformations it is clear how the sum of optimal set of pairs is changed. But very special consideration is required when several events occur for the same direction. Currently, please, see tester’s solution as a reference. It currently has no explanatory comments. It will be improved soon, as well as the required explanation will be added here.
Now let’s discuss how to deal with real numbers in this problem with 13 digits after the decimal point. We are very disappointed that nobody figured it out the beautiful way we invent for this. All accepted solutions have cumbersome calculation of square root using some usual algorithm. Common guys! Are you really thought that we want you to add a ton of ugly code just to print crappy 13 digits after the decimal point? If we would really decide to see how you could calculate the square root we would ask to print, say, 100000 digits.
In this problem the only thing that is needed to deal with square roots with required precision is the following simple formula:
**sqrt(A) = M + (A - M * M) / (M + sqrt(A))**,where M = floor(sqrt(A)). Now the idea is to deal with real numbers here as with pairs (I, F), where I is an integer part and F is a fractional part, so 0 <= F < 1. It is clear how to compare and add such numbers (see tester’s solution as a reference). And the mentioned formula shows that for sqrt(A) we have I = M and F = (A - M * M) / (M + sqrt(A)). Here F should be calculated directly using built-in function sqrt() since we only need F with 13 significant digits.
You should use this not only for printing the answer but also for comparing current answer with some candidate answer. We have tests where very close different choices exist. The simplest such test is the following:
1 1
0 1000000000
1000000000 -1000000000
-1000000000 -1000000000
1 0
We can either use the spell card or not. The answers in these cases are:
6472135954.9995793928183
6472135954.9995793931761
So their relative error is 5.5279e-20 and can’t be distinguished by double.
Well, we even have several such tests involving large tests with many different sectors. So we had no doubt that solution with comparison in double will fail. But contestants as usually cheated on us. Check it out this solution and also this one.
Also in the first solution instead of smart iterating over critical directions simple heuristic to try many random directions was used. In fact, we have tests where this heuristic used with correct maximization of f(X, Y) will fail even after 1000000 of random tries. These very special tests were designed by the tester and the probability to guess the direction is lower then 1/11,000,000 in these tests. But look, in this solution for the given direction all prefix sets are tried instead of choosing the optimal set. In fact, after creating those tough test data we develop this hacking solution as well and saw that it easily passes, but we had no time at all to think how to fail it. It was one of the reasons of such delay with this problem. In the second solution some much smarter scheme to crack our test data was used.
After some checks we now realize that by adding the following test
1 1
1 1000000000
1000000000 -1000000000
-1000000000 -1000000000
-1 0
along with the mentioned small test (which presents in test data) we would fail comparison in double.
Finally, let us say a few words about tricky overflows in this problem. The first one can occur when the coordinates of the shifted point A are subtracted from the coordinates of the points B and C, if this calculation is performed in 32-bit signed integer type. This will happen on the following test
500 500
1000000000 1000000000
-1000000000 -1000000000
-1000000000 -1000000000
1000000 1000000 (500 times)
The difference in this case can be 2500000000 = 2.5e9 > 231.
The same test also causes overflow during calculation of the sum of the squares (Bx - X)2 + (By - Y)2, where (X, Y) are the coordinates of the shifted point A, if this calculation performed in 64-bit signed integer type. This sum of squares can equal to 1.25e19 > 263. But using unsigned 64-bit integer type will resolve the issue. We can proudly said that even ACRush was caught by this overflow for some time 
SETTER’S AND TESTER’S SOLUTIONS:
Setter’s solution will be updated soon.
Tester’s solution can be found here.
RELATED PROBLEMS:
CPOINT - Codechef April 2011 Long Contest - Central Point.
The tester has invented this hack with square root exactly during participation at that contest (see his solution here).
