Skip to content

CDQ Divide and Conquer

This page introduces CDQ divide and conquer.

Introduction

CDQ divide and conquer is a concept rather than a specific algorithm, similar to dynamic programming. Currently, the extensions of this idea are quite broad, and based on different principles and implementations, they can be roughly categorized into three types:

  • Solving problems related to point pairs.
  • Optimization and transition of 1D dynamic programming.
  • Transforming some dynamic problems into static problems through CDQ divide and conquer.

The idea of CDQ divide and conquer was first organized and summarized by Chen Danqi, a gold medalist of IOI2008, during high school, and it is named after her.1

These problems are mostly similar to "Given a sequence of length \(n\), count the number of point pairs \((i,j)\) with certain characteristics" or "Given a sequence of length \(n\), find a point pair \((i,j)\) that maximizes some function."

The algorithm flow for solving such problems with CDQ divide and conquer is as follows:

  1. Find the midpoint \(mid\) of the sequence;

  2. Divide all point pairs \((i,j)\) into three categories:

    1. Point pairs where \(1 \leq i \leq mid, 1 \leq j \leq mid\);
    2. Point pairs where \(1 \leq i \leq mid, mid+1 \leq j \leq n\);
    3. Point pairs where \(mid+1 \leq i \leq n, mid+1 \leq j \leq n\).
  3. Split the sequence \((1,n)\) into two sequences \((1,mid)\) and \((mid+1,n)\). At this point, the first and third categories of point pairs are within these two sequences;

  4. Recursively process these two categories of point pairs;

  5. Handle the second category of point pairs.

It can be seen that the idea of CDQ divide and conquer is to continuously assign point pairs to the left and right intervals through recursion.

In practical applications, we usually use a function solve(l,r) to process point pairs where \(l \leq i \leq r, l \leq j \leq r\). The recursive part of the above algorithm flow is implemented through solve(l,mid) and solve(mid,r). The remaining second category of point pairs requires additional algorithmic design to solve.

Example Problems

3D Partial Order

Given a sequence, each point has three attributes \(a_i, b_i, c_i\). Find: how many point pairs \((i,j)\) in this sequence satisfy \(a_j \leq a_i\), \(b_j \leq b_i\), \(c_j \leq c_i\), and \(j \ne i\).

Solution Idea

3D partial order is a classic problem for CDQ divide and conquer.

The problem requires counting the number of point pairs in the sequence, so try using CDQ divide and conquer.

First, sort the sequence by \(a\).

Suppose we have implemented solve(l,r) and recursively handled solve(l,mid) and solve(mid+1,r). Now we need to count, among the point pairs \((i,j)\) satisfying \(l \leq i \leq mid\), \(mid+1 \leq j \leq r\), how many also satisfy \(a_{i} \leq a_{j}\), \(b_{i} \leq b_{j}\), \(c_{i} \leq c_{j}\).

A little thought reveals that the condition \(a_{i} \leq a_{j}\) is redundant: since \(i\) is less than \(mid\) and \(j\) is greater than \(mid\), \(i\) must be smaller than \(j\); and since the sequence is sorted by \(a\), we must have \(a_{i} \leq a_{j}\). Now we are left with two conditions: \(b_{i} \leq b_{j}\) and \(c_{i} \leq c_{j}\). Based on these conditions, we can enumerate \(j\) and find how many \(i\) satisfy the conditions.

For convenience, we sort all points in \((l,mid)\) and \((mid+1,r)\) by their \(b\) values in ascending order. Then we enumerate each \(j\) and insert all points \(i\) with \(b_{i} \leq b_{j}\) into a data structure (here we choose a Fenwick Tree). At this point, we only need to query how many points in the Fenwick Tree have a \(c\) value less than or equal to \(c_{j}\) to determine how many \(i\) can legally match this \(j\).

When inserting a point with \(c\) value equal to \(x\), we increment the position \(x\) in the Fenwick Tree by one. Querying how many points are less than \(x\) in the Fenwick Tree is essentially finding a prefix sum. As long as we have discretized all \(c\) values in advance, our complexity is correct.

For each \(j\), we need to insert all points \(i\) with \(b_{i} \leq b_{j}\) into the Fenwick Tree. Since all \(i\) and \(j\) are already sorted by \(b\) values, we can insert points into the Fenwick Tree using a two-pointer approach, reducing the number of insertion operations from \(O(n^2)\) to \(O(n)\).

Through this algorithm flow, we process information about the second category of point pairs in \(O(n\log n)\) time. The overall time complexity of the algorithm is \(T(n)=T(\lfloor \frac{n}{2} \rfloor)+T(\lceil \frac{n}{2} \rceil)+O(n\log n)=O(n\log^2n)\).

Example Code
  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
#include <algorithm>
#include <iostream>

constexpr int MAXN = 1e5 + 10;
constexpr int MAXK = 2e5 + 10;

int n, k;

struct Element {
  int a, b, c;
  int cnt;
  int res;

  bool operator!=(Element other) const {
    if (a != other.a) return true;
    if (b != other.b) return true;
    if (c != other.c) return true;
    return false;
  }
};

Element e[MAXN];
Element ue[MAXN];
int m, t;
int res[MAXN];

struct BinaryIndexedTree {
  int node[MAXK];

  int lowbit(int x) { return x & -x; }

  void Add(int pos, int val) {
    while (pos <= k) {
      node[pos] += val;
      pos += lowbit(pos);
    }
    return;
  }

  int Ask(int pos) {
    int res = 0;
    while (pos) {
      res += node[pos];
      pos -= lowbit(pos);
    }
    return res;
  }
} BIT;

bool cmpA(Element x, Element y) {
  if (x.a != y.a) return x.a < y.a;
  if (x.b != y.b) return x.b < y.b;
  return x.c < y.c;
}

bool cmpB(Element x, Element y) {
  if (x.b != y.b) return x.b < y.b;
  return x.c < y.c;
}

void CDQ(int l, int r) {
  if (l == r) return;
  int mid = (l + r) / 2;
  CDQ(l, mid);
  CDQ(mid + 1, r);
  std::sort(ue + l, ue + mid + 1, cmpB);
  std::sort(ue + mid + 1, ue + r + 1, cmpB);
  int i = l;
  int j = mid + 1;
  while (j <= r) {
    while (i <= mid && ue[i].b <= ue[j].b) {
      BIT.Add(ue[i].c, ue[i].cnt);
      i++;
    }
    ue[j].res += BIT.Ask(ue[j].c);
    j++;
  }
  for (int k = l; k < i; k++) BIT.Add(ue[k].c, -ue[k].cnt);
  return;
}

using std::cin;
using std::cout;

int main() {
  cin.tie(nullptr)->sync_with_stdio(false);
  cin >> n >> k;
  for (int i = 1; i <= n; i++) cin >> e[i].a >> e[i].b >> e[i].c;
  std::sort(e + 1, e + n + 1, cmpA);
  for (int i = 1; i <= n; i++) {
    t++;
    if (e[i] != e[i + 1]) {
      m++;
      ue[m].a = e[i].a;
      ue[m].b = e[i].b;
      ue[m].c = e[i].c;
      ue[m].cnt = t;
      t = 0;
    }
  }
  CDQ(1, m);
  for (int i = 1; i <= m; i++) res[ue[i].res + ue[i].cnt - 1] += ue[i].cnt;
  for (int i = 0; i < n; i++) cout << res[i] << '\n';
  return 0;
}
CQOI2011 Dynamic Inversion Pairs

For a sequence \(a\), its number of inversion pairs is defined as the number of elements in the set \(\{(i,j)| i < j \wedge a_i > a_j \}\).

Now given a permutation of \(1\sim n\), delete \(m\) elements in some order. Your task is to count the total number of inversion pairs in the sequence before each deletion.

Example Code
  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
// After careful deduction, the formula is similar to three-dimensional partial order, essentially a template for three-dimensional partial order.
#include <algorithm>
#include <iostream>
using namespace std;
using ll = long long;
int n;
int m;

struct treearray {
  int ta[200010];

  void ub(int& x) { x += x & (-x); }

  void db(int& x) { x -= x & (-x); }

  void c(int x, int t) {
    for (; x <= n + 1; ub(x)) ta[x] += t;
  }

  int sum(int x) {
    int r = 0;
    for (; x > 0; db(x)) r += ta[x];
    return r;
  }
} ta;

struct data_ {
  int val;
  int del;
  int ans;
} a[100010];

int rv[100010];
ll res;

// Override two comparison functions
bool cmp1(const data_& a, const data_& b) { return a.val < b.val; }

bool cmp2(const data_& a, const data_& b) { return a.del < b.del; }

void solve(int l, int r) {  // Below are the specific formulas, applied directly
  if (r - l == 1) {
    return;
  }
  int mid = (l + r) / 2;
  solve(l, mid);
  solve(mid, r);
  int i = l + 1;
  int j = mid + 1;
  while (i <= mid) {
    while (a[i].val > a[j].val && j <= r) {
      ta.c(a[j].del, 1);
      j++;
    }
    a[i].ans += ta.sum(m + 1) - ta.sum(a[i].del);
    i++;
  }
  i = l + 1;
  j = mid + 1;
  while (i <= mid) {
    while (a[i].val > a[j].val && j <= r) {
      ta.c(a[j].del, -1);
      j++;
    }
    i++;
  }
  i = mid;
  j = r;
  while (j > mid) {
    while (a[j].val < a[i].val && i > l) {
      ta.c(a[i].del, 1);
      i--;
    }
    a[j].ans += ta.sum(m + 1) - ta.sum(a[j].del);
    j--;
  }
  i = mid;
  j = r;
  while (j > mid) {
    while (a[j].val < a[i].val && i > l) {
      ta.c(a[i].del, -1);
      i--;
    }
    j--;
  }
  sort(a + l + 1, a + r + 1, cmp1);
  return;
}

int main() {
  cin.tie(nullptr)->sync_with_stdio(false);
  cin >> n >> m;
  for (int i = 1; i <= n; i++) {
    cin >> a[i].val;
    rv[a[i].val] = i;
  }
  for (int i = 1; i <= m; i++) {
    int p;
    cin >> p;
    a[rv[p]].del = i;
  }
  for (int i = 1; i <= n; i++) {
    if (a[i].del == 0) a[i].del = m + 1;
  }
  for (int i = 1; i <= n; i++) {
    res += ta.sum(n + 1) - ta.sum(a[i].val);
    ta.c(a[i].val, 1);
  }
  for (int i = 1; i <= n; i++) {
    ta.c(a[i].val, -1);
  }
  solve(0, n);
  sort(a + 1, a + n + 1, cmp2);
  for (int i = 1; i <= m; i++) {
    cout << res << '\n';
    res -= a[i].ans;
  }
  return 0;
}

CDQ Divide and Conquer for Optimizing 1D/1D Dynamic Programming Transitions

Related content: CDQ Divide and Conquer for Optimizing DP

1D/1D dynamic programming refers to a specific class of DP problems where the DP array is one-dimensional and transitions are \(O(n)\). Under favorable conditions, CDQ divide and conquer can sometimes reduce the time complexity from \(O(n^2)\) to \(O(n\log^2n)\).

For example, given a sequence where each element has two attributes \(a\) and \(b\), we want to compute a DP value with the following transition equation:

\(dp_{i}=1+ \max_{j=1}^{i-1}dp_{j}[a_{j} < a_{i}][b_{j} < b_{i}]\)

This is a DP equation for a 2D longest increasing subsequence, meaning only points \(j\) with \(j < i\), \(a_{j} < a_{i}\), and \(b_{j} < b_{i}\) can update the DP value of point \(i\).

Direct transition is clearly \(O(n^2)\). The following explains how to use CDQ divide and conquer to optimize the transition process.

We observe that the transition from \(dp_{j}\) to \(dp_{i}\) is also a type of point-pair relationship, so we handle it similarly to how CDQ divide and conquer processes point-pair relationships.

This transition process is relatively routine. Suppose we are currently processing the interval \((l,r)\). The algorithm flow is roughly as follows:

  1. If \(l=r\), it means \(dp_{r}\) has already been computed. Simply increment \(dp_{r}\) by one and return;
  2. Recursively call solve(l,mid);
  3. Process all transition relationships where \(l \leq j \leq mid\) and \(mid+1 \leq i \leq r\);
  4. Recursively call solve(mid+1,r).

The third step is similar to CDQ divide and conquer for 3D partial order. When processing transition relationships where \(l \leq j \leq mid\) and \(mid+1 \leq i \leq r\), we no longer need to consider the \(j < i\) constraint. Therefore, we still sort all points \(i\) and \(j\) by their \(a\) values, insert point \(j\) into a Fenwick Tree using a two-pointer approach, and finally query the prefix maximum to update \(dp_{i}\).

Proof of Correctness for the Transition Process

The key difference between this CDQ implementation and the one for handling point-pair relationships lies in processing the part where \(l \leq j \leq mid\) and \(mid+1 \leq i \leq r\). In the point-pair CDQ implementation, this part can be placed anywhere. However, when using CDQ divide and conquer to optimize DP, this step must be placed between solve(l,mid) and solve(mid+1,r). The reason is that DP transitions are ordered and must satisfy two conditions; otherwise, they are incorrect:

  1. All \(dp_{j}\) values used to compute \(dp_{i}\) must be fully computed; no "half-finished" values are allowed.
  2. All \(dp_{j}\) values used to compute \(dp_{i}\) must be able to update \(dp_{i}\); no \(dp_{j}\) values should be left out.

These two conditions are easily satisfied in an \(O(n^2)\) brute-force approach. However, after using CDQ divide and conquer, the transition order is clearly disrupted, so it is necessary to verify the correctness of the transitions.

The recursion tree of CDQ divide and conquer is shown below.

Recursion tree of CDQ divide and conquer

If we execute the algorithm flow described above, take point \(8\) as an example: its DP value is updated in the functions solve(1,8), solve(5,8), and solve(7,8). The three intervals used to update it are \((1,4)\), \((5,6)\), and \((7,7)\), which are disjoint. Similarly, for point \(5\), its DP value is resolved in solve(1,4), and the updating interval is \((1,4)\). Careful observation reveals that the DP value of a point \(i\) is updated \(\log\) times, and the intervals updating it are exactly the \(\log\) intervals into which \((1,i)\) is split on a segment tree. Thus, we indeed ensure that all valid \(j\) have updated point \(i\), satisfying condition 2.

Next, analyze the execution flow of our algorithm:

  1. The first function to finish is solve(1,1). At this point, we find that \(dp_{1}\) has been computed.
  2. The first function to execute the transition process is solve(1,2). At this point, we find that \(dp_{2}\) has been updated.
  3. The second function to finish is solve(2,2). At this point, we find that \(dp_{2}\) has been computed.
  4. Next, solve(1,2) finishes, and the DP values for the interval \((1,2)\) are all computed.
  5. The next function to execute the transition process is solve(1,4). After this transition, we find that \(dp_{3}\) has been updated.
  6. The next function to finish is solve(3,3). We find that \(dp_{3}\) has been computed.
  7. The next transition executed is solve(2,4). At this point, \(dp_{4}\) was updated once by \((1,2)\) in solve(1,4) and again by \((3,3)\), so \(dp_{4}\) is also updated.
  8. solve(4,4) finishes, and \(dp_{4}\) is computed.
  9. solve(3,4) finishes, and the values for \((3,4)\) are computed.
  10. solve(1,4) finishes, and the values for \((1,4)\) are computed.
  11. ……

By simulating the function flow, we observe that whenever solve(l,r) finishes, the DP values for the interval \((l,r)\) are fully computed. Since each time we execute the transition function, solve(l,mid) has already finished, every transition process we execute is valid, satisfying condition 1.

In the process above, we see that if we view the recursion tree of CDQ divide and conquer as a segment tree, then CDQ divide and conquer is the in-order traversal function of this segment tree. Thus, we process all DP values in order, only with the transition order split up, so the algorithm is correct.

Example Problems

SDOI2011 Intercept Missiles

A country has developed a missile interception system to defend against enemy missile attacks. However, this system has a flaw: although its first shot can reach any height and intercept missiles at any speed, each subsequent shot cannot exceed the height of the previous shot, and the speed of the intercepted missile cannot exceed that of the previous shot. One day, radar detects incoming enemy missiles. Since the system is still in the trial phase, only one set is available, so it may not intercept all missiles.

When not all missiles can be intercepted, we naturally choose the plan that minimizes national loss, i.e., intercepts the maximum number of missiles. However, there may be multiple optimal plans for intercepting the maximum number of missiles. If there are multiple optimal plans, we randomly select one as the final interception blueprint.

Our spies have obtained the heights and speeds of all enemy missiles. Your task is to calculate the probability that each missile is intercepted when executing the above decision.

Example Code
  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
// A problem of 2D longest increasing subsequence
// To determine whether an element is in the longest increasing subsequence, run CDQ twice (forward and backward)
#include <algorithm>
#include <iomanip>
#include <iostream>
using namespace std;
using db = double;
constexpr int N = 1e6 + 10;

struct data_ {
  int h;
  int v;
  int p;
  int ma;
  db ca;
} a[2][N];

int n;
bool tr;

// Below are rewritten comparison functions
bool cmp1(const data_& a, const data_& b) {
  if (tr)
    return a.h > b.h;
  else
    return a.h < b.h;
}

bool cmp2(const data_& a, const data_& b) {
  if (tr)
    return a.v > b.v;
  else
    return a.v < b.v;
}

bool cmp3(const data_& a, const data_& b) {
  if (tr)
    return a.p < b.p;
  else
    return a.p > b.p;
}

bool cmp4(const data_& a, const data_& b) { return a.v == b.v; }

struct treearray {
  int ma[2 * N];
  db ca[2 * N];

  void c(int x, int t, db c) {
    for (; x <= n; x += x & (-x)) {
      if (ma[x] == t) {
        ca[x] += c;
      } else if (ma[x] < t) {
        ca[x] = c;
        ma[x] = t;
      }
    }
  }

  void d(int x) {
    for (; x <= n; x += x & (-x)) {
      ma[x] = 0;
      ca[x] = 0;
    }
  }

  void q(int x, int& m, db& c) {
    for (; x > 0; x -= x & (-x)) {
      if (ma[x] == m) {
        c += ca[x];
      } else if (m < ma[x]) {
        c = ca[x];
        m = ma[x];
      }
    }
  }
} ta;

int rk[2][N];

void solve(int l, int r, int t) {  // Recursive processing
  if (r - l == 1) {
    return;
  }
  int mid = (l + r) / 2;
  solve(l, mid, t);
  sort(a[t] + mid + 1, a[t] + r + 1, cmp1);
  int p = l + 1;
  for (int i = mid + 1; i <= r; i++) {
    for (; (cmp1(a[t][p], a[t][i]) || a[t][p].h == a[t][i].h) && p <= mid;
         p++) {
      ta.c(a[t][p].v, a[t][p].ma, a[t][p].ca);
    }
    db c = 0;
    int m = 0;
    ta.q(a[t][i].v, m, c);
    if (a[t][i].ma < m + 1) {
      a[t][i].ma = m + 1;
      a[t][i].ca = c;
    } else if (a[t][i].ma == m + 1) {
      a[t][i].ca += c;
    }
  }
  for (int i = l + 1; i <= mid; i++) {
    ta.d(a[t][i].v);
  }
  sort(a[t] + mid, a[t] + r + 1, cmp3);
  solve(mid, r, t);
  sort(a[t] + l + 1, a[t] + r + 1, cmp1);
}

void ih(int t) {
  sort(a[t] + 1, a[t] + n + 1, cmp2);
  rk[t][1] = 1;
  for (int i = 2; i <= n; i++) {
    rk[t][i] = (cmp4(a[t][i], a[t][i - 1])) ? rk[t][i - 1] : i;
  }
  for (int i = 1; i <= n; i++) {
    a[t][i].v = rk[t][i];
  }
  sort(a[t] + 1, a[t] + n + 1, cmp3);
  for (int i = 1; i <= n; i++) {
    a[t][i].ma = 1;
    a[t][i].ca = 1;
  }
}

int len;
db ans;

int main() {
  cin.tie(nullptr)->sync_with_stdio(false);
  cin >> n;
  for (int i = 1; i <= n; i++) {
    cin >> a[0][i].h >> a[0][i].v;
    a[0][i].p = i;
    a[1][i].h = a[0][i].h;
    a[1][i].v = a[0][i].v;
    a[1][i].p = i;
  }
  ih(0);
  solve(0, n, 0);
  tr = true;
  ih(1);
  solve(0, n, 1);
  tr = true;
  sort(a[0] + 1, a[0] + n + 1, cmp3);
  sort(a[1] + 1, a[1] + n + 1, cmp3);
  for (int i = 1; i <= n; i++) {
    len = max(len, a[0][i].ma);
  }
  cout << len << '\n';
  for (int i = 1; i <= n; i++) {
    if (a[0][i].ma == len) {
      ans += a[0][i].ca;
    }
  }
  cout << fixed << setprecision(5);
  for (int i = 1; i <= n; i++) {
    if (a[0][i].ma + a[1][i].ma - 1 == len) {
      cout << (a[0][i].ca * a[1][i].ca) / ans << ' ';
    } else {
      cout << "0.00000 ";
    }
  }
  return 0;
}

Transforming Dynamic Problems into Static Problems

In the first two cases, the purpose of using CDQ divide and conquer is to recursively process point-pair relationships after halving the sequence to achieve good complexity. However, in this section, what is halved is not a general sequence but a time sequence.

This approach is suitable for some data structure problems that "require supporting xxx modifications and then xxx queries." Such problems have two characteristics:

  • If the queries are offline, all operations naturally form a sequence in time order.
  • Each modification is closely related to subsequent query operations. There will be \(O(n^2)\) such "modification-query" pairs.

We can use CDQ divide and conquer to partition this operation sequence and handle the relationships between modifications and queries.

Similar to CDQ divide and conquer for handling point-pair relationships, assume we are partitioning the sequence \((l,r)\). We first recursively handle the modification-query relationships between \((l,mid)\) and \((mid,r)\), then process all modification-query relationships where \(l \leq i \leq mid\) and \(mid+1 \leq j \leq r\), where \(i\) is a modification and \(j\) is a query.

Note that if the modifications are independent of each other, we do not need to handle the temporal order between \(l \leq i \leq mid\) and \(mid+1 \leq j \leq r\), or between solve(l,mid) and solve(mid+1,r) (e.g., ordinary addition/subtraction problems). However, if the modifications are not independent (e.g., assignment operations), the state of the sequence after a modification may depend on the previous sequence. In this case, the step of processing all modification-query relationships that cross \(mid\) must be placed between solve(l,mid) and solve(mid+1,r). The reason is the same as for CDQ divide and conquer optimizing 1D/1D dynamic programming: partitioning in in-order traversal ensures that each modification is executed strictly in time order.

Example Problems

Rectangle Addition and Rectangle Sum Query

Maintain a 2D array, supporting adding a number to a rectangular region and querying the sum of a rectangular region.

Solution Idea

For the non-modification version of this problem, i.e., "Given a 2D array, query the sum of a rectangular region multiple times," there is a classic approach using scanline with a segment tree. Specifically, each rectangle is split into insertion and deletion operations, each query is split into the difference of 2D prefix sums, and then processed offline. However, the original problem includes modifications, so this approach cannot be directly applied.

Try using CDQ divide and conquer. We offline all query and modification operations. These operations form a sequence with \(O(N^2)\) modification-query pairs. Following the general CDQ divide and conquer flow, we divide all relationships into three categories. In this layer of partitioning, we only handle modification-query relationships that cross \(mid\), and the remaining modification-query relationships are solved recursively.

We find that all modifications are completed before the queries. At this point, the original problem is equivalent to "There are static rectangles on the plane, and we continuously query the sum of a rectangular region."

Use a scanline to handle all modification-query relationships crossing \(mid\) in \(O(n\log n)\) time. The remaining task is to recursively partition the modification-query relationships on the left and right sides.

In this CDQ divide and conquer implementation, each query is processed \(O(\log n)\) times. However, this is fine because the modifications contributing to this query are disjoint each time. The overall time complexity is \(T(n)=T(\lfloor \frac{n}{2} \rfloor)+T(\lceil \frac{n}{2} \rceil)+ O(n\log n)=O(n\log^2n)\).

Observing the algorithm flow, we see that initially we could only solve the static rectangle addition and sum query problem. However, by simply using CDQ divide and conquer, we can offline solve a dynamic rectangle addition and sum query problem. The essence of transforming dynamic problems into static problems lies in CDQ divide and conquer processing only the modification-query relationships crossing a certain point each time. This way, we only need to consider the simple problem of "all queries occur after modifications." It is precisely because of this that CDQ divide and conquer is called a "tool for transforming dynamic problems into static problems."

[Ynoi2016] Insects in the Mirror

Maintain a sequence \(a_i\) of length \(n\), with \(m\) operations.

  1. Change the value of the interval \([l,r]\) to \(x\);
  2. Query how many distinct numbers appear in the interval \([l,r]\), i.e., multiple occurrences of the same number count as one.

In one sentence: interval assignment and interval color counting.

Solution Idea

Maintain the position of the first same-colored point to the left of each position, denoted as \(pre_{i}\). Then interval color counting is transformed into a classic 2D point counting problem.

By treating consecutive segments of the same color as a single point, it can be proven that the total change in \(pre\) is \(O(n+m)\), i.e., each operation only causes \(O(1)\) changes in \(pre\) values. Thus, we can use CDQ divide and conquer to solve the dynamic single-point addition and rectangle sum query problem.

The specific changes in the \(pre\) array can be handled using std::set. This technique of using set to maintain consecutive intervals is also known as old driver tree.

Example Code

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
```cpp
#include <algorithm>
#include <iostream>
#include <map>
#include <set>
#define SNI set<nod>::iterator
#define SDI set<data>::iterator
using namespace std;
constexpr int N = 1e5 + 10;
int n;
int m;
int pre[N];
int npre[N];
int a[N];
int tp[N];
int lf[N];
int rt[N];
int co[N];

struct modi {
  int t;
  int pos;
  int pre;
  int va;

  friend bool operator<(modi a, modi b) { return a.pre < b.pre; }
} md[10 * N];

int tp1;

struct qry {
  int t;
  int l;
  int r;
  int ans;

  friend bool operator<(qry a, qry b) { return a.l < b.l; }
} qr[N];

int tp2;
int cnt;

bool cmp(const qry& a, const qry& b) { return a.t < b.t; }

void modify(int pos, int co)  // Modification function
{
  if (npre[pos] == co) return;
  md[++tp1] = modi{++cnt, pos, npre[pos], -1};
  md[++tp1] = modi{++cnt, pos, npre[pos] = co, 1};
}

namespace prew {
int lst[2 * N];
map<int, int> mp;  // Pre-discretization

void prew() {
  cin.tie(nullptr)->sync_with_stdio(false);
  cin >> n >> m;
  for (int i = 1; i <= n; i++) cin >> a[i], mp[a[i]] = 1;
  for (int i = 1; i <= m; i++) {
    cin >> tp[i] >> lf[i] >> rt[i];
    if (tp[i] == 1) cin >> co[i], mp[co[i]] = 1;
  }
  map<int, int>::iterator it, it1;
  for (it = mp.begin(), it1 = it, ++it1; it1 != mp.end(); ++it, ++it1)
    it1->second += it->second;
  for (int i = 1; i <= n; i++) a[i] = mp[a[i]];
  for (int i = 1; i <= n; i++)
    if (tp[i] == 1) co[i] = mp[co[i]];
  for (int i = 1; i <= n; i++) pre[i] = lst[a[i]], lst[a[i]] = i;
  for (int i = 1; i <= n; i++) npre[i] = pre[i];
}
}  // namespace prew

namespace colist {
struct data {
  int l;
  int r;
  int x;

  friend bool operator<(data a, data b) { return a.r < b.r; }
};

set<data> s;

struct nod {
  int l;
  int r;

  friend bool operator<(nod a, nod b) { return a.r < b.r; }
};

set<nod> c[2 * N];
set<int> bd;

void split(int mid) {  // Split a node into two nodes
  SDI it = s.lower_bound(data{0, mid, 0});
  data p = *it;
  if (mid == p.r) return;
  s.erase(p);
  s.insert(data{p.l, mid, p.x});
  s.insert(data{mid + 1, p.r, p.x});
  c[p.x].erase(nod{p.l, p.r});
  c[p.x].insert(nod{p.l, mid});
  c[p.x].insert(nod{mid + 1, p.r});
}

void del(set<data>::iterator it) {  // Delete an iterator
  bd.insert(it->l);
  SNI it1, it2;
  it1 = it2 = c[it->x].find(nod{it->l, it->r});
  ++it2;
  if (it2 != c[it->x].end()) bd.insert(it2->l);
  c[it->x].erase(it1);
  s.erase(it);
}

void ins(data p) {  // Insert a node
  s.insert(p);
  SNI it = c[p.x].insert(nod{p.l, p.r}).first;
  ++it;
  if (it != c[p.x].end()) {
    bd.insert(it->l);
  }
}

void stv(int l, int r, int x) {  // Interval assignment
  if (l != 1) split(l - 1);
  split(r);
  int p = l;  // After splitting twice, delete all intervals
  while (p != r + 1) {
    SDI it = s.lower_bound(data{0, p, 0});
    p = it->r + 1;
    del(it);
  }
  ins(data{l, r, x});  // Scan the set to process all changed pre values
  for (set<int>::iterator it = bd.begin(); it != bd.end(); ++it) {
    SDI it1 = s.lower_bound(data{0, *it, 0});
    if (*it != it1->l)
      modify(*it, *it - 1);
    else {
      SNI it2 = c[it1->x].lower_bound(nod{0, *it});
      if (it2 != c[it1->x].begin())
        --it2, modify(*it, it2->r);
      else
        modify(*it, 0);
    }
  }
  bd.clear();
}

void ih() {
  int nc = a[1];
  int ccnt = 1;  // Insert a continuous segment into the set
  for (int i = 2; i <= n; i++)
    if (nc != a[i]) {
      s.insert(data{i - ccnt, i - 1, nc}), c[nc].insert(nod{i - ccnt, i - 1});
      nc = a[i];
      ccnt = 1;
    } else {
      ccnt++;
    }
  s.insert(data{n - ccnt + 1, n, a[n]}), c[a[n]].insert(nod{n - ccnt + 1, n});
}
}  // namespace colist

namespace CDQ {
struct treearray  // Fenwick Tree (Binary Indexed Tree)
{
  int ta[N];

  void c(int x, int t) {
    for (; x <= n; x += x & (-x)) ta[x] += t;
  }

  void d(int x) {
    for (; x <= n; x += x & (-x)) ta[x] = 0;
  }

  int q(int x) {
    int r = 0;
    for (; x; x -= x & (-x)) r += ta[x];
    return r;
  }

  void clear() {
    for (int i = 1; i <= n; i++) ta[i] = 0;
  }
} ta;

int srt[N];

bool cmp1(const int& a, const int& b) { return pre[a] < pre[b]; }

void solve(int l1, int r1, int l2, int r2, int L, int R) {  // CDQ
  if (l1 == r1 || l2 == r2) return;
  int mid = (L + R) / 2;
  int mid1 = l1;
  while (mid1 != r1 && md[mid1 + 1].t <= mid) mid1++;
  int mid2 = l2;
  while (mid2 != r2 && qr[mid2 + 1].t <= mid) mid2++;
  solve(l1, mid1, l2, mid2, L, mid);
  solve(mid1, r1, mid2, r2, mid, R);
  if (l1 != mid1 && mid2 != r2) {
    sort(md + l1 + 1, md + mid1 + 1);
    sort(qr + mid2 + 1, qr + r2 + 1);
    for (int i = mid2 + 1, j = l1 + 1; i <= r2; i++) {  // Consider contribution from left to right
      while (j <= mid1 && md[j].pre < qr[i].l) ta.c(md[j].pos, md[j].va), j++;
      qr[i].ans += ta.q(qr[i].r) - ta.q(qr[i].l - 1);
    }
    for (int i = l1 + 1; i <= mid1; i++) ta.d(md[i].pos);
  }
}

void mainsolve() {
  colist::ih();
  for (int i = 1; i <= m; i++)
    if (tp[i] == 1)
      colist::stv(lf[i], rt[i], co[i]);
    else
      qr[++tp2] = qry{++cnt, lf[i], rt[i], 0};
  sort(qr + 1, qr + tp2 + 1);
  for (int i = 1; i <= n; i++) srt[i] = i;
  sort(srt + 1, srt + n + 1, cmp1);
  for (int i = 1, j = 1; i <= tp2; i++) {  // Initialize the value for each query
    while (j <= n && pre[srt[j]] < qr[i].l) ta.c(srt[j], 1), j++;
    qr[i].ans += ta.q(qr[i].r) - ta.q(qr[i].l - 1);
  }
  ta.clear();
  sort(qr + 1, qr + tp2 + 1, cmp);
  solve(0, tp1, 0, tp2, 0, cnt);
  sort(qr + 1, qr + tp2 + 1, cmp);
  for (int i = 1; i <= tp2; i++) cout << qr[i].ans << '\n';
}
}  // namespace CDQ

int main() {
  prew::prew();
  CDQ::mainsolve();
  return 0;
}
```

[HNOI2010] City Construction

PS Country is a large country with many cities. King Louis has racked his brains for urban transportation construction. Louis can build roads between certain cities, with different costs for building roads between different cities.

Louis hopes to build the minimum number of roads to connect all cities. However, due to certain factors, the cost of building roads between cities changes over time. Louis continuously receives messages about changes in road construction costs. He wants to know the minimum total cost to connect the cities immediately after each message. Louis decides to ask for your help to complete this task.

In one sentence: Given a graph supporting dynamic modifications of edge weights, output the minimum cost of the minimum spanning tree after each edge weight modification.

Solution Idea

In fact, there is an approach using segment tree divide and conquer with LCT to solve this problem, but this implementation has a large constant and may require sophisticated constant optimization to pass. Therefore, consider using CDQ divide and conquer to solve this problem.

Unlike problems typically solved by CDQ divide and conquer, there is no modification-query relationship to partition here because we cannot separately consider "the contribution of modifying an edge to the entire graph's minimum spanning tree." The traditional CDQ divide and conquer approach does not seem to work well.

From the previous example, we see that general CDQ divide and conquer has a special connection with segment trees: during CDQ divide and conquer, we implicitly build a segment tree (because the recursion tree of CDQ divide and conquer is a segment tree). Usually, CDQ considers the relationship between left and right children of the segment tree. For this problem, we need to consider the relationship between parent and child; in other words, when processing the interval $solve(l,r)$, if we can reduce the graph's scale to a variable related to the interval length, we can solve the problem.

So how do we design the algorithm specifically?

Suppose we are constructing the edge set of the minimum spanning tree for the interval \((l,r)\), and we know the edge set of its parent's minimum spanning tree. We assign edge weights of \(+ \infty\) and \(-\infty\) to the edges that change in the interval \((l,r)\), respectively, and run Kruskal once for each to find the edges that appear in the minimum spanning tree.

For an edge:

  • If all modified edge weights in the minimum spanning tree are assigned \(+\infty\), and this edge does not appear in the tree, then it cannot appear in the minimum spanning tree for the queries in \((l,r)\). Thus, we only add the tree edges of the minimum spanning tree to the edge set for \((l,r)\).
  • If all modified edge weights in the minimum spanning tree are assigned \(-\infty\), and this edge appears in the tree, then it must appear in the minimum spanning tree for the interval \((l,r)\). In this case, we can use a union-find set to contract the points connected by these edges and add the edge weights to the answer.

This way, we construct the edge set for the interval \((l,r)\). The minimum spanning tree computed with these edges is equivalent to directly computing the original graph's minimum spanning tree.

So why is our complexity correct?

First, modified edges are always added to our edge set, and the number of such edges is \(O(len)\).

Next, we need to prove that there are not too many unmodified edges in the edge set. We only add tree edges from the minimum spanning tree with all edge weights set to \(+\infty\), so the number of edges we add does not exceed the number of vertices in the current graph.

Now we only need to prove that the number of vertices per recursion level is \(O(len)\) to show that the number of edges is \(O(len)\).

Proving that the number of vertices is \(O(len)\) becomes straightforward. Each time we recursively descend, we contract edges that are unmodified edges appearing in the \(-\infty\) spanning tree. Conversely, we cut all modified edges that appear in the \(-\infty\) spanning tree. Clearly, we cut at most \(len\) edges, so the graph splits into at most \(O(len)\) connected components, meaning the new graph has \(O(len)\) vertices. Thus, we prove that the graph we run Kruskal on each time is \(O(len)\) in size, so the time complexity per level is \(O(n\log n)\).

The overall time complexity is \(T(n)=T(\lfloor \frac{n}{2} \rfloor)+T(\lceil \frac{n}{2} \rceil)+ O(n\log n)=O(n\log^2n)\).

The code implementation may be somewhat challenging. Note that the union-find set cannot use path compression, as it would not support rollback operations. Also, when performing contraction operations, there is no need to actually contract; instead, each level's Kruskal can be directly performed on the previous level's union-find set.

Example Code

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
```cpp
#include <algorithm>
#include <iostream>
#include <stack>
#include <vector>
using namespace std;
using ll = long long;
int n;
int m;
int ask;

struct bcj {
  int fa[20010];
  int size[20010];

  struct opt {
    int u;
    int v;
  };

  stack<opt> st;

  void ih() {
    for (int i = 1; i <= n; i++) fa[i] = i, size[i] = 1;
  }

  int f(int x) { return (fa[x] == x) ? x : f(fa[x]); }

  void u(int x, int y) {  // Union with rollback capability
    int u = f(x);
    int v = f(y);
    if (u == v) return;
    if (size[u] < size[v]) swap(u, v);
    size[u] += size[v];
    fa[v] = u;
    opt o;
    o.u = u;
    o.v = v;
    st.push(o);
  }

  void undo() {
    opt o = st.top();
    st.pop();
    fa[o.v] = o.v;
    size[o.u] -= size[o.v];
  }

  void clear(int tim) {
    while (st.size() > tim) {
      undo();
    }
  }
} s, s1;

struct edge  // Static edge
{
  int u;
  int v;
  ll val;
  int mrk;

  friend bool operator<(edge a, edge b) { return a.val < b.val; }
} e[50010];

struct moved {
  int u;
  int v;
};  // Dynamic edge

struct query {
  int num;
  ll val;
  ll ans;
} q[50010];

bool book[50010];  // Query marker
vector<edge> ve[30];
vector<moved> vq;
vector<edge> tr;
ll res[30];
int tim[30];

void pushdown(int dep)  // Edge contraction
{
  tr.clear();  // Make a copy here to avoid rollback issues
  for (int i = 0; i < ve[dep].size(); i++) {
    tr.push_back(ve[dep][i]);
  }
  sort(tr.begin(), tr.end());
  for (int i = 0; i < tr.size(); i++) {  // Useless edges
    if (s1.f(tr[i].u) == s1.f(tr[i].v)) {
      tr[i].mrk = -1;
      continue;
    }
    s1.u(tr[i].u, tr[i].v);
  }
  s1.clear(0);
  res[dep + 1] = res[dep];
  for (int i = 0; i < vq.size(); i++) {
    s1.u(vq[i].u, vq[i].v);
  }
  vq.clear();
  for (int i = 0; i < tr.size(); i++) {  // Essential edges
    if (tr[i].mrk == -1 || s1.f(tr[i].u) == s1.f(tr[i].v)) continue;
    tr[i].mrk = 1;
    s1.u(tr[i].u, tr[i].v);
    s.u(tr[i].u, tr[i].v);
    res[dep + 1] += tr[i].val;
  }
  s1.clear(0);
  ve[dep + 1].clear();
  for (int i = 0; i < tr.size(); i++) {  // Contract edges
    if (tr[i].mrk != 0) continue;
    edge p;
    p.u = s.f(tr[i].u);
    p.v = s.f(tr[i].v);
    if (p.u == p.v) continue;
    p.val = tr[i].val;
    p.mrk = 0;
    ve[dep + 1].push_back(p);
  }
  return;
}

void solve(int l, int r, int dep) {
  tim[dep] = s.st.size();
  int mid = (l + r) / 2;
  if (r - l == 1) {  // Termination condition
    edge p;
    p.u = s.f(e[q[r].num].u);
    p.v = s.f(e[q[r].num].v);
    p.val = q[r].val;
    e[q[r].num].val = q[r].val;
    p.mrk = 0;
    ve[dep].push_back(p);
    pushdown(dep);
    q[r].ans = res[dep + 1];
    s.clear(tim[dep - 1]);
    return;
  }
  for (int i = l + 1; i <= mid; i++) {
    book[q[i].num] = true;
  }
  for (int i = mid + 1; i <= r; i++) {  // Dynamic to static
    if (book[q[i].num]) continue;
    edge p;
    p.u = s.f(e[q[i].num].u);
    p.v = s.f(e[q[i].num].v);
    p.val = e[q[i].num].val;
    p.mrk = 0;
    ve[dep].push_back(p);
  }
  for (int i = l + 1; i <= mid; i++) {  // Query to dynamic
    moved p;
    p.u = s.f(e[q[i].num].u);
    p.v = s.f(e[q[i].num].v);
    vq.push_back(p);
  }
  pushdown(dep);  // Below is rollback
  for (int i = mid + 1; i <= r; i++) {
    if (book[q[i].num]) continue;
    ve[dep].pop_back();
  }
  for (int i = l + 1; i <= mid; i++) {
    book[q[i].num] = false;
  }
  solve(l, mid, dep + 1);
  for (int i = 0; i < ve[dep].size(); i++) {
    ve[dep][i].mrk = 0;
  }
  for (int i = mid + 1; i <= r; i++) {
    book[q[i].num] = true;
  }
  for (int i = l + 1; i <= mid; i++) {  // Dynamic to static
    if (book[q[i].num]) continue;
    edge p;
    p.u = s.f(e[q[i].num].u);
    p.v = s.f(e[q[i].num].v);
    p.val = e[q[i].num].val;
    p.mrk = 0;
    ve[dep].push_back(p);
  }
  for (int i = mid + 1; i <= r; i++) {  // Query to dynamic
    book[q[i].num] = false;
    moved p;
    p.u = s.f(e[q[i].num].u);
    p.v = s.f(e[q[i].num].v);
    vq.push_back(p);
  }
  pushdown(dep);
  solve(mid, r, dep + 1);
  s.clear(tim[dep - 1]);
  return;  // Time travel back to previous layer
}

int main() {
  cin.tie(nullptr)->sync_with_stdio(false);
  cin >> n >> m >> ask;
  s.ih();
  s1.ih();
  for (int i = 1; i <= m; i++) {
    cin >> e[i].u >> e[i].v >> e[i].val;
  }
  for (int i = 1; i <= ask; i++) {
    cin >> q[i].num >> q[i].val;
  }
  for (int i = 1; i <= ask; i++) {  // Initial dynamic edges
    book[q[i].num] = true;
    moved p;
    p.u = e[q[i].num].u;
    p.v = e[q[i].num].v;
    vq.push_back(p);
  }
  for (int i = 1; i <= m; i++) {  // Initial static edges
    if (book[i]) continue;
    ve[1].push_back(e[i]);
  }
  for (int i = 1; i <= ask; i++) {
    book[q[i].num] = false;
  }
  solve(0, ask, 1);
  for (int i = 1; i <= ask; i++) {
    cout << q[i].ans << '\n';
  }
  return 0;
}
```

References and Notes