PROBLEM LINK:
Practice
Contest: Division 1
Contest: Division 2
Contest: Division 3
Author: Agnimandur
Tester: Danny Mittal
Editorialist: Nishank Suresh
DIFFICULTY:
Medium - Hard
PREREQUISITES:
Mo’s algorithm
PROBLEM:
You have N distinct points P_1, P_2, \dots, P_N on the plane, and Q queries of the form \left(L, R, X\right).
For each query, consider the polygon formed by the points (0, 0), (X, 0), and points P_L, P_{L+1}, \dots, P_R.
If this polygon is not convex, output -1.
If it is convex, triangulate it by treating (X, 0) as a pivot, and print (twice) the maximum area among the formed triangles.
QUICK EXPLANATION
- Use Mo’s algorithm to process the queries
- Maintain a set of currently active points, sorted by angle from the origin. Also maintain the number of bad turns of the current polygon, to easily tell whether it is convex or not.
- Adding a new point translates to breaking an edge of the current polygon and inserting two more, while removing a point translates to breaking two edges and inserting one
- For a given (X, 0), the area of each triangle looks like AX + B for some constants A and B, where A is the difference of two y-coordinates of consecutive points and B is their cross product. There can be at most \sqrt{2C} distinct values of A (because their total sum is 2C), so iterate over all of these for each query and find the best value of B.
EXPLANATION:
As the constraints and time limit might suggest, we will use Mo’s algorithm to solve this problem.
To be able to use Mo’s algorithm, we need to know what information to maintain about the current set of points, and how to update that information when adding/removing a point.
We require two pieces of information - one to determine whether the current polygon is convex, and another to actually compute the answer when we know that it is convex.
Checking convexity
This part is relatively simple - we maintain the number of angles which are strictly larger than 180^\circ.
To do so, consider the current set of points S to be sorted in lexicographical order.
- When adding a new point p, let q_1 < p < p_1 be the points immediately before and after p in S.
Inserting p is equivalent to erasing the edge between p_1 and q_1 and adding edges (p_1, p) and (p, q_1), which in turn removes two old angles and adds 3 new ones - all of which can be processed in \mathcal{O}(1) time. - Similarly, deleting a point is equivalent to deleting two edges and adding a new one, which means three old angles are removed and two new ones are added.
Given that this information is maintained, checking convexity is trivial.
Computing the maximum area
Now we only care about the case where the polygon is convex.
Consider any triangle in the triangulation we obtain, with adjacent vertices (x_1, y_1), (x_2, y_2), (X, 0) where x_1 < x_2.
The (doubled) area of this triangle is known to be
\begin{vmatrix}
X & 0 & 1 \\
x_1 & y_1 & 1 \\
x_2 & y_2 & 1
\end{vmatrix}
= X\cdot (y_1 - y_2) + (x_1y_2 - x_2y_1)
Now, note that because of convexity, the sum of |y_1 - y_2| over all pairs of adjacent vertices in the polygon equals 2M, where M is the maximum y-coordinate of any point.
Thus, there are \mathcal{O}(\sqrt{M}) distinct values of |y_1 - y_2|.
Suppose we fix a value of y_1 - y_2 > 0. Then, X\cdot (y_1 - y_2) + (x_1y_2 - x_2y_1) is maximized when x_1y_2 - x_2y_1 is maximized.
Further, if y_1 - y_2 < 0, we can swap the second and third rows of the determinant to see that the area is also equal to X\cdot (y_2 - y_1) + (x_2y_1 - x_1y_2), which is maximized when x_2y_1 - x_1y_2 is maximized, i.e, x_1y_2 - x_2y_1 is minimized.
So, once we fix a value of |y_1 - y_2|, it is enough to test the maximum and minimum value of x_1y_2 - x_2y_1 among all edges with this slope.
Since there are only \mathcal{O}(\sqrt{C}) distinct values of |y_1 - y_2|, it suffices to iterate over all of them and take the best answer.
This also tells us the information which needs to be maintained - the set of active values of |y_1 - y_2| among adjacent vertices of the polygon, and a multiset of their x_1y_2 - x_2y_1 values.
Updating these when adding/removing an edge is easy.
TIME COMPLEXITY:
\mathcal{O}(N\sqrt{N}\log N + Q(\sqrt{M} + \log N))
CODE:
Setter (C++)
#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <algorithm>
#include <cmath>
#include <vector>
#include <set>
#include <map>
#include <unordered_set>
#include <unordered_map>
#include <queue>
#include <ctime>
#include <cassert>
#include <complex>
#include <string>
#include <cstring>
#include <chrono>
#include <random>
#include <bitset>
#include <array>
using namespace std;
#ifdef LOCAL
#define eprintf(...) {fprintf(stderr, __VA_ARGS__);fflush(stderr);}
#else
#define eprintf(...) 42
#endif
using ll = long long;
using ld = long double;
using uint = unsigned int;
using ull = unsigned long long;
template<typename T>
using pair2 = pair<T, T>;
using pii = pair<int, int>;
using pli = pair<ll, int>;
using pll = pair<ll, ll>;
mt19937_64 rng(chrono::steady_clock::now().time_since_epoch().count());
ll myRand(ll B) {
return (ull)rng() % B;
}
#define pb push_back
#define mp make_pair
#define all(x) (x).begin(),(x).end()
#define fi first
#define se second
clock_t startTime;
double getCurrentTime() {
return (double)(clock() - startTime) / CLOCKS_PER_SEC;
}
struct Point {
int x, y;
Point() : x(), y() {}
Point(int _x, int _y) : x(_x), y(_y) {}
void scan() {
scanf("%d%d", &x, &y);
}
Point operator + (const Point &a) const {
return Point(x + a.x, y + a.y);
}
Point operator - (const Point &a) const {
return Point(x - a.x, y - a.y);
}
ll operator % (const Point &a) const {
return (ll)x * a.x + (ll)y * a.y;
}
ll operator * (const Point &a) const {
return (ll)x * a.y - (ll)y * a.x;
}
bool operator < (const Point &a) const {
ll z = *this * a;
if (z != 0) return z < 0;
return x < a.x;
}
};
const ll INF = (ll)1e18;
const ll D = (ll)2e12;
int B;
const int N = 10010;
const int QQ = 200200;
const int C = 100100;
int n;
Point a[N];
pair<Point, int> ord[N];
int id[N];
set<int> setik;
int bad = 0;
int L, R;
ll ANS[QQ];
set<ll> vals[2 * C];
ll maxVal[2 * C];
set<int> activeY;
inline bool getBad(int x, int y, int z) {
if (x == -1 || y == -1 || z == -1) return false;
return (ord[y].first - ord[x].first) * (ord[z].first - ord[y].first) > 0;
}
inline void addSide(int p, int q) {
if (p == -1 || q == -1) return;
int k = ord[q].first.y - ord[p].first.y;
ll v = (D + ord[q].first * ord[p].first) * N + p;
activeY.insert(k);
vals[C + k].insert(-v);
maxVal[C + k] = -(*vals[C + k].begin());
}
inline void remSide(int p, int q) {
if (p == -1 || q == -1) return;
int k = ord[q].first.y - ord[p].first.y;
ll v = (D + ord[q].first * ord[p].first) * N + p;
vals[C + k].erase(-v);
if (vals[C + k].empty()) {
activeY.erase(k);
} else {
maxVal[C + k] = -(*vals[C + k].begin());
}
}
void myAdd(int p) {
//eprintf("myAdd %d : (%d %d)\n", p, ord[p].first.x, ord[p].first.y);
auto it = setik.upper_bound(p);
int q1 = -1, q2 = -1, p1 = -1, p2 = -1;
if (it != setik.end()) {
q1 = *it;
it++;
if (it != setik.end()) q2 = *it;
it--;
}
assert(it != setik.begin());
it--;
p1 = *it;
assert(p1 < p);
if (it != setik.begin()) p2 = *prev(it);
bad -= getBad(p2, p1, q1);
bad -= getBad(p1, q1, q2);
bad += getBad(p2, p1, p);
bad += getBad(p1, p, q1);
bad += getBad(p, q1, q2);
remSide(p1, q1);
addSide(p1, p);
addSide(p, q1);
setik.insert(p);
}
void myRem(int p) {
//eprintf("myRem %d : (%d %d)\n", p, ord[p].first.x, ord[p].first.y);
auto it = setik.lower_bound(p);
assert(it != setik.begin() && it != setik.end() && *it == p);
int q1 = -1, q2 = -1, p1 = -1, p2 = -1;
it++;
if (it != setik.end()) {
q1 = *it;
it++;
if (it != setik.end()) q2 = *it;
it--;
}
it--;
it--;
p1 = *it;
if (it != setik.begin()) p2 = *prev(it);
bad += getBad(p2, p1, q1);
bad += getBad(p1, q1, q2);
bad -= getBad(p2, p1, p);
bad -= getBad(p1, p, q1);
bad -= getBad(p, q1, q2);
remSide(p1, p);
remSide(p, q1);
addSide(p1, q1);
setik.erase(p);
}
ll solve(int X) {
if (bad > 0) return -1;
Point P, Q, R = Point(X, 0);
assert((int)setik.size() >= 2);
auto it = setik.end();
it--;
Q = ord[*it].first;
it--;
P = ord[*it].first;
if ((Q - P) * (R - Q) > 0) return -1;
ll ans = 0;
for (int k : activeY) {
ll w = maxVal[C + k] / N - D;
ans = max(ans, (ll)k * X + w);
}
return ans;
}
struct Query {
int l, r, x, id;
Query() : l(), r(), x(), id() {}
void scan(int _id) {
id = _id;
scanf("%d%d%d", &l, &r, &x);
r++;
}
bool operator < (const Query &Q) const {
int id1 = l / B, id2 = Q.l / B;
if (id1 != id2) return id1 < id2;
return (id1 & 1) ^ (r < Q.r);
}
};
Query Q[QQ];
int main()
{
startTime = clock();
// freopen("input.txt", "r", stdin);
// freopen("output.txt", "w", stdout);
scanf("%d", &n);
n++;
for (int i = 1; i < n; i++) {
a[i].scan();
ord[i] = mp(a[i], i);
}
sort(ord, ord + n);
for (int i = 0; i < n; i++)
id[ord[i].second] = i;
int q;
scanf("%d", &q);
for (int i = 0; i < q; i++)
Q[i].scan(i);
B = max(1, (int)(n / sqrt(q) + 0.5));
sort(Q, Q + q);
setik.insert(0);
L = R = 1;
for (int i = 0; i < q; i++) {
while(R < Q[i].r) myAdd(id[R++]);
while(L > Q[i].l) myAdd(id[--L]);
while(R > Q[i].r) myRem(id[--R]);
while(L < Q[i].l) myRem(id[L++]);
ANS[Q[i].id] = solve(Q[i].x);
}
for (int i = 0; i < q; i++)
printf("%lld\n", ANS[i]);
return 0;
}
Tester (Kotlin)
import java.io.BufferedInputStream
import java.util.*
import kotlin.math.max
const val BILLION = 1000000000
const val X_LIMIT = 100000
fun main(omkar: Array<String>) {
val jin = FastScanner()
val n = jin.nextInt(10000)
val psOriginal = arrayOf(Point(0L, 0L)) + Array(n) {
val x = jin.nextInt(X_LIMIT, false).toLong()
val y = jin.nextInt(X_LIMIT).toLong()
Point(x, y)
}
val indexes = (0..n).sortedWith(object: Comparator<Int> {
override fun compare(j: Int, k: Int): Int {
val comparison = (if (j == 0) .0 else (psOriginal[j].x.toDouble() / psOriginal[j].y.toDouble())).compareTo(if (k == 0) .0 else (psOriginal[k].x.toDouble() / psOriginal[k].y.toDouble()))
if (comparison != 0) {
return comparison
}
return psOriginal[j].x.compareTo(psOriginal[k].x)
}
})
val ps = Array(n + 1) { psOriginal[indexes[it]] }
val where = IntArray(n + 1)
for (j in 0..n) {
where[indexes[j]] = j
}
val q = jin.nextInt(200000)
val queries = Array(q) {
val l = jin.nextInt(1, n, false)
val r = jin.nextInt(l, n, false)
val x = jin.nextInt(BILLION).toLong()
Query(l, r, x, it)
}
queries.sortWith(compareBy({ it.from / 22 }, { it.to }))
val points = TreeSet<Int>()
val left = IntArray(n + 1)
val right = IntArray(n + 1)
points.add(0)
var amtNotConvex = 0
var from = 1
var to = 0
val lines = mutableMapOf<Long, TreeSet<Pair<Long, Int>>>()
fun accountForNonConvexity(center: Int, delta: Int) {
if (center != 0 && right[center] != 0 && !areConvex(ps[left[center]], ps[center], ps[right[center]])) {
amtNotConvex += delta
}
}
var position = 0
fun addLine(j: Int) {
val (a, b) = lineFrom(ps[j], ps[right[j]])
lines.computeIfAbsent(a) { TreeSet(compareBy({ it.first }, { it.second }))}.add(Pair(b, j))
}
fun removeLine(j: Int) {
val (a, b) = lineFrom(ps[j], ps[right[j]])
val treeSet = lines[a]!!
if (treeSet.size == 1) {
lines.remove(a)
} else {
treeSet.remove(Pair(b, j))
}
}
fun addPoint(j: Int) {
val lower = points.lower(j)!!
val higher = right[lower]
if (higher != 0) {
accountForNonConvexity(lower, -1)
accountForNonConvexity(higher, -1)
removeLine(lower)
}
points.add(j)
left[j] = lower
right[j] = higher
right[lower] = j
left[higher] = j
accountForNonConvexity(lower, 1)
addLine(lower)
if (higher != 0) {
accountForNonConvexity(j, 1)
accountForNonConvexity(higher, 1)
addLine(j)
}
}
fun removePoint(j: Int) {
val lower = left[j]
val higher = right[j]
accountForNonConvexity(lower, -1)
removeLine(lower)
if (higher != 0) {
accountForNonConvexity(j, -1)
accountForNonConvexity(higher, -1)
removeLine(j)
}
points.remove(j)
right[lower] = higher
left[higher] = lower
if (higher != 0) {
accountForNonConvexity(lower, 1)
accountForNonConvexity(higher, 1)
addLine(lower)
}
}
val answer = LongArray(q)
for (query in queries) {
position++
while (to < query.to) {
to++
addPoint(where[to])
}
while (from > query.from) {
from--
addPoint(where[from])
}
while (to > query.to) {
removePoint(where[to])
to--
}
while (from < query.from) {
removePoint(where[from])
from++
}
if (amtNotConvex == 0 && areConvex(ps[points.lower(points.last()!!)!!], ps[points.last()!!], Point(query.x, 0L))) {
for ((a, treeSet) in lines) {
val b = treeSet.last().first
answer[query.index] = max(answer[query.index], (a * query.x) + b)
}
} else {
answer[query.index] = -1L
}
}
println(answer.joinToString("\n"))
jin.assureInputDone()
}
fun sideOfLine(a: Point, b: Point, p: Point): Long {
val left = (b.x - a.x) * (p.y - a.y)
val right = (b.y - a.y) * (p.x - a.x)
return left - right
}
fun areConvex(p1: Point, p2: Point, p3: Point) = sideOfLine(p1, p2, p3) <= 0L
data class Point(val x: Long, val y: Long)
fun lineFrom(left: Point, right: Point) = Point(right.y - left.y, (right.x * left.y) - (left.x * right.y))
data class Query(val from: Int, val to: Int, val x: Long, val index: Int) {
var position: Int = 0
}
class FastScanner {
private val BS = 1 shl 16
private val NC = 0.toChar()
private val buf = ByteArray(BS)
private var bId = 0
private var size = 0
private var c = NC
private var `in`: BufferedInputStream? = null
constructor() {
`in` = BufferedInputStream(System.`in`, BS)
}
private val char: Char
private get() {
while (bId == size) {
size = try {
`in`!!.read(buf)
} catch (e: Exception) {
return NC
}
if (size == -1) return NC
bId = 0
}
return buf[bId++].toChar()
}
fun assureInputDone() {
if (char != NC) {
throw IllegalArgumentException("excessive input")
}
}
fun nextInt(endsLine: Boolean): Int {
var neg = false
c = char
if (c !in '0'..'9' && c != '-' && c != ' ' && c != '\n') {
throw IllegalArgumentException("found character other than digit, negative sign, space, and newline")
}
if (c == '-') {
neg = true
c = char
}
var res = 0
while (c in '0'..'9') {
res = (res shl 3) + (res shl 1) + (c - '0')
c = char
}
if (endsLine) {
if (c != '\n') {
throw IllegalArgumentException("found character other than newline, character code = ${c.toInt()}")
}
} else {
if (c != ' ') {
throw IllegalArgumentException("found character other than space, character code = ${c.toInt()}")
}
}
return if (neg) -res else res
}
fun nextInt(from: Int, to: Int, endsLine: Boolean = true): Int {
val res = nextInt(endsLine)
if (res !in from..to) {
throw IllegalArgumentException("$res not in range $from..$to")
}
return res
}
fun nextInt(to: Int, endsLine: Boolean = true) = nextInt(1, to, endsLine)
}