1.1.5 Counting Inversions

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

/*

The number of inversions for a sequence $a$ is the number of ordered pairs $(i, j)$ such that
$i < j$ and $a[i] > a[j]$. This is roughly how "close" an array is to being sorted, but is *not* the
minimum number of swaps required to sort the array. If the array is sorted, then the inversion count
is $0$. If the array is sorted in decreasing order, then the inversion count is maximal. The
following two functions are each techniques to efficiently count inversions. In the merge sort
approach, whenever the merge step emits an element from the right half, that element jumps ahead of
every unmerged left-half element, and exactly that many inversions are added to the count.

- `inversions(lo, hi, comp = std::less<>)` uses merge sort to return the number of inversions given
  two random-access iterators as a range $[`lo`, `hi`)$. The input range will be sorted after the
  function call. Optionally, a comparison function object specifying a strict weak ordering may be
  specified to replace the default `operator<`.
- `inversions(a)` uses coordinate compression and a Fenwick tree to return the number of inversions
  in an integer vector without modifying it.

Time Complexity:
- O(n log n) per call to `inversions(lo, hi)`, where $n$ is the distance between `lo` and `hi`.
- O(n log n) per call to `inversions(a)`.

Space Complexity:
- O(log n) auxiliary stack space and O(n) auxiliary heap space for `inversions(lo, hi)`.
- O(n) auxiliary heap space for `inversions(a)`.

*/

#include <algorithm>
#include <cstdint>
#include <functional>
#include <iterator>
#include <vector>

template<typename It, typename Compare = std::less<>>
int64_t inversions(It lo, It hi, Compare comp = Compare{}) {
  int n = static_cast<int>(hi - lo);
  if (n < 2) {
    return 0;
  }
  It mid = lo + (n - 1) / 2, a = lo, c = mid + 1;
  int64_t res = inversions(lo, mid + 1, comp) + inversions(mid + 1, hi, comp);
  using T = typename std::iterator_traits<It>::value_type;
  std::vector<T> merged;
  merged.reserve(n);
  while (a <= mid && c < hi) {
    if (comp(*c, *a)) {
      merged.push_back(*(c++));
      res += (mid - a) + 1;
    } else {
      merged.push_back(*(a++));
    }
  }
  if (a > mid) {
    for (It it = c; it != hi; ++it) {
      merged.push_back(*it);
    }
  } else {
    for (It it = a; it <= mid; ++it) {
      merged.push_back(*it);
    }
  }
  for (It it = lo; it != hi; ++it) {
    *it = merged[it - lo];
  }
  return res;
}

int64_t inversions(const std::vector<int> &a) {
  int n = static_cast<int>(a.size());
  std::vector<int> values(a);
  std::sort(values.begin(), values.end());
  values.resize(std::unique(values.begin(), values.end()) - values.begin());
  std::vector<int> bit(values.size() + 1);
  int64_t res = 0;
  for (int i = n - 1; i >= 0; i--) {
    int id =
        static_cast<int>(std::lower_bound(values.begin(), values.end(), a[i]) - values.begin()) + 1;
    for (int j = id - 1; j > 0; j -= j & -j) {
      res += bit[j];
    }
    for (int j = id; j < static_cast<int>(bit.size()); j += j & -j) {
      bit[j]++;
    }
  }
  return res;
}

/*** Example Usage ***/

#include <cassert>
using namespace std;

int main() {
  {
    vector<int> a{6, 9, 1, 14, 8, 12, 3, 2};
    assert(inversions(a.begin(), a.end()) == 16);
    assert(is_sorted(a.begin(), a.end()));
  }
  {
    vector<int> a{6, 9, 1, 14, 8, 12, 3, 2};
    assert(inversions(a) == 16);
    assert((a == vector<int>{6, 9, 1, 14, 8, 12, 3, 2}));
  }
  {
    vector<int> a{2, 2, 1};
    assert(inversions(a) == 2);  // Equal elements are not inversions.
  }
  return 0;
}