康託展開+逆展開（Cantor expension）詳解+優化

chesium發表於2021-10-26
優化

康託展開

引入

康託展開（Cantor expansion）用於將排列轉換為字典序的索引（逆展開則相反）
百度百科
維基百科

方法

假設我們要求排列 5 2 4 1 3 的字典序索引

逐位處理：

  • 第一位：5 2 4 1 3，如果一個排列的第一位比 \(5\) 小（有 \(4\) 種情況）
    則不管其後 \(4\) 位如何（有 \(4!\) 種情況），其字典序都更小
    所以，至少有 \(4\times 4!\) 個排列字典序更小。
  • 第二位：5 2 4 1 3，如果另一個排列的第一位就是 \(5\) ，但第二位比 \(2\) 小（有 \(1\) 種情況）
    則不管其後 \(3\) 位如何（有 \(3!\) 種情況），其字典序都更小
    所以， 至少有 \(4\times 4!+1\times 3!\) 個排列字典序更小。
  • 第三位：5 2 4 1 3，如果另一個排列的前兩位與我們的相同，但第三位比 \(4\) 小（\(2\) 不能選了，從右往左看，有 \(2\) 種情況）
    則不管其後 \(2\) 位如何（有 \(2!\) 種情況），其字典序都更小
    所以， 至少有 \(4\times 4!+1\times 3!+2\times 2!\) 個排列字典序更小。
  • 第四位：5 2 4 1 3，如果另一個排列的前三位與我們的相同，但第四位比 \(1\) 小（不可能， 有 \(0\) 種情況）
    則不管其後 \(1\) 位如何（有 \(1!\) 種情況），其字典序都更小
    所以，至少有 \(4\times 4!+1\times 3!+2\times 2!+0\times1!\) 個排列字典序更小。
  • 第五位：5 2 4 1 3，按照上面的方法操作，很顯然， 有 \(4\times 4!+1\times 3!+2\times 2!+0\times1!+0\times0!\) 個排列字典序更小。

因此，若索引從 \(1\) 開始，則 5 2 4 1 3 的索引是 \(4\times 4!+1\times 3!+2\times 2!+0\times1!+0\times0!+1=107\) 。

演算法

總結上述方法，可以歸納出以下演算法：

  • 列舉排列的每一位，對值為 \(p_k\) 的第 \(k\) 位：

    • 找出後面所有位（ \(k+1\) 至 \(n\) ）中小於 \(p_k\) 的位數 \(a_k\) （也就是 \(p_k\) 是第 \(k\) 至 \(n\) 位第 \(a_k+1\) 小的）

    若用這些位中的某一位替換第 \(k\) 位，則無論後面 \(n-k\) 位如何排列（總共有 \(a_k(n-k)!\) 種情況），最終的字典序肯定更小

    • 把這些字典序更小的排列數加起來

  • 再加 \(1\) 即為該排列的字典序索引

公式

用公式來表示即為：

\[\mathrm{Index}(p)=1+\sum^n_{k=1} (n-k)! |\{p_i\mid p_i<p_k,\ k+1\leq i\leq n\}| \]

優化

其中 \(a_k\) 的求值過程可以進行優化，設一序列 \(P=\underbrace{[1,1,1,\dots,1]}_{n\text{個}1}\)
我們每處理一位就置 \(P_{p_k}\) 為 \(0\) 。
這樣在排列 \(p\) 中，我們沒處理過的值就可以被表示為序列 \(P\) 中值為 \(1\) 的索引，即：

\[\{p_i\mid k+1\leq i\leq n\}=\{i\mid P_i=1\} \]

我們可以用線段樹、樹狀陣列這樣的資料結構來維護 \(P\) ，要求小於 \(p_k\) 的位數，即是求序列 \(P\) 區間 \([1,p_k-1]\) 中 \(1\) 的個數，這不就是區間求和嗎？
演算法可以改進為以下這樣：

優化演算法

  • 初始化長度為 \(n\)（索引從 \(1\) 開始），值全為 \(1\) 的序列 \(P\)
  • 列舉排列的每一位，對值為 \(p_k\) 的第 \(k\) 位：
    • 修改 \(P_{p_k}=0\)
    • 求出 \(a=(n-k)!\sum^{p_k-1}_{i=1}P_i\) ，和式為序列 \(P\) 區間 \([1,p_k-1]\) 的元素和
  • 把所有 \(a\) 相加，最後加 \(1\) 即為索引。

時間複雜度為 \(\mathrm{O}(n\log n)\)

參考程式碼

題目：P5367 【模板】康託展開 - 洛谷
（線段樹版本）

#include <cctype>
#include <cstdio>
#include <cstring>
#include <iostream>
//---//
#include <algorithm>
#include <cmath>
#include <vector>
using namespace std;

typedef unsigned int u;
typedef long long ll;
typedef unsigned long long llu;

#define rep(i, a, b) for (ll i = a; i < b; i++)
#define REP(i, a, b) for (ll i = a; i <= b; i++)
#define per(i, b, a) for (ll i = b; i >= a; i--)

const ll N = 1000005;
const ll mod = 998244353;

ll p[N], f[N], d[4 * N];
bool lz[4 * N];

#define mid ll m = s + ((t - s) >> 1)

void bd(ll s, ll t, ll i) {
  if (s == t) {
    d[i] = 1;
    return;
  }
  mid;
  bd(s, m, i * 2);
  bd(m + 1, t, i * 2 + 1);
  d[i] = (d[i * 2] + d[i * 2 + 1]) % mod;
}

void spr(ll s, ll t, ll i) {
  if (lz[i]) {
    d[2 * i] = 0;
    d[2 * i + 1] = 0;
    lz[2 * i] = true;
    lz[2 * i + 1] = true;
    lz[i] = false;
  }
}

void upd(ll l, ll r, ll s, ll t, ll i) {
  if (l <= s && t <= r) {
    d[i] = 0;
    lz[i] = true;
    return;
  }
  mid;
  spr(s, t, i);
  if (l <= m) upd(l, r, s, m, 2 * i);
  if (r > m) upd(l, r, m + 1, t, 2 * i + 1);
  d[i] = (d[i * 2] + d[i * 2 + 1]) % mod;
}

ll get(ll l, ll r, ll s, ll t, ll i) {
  if (l > r) return 0;
  if (l <= s && t <= r) return d[i];
  mid;
  spr(s, t, i);
  ll sum;
  if (l <= m) sum = get(l, r, s, m, 2 * i) % mod;
  if (r > m) sum = (sum + get(l, r, m + 1, t, 2 * i + 1)) % mod;
  return sum;
}

signed main() {
  ll n, r = 1;
  f[0] = 1;
  scanf("%lld", &n);
  bd(1, n, 1);
  REP(i, 1, n) scanf("%lld", &p[i]);
  REP(i, 1, n - 1) f[i] = f[i - 1] * i % mod;
  REP(k, 1, n) {
    upd(p[k], p[k], 1, n, 1);
    r = (r + get(1, p[k] - 1, 1, n, 1) * f[n - k] % mod) % mod;
  }
  printf("%lld", r);
}

// https://www.luogu.com.cn/record/60977896

康託逆展開

引入

康託逆展開用於通過一個已知長度排列的字典序索引反求出該排列

推導

剛剛，我們知道了：

\[\mathrm{Index}(p)=1+\sum^n_{k=1}(n-k)!a_k \]

其中：

\[a_k=|\{p_i\mid p_i<p_k,\ k+1\leq i\leq n\}| \]

已知 \(n\) 位排列 \(p\) 的字典序索引為 \(\mathrm{Index}(p)\) ，通過康託逆展開，我們可以求出各個 \(a_i\) ，從而求出排列 \(p\) 。
首先，我們把右側的 \(1\) 移至左側，再提出和式中的第一項：

\[\mathrm{Index}(p)-1=(n-1)!a_1+\sum^n_{k=2}(n-k)!a_k \quad(1) \]

好像我們能用整除直接求出 \(a_1\) ，但我們得先證明後面那個和式不影響結果。
由 \(a_k\) 的定義，我們有 \(a_k\leq n-k\)，所以：

\[\begin{aligned} &\ \sum^n_{k=2}(n-k)!a_k\\ \leq &\ \sum^n_{k=2}(n-k)!(n-k)\\ = &\ \sum^n_{k=2}(n-k)!(n-k+1-1)\\ = &\ \sum^n_{k=2}(n-k+1)!-(n-k)!\\ = &\ \sum^{n-2}_{k=0}(k+1)!-k!\\ = &\ (n-1)!-1\\ < &\ (n-1)! \end{aligned} \]

把式 \((1)\) 右側的和式看成帶餘除法的餘數 \(R_1\) ，原式寫為：

\[\mathrm{Index}(p)-1=(n-1)!a_1+R_1\quad(R_1<(n-1)!) \]

把 \((n-1)!\) 看作除數，\(a_1\) 看作商，我們就可以表示出 \(a_1\) 了：

\[a_1=\lfloor\frac{\mathrm{Index}(p)-1}{(n-1)!}\rfloor \]

另外有

\[R_1=(\mathrm{Index}(p)-1)\ \%\ (n-1)! \]

我們繼續拆出餘下和式的第一項：

\[R_1=(n-2)!a_2+\sum^n_{k=3}(n-k)!a_k \]

這個式子和剛剛的是一模一樣的，於是我們能同樣地求出 \(a_2\) 以及其他所有 \(a\) 值，按照以下遞推式：

\[\begin{aligned} a_i&=\lfloor\frac{R_{i+1}}{(n-i)!}\rfloor\\ R_i&=R_{i-1}\ \%\ (n-i)! \end{aligned} \]

這樣，我們就求出了所有 \(a\) ，回想一下， \(a_k\) 是指排列的第 \(k\) 位是排列第 \(k\) 位至 \(n\) 位（未處理的所有位）中第 \(a_k+1\) 小的，我們只需從第 \(1\) 位開始，對第 \(i\) 位，從未選中的數中選擇第 \(a_i+1\) 小的新增到排列中，最後就能形成對應字典序的排列了。
我們可以用一個初始化為[1,2,3,...,n]的vector，每處理一位則erase()一位，每次取索引為 \(a_i\) 的元素即可。

演算法

歸納上述推導過程就有了以下的演算法：

  • 初始化一個 \(1\) 至 \(n\) 的vector<int>vec
  • 求出 \(1\) 至 \(n-1\) 的階乘，備用
  • 初始化 \(R=\mathrm{Index}(p)-1\)
  • 對 \(1\) 至 \(n\) 的 \(i\) 執行：
    • 求出\(a=\lfloor\frac{R_{i+1}}{(n-i)!}\rfloor\)
    • 排列的 \(p\) 的第 \(i\) 位即為vec[a]
    • 移除vec[a]
    • 更新：\(R\leftarrow R\ \%\ (n-i)!\)
  • 排列 \(p\) 即為答案

參考程式碼

題目：P3014 [USACO11FEB]Cow Line S - 洛谷

#include <cctype>
#include <cstdio>
#include <cstring>
#include <iostream>
//---//
#include <algorithm>
#include <cmath>
#include <vector>
using namespace std;

typedef unsigned int u;
typedef long long ll;
typedef unsigned long long llu;

#define rep(i, a, b) for (ll i = a; i < b; i++)
#define REP(i, a, b) for (ll i = a; i <= b; i++)
#define per(i, b, a) for (ll i = b; i >= a; i--)

const ll N = 21;

ll p[N], f[N], d[4 * N];
vector<ll> vec;
bool lz[4 * N];

#define mid ll m = s + ((t - s) >> 1)

void bd(ll s, ll t, ll i) {
  lz[i] = false;
  if (s == t) {
    d[i] = 1;
    return;
  }
  mid;
  bd(s, m, i * 2);
  bd(m + 1, t, i * 2 + 1);
  d[i] = d[i * 2] + d[i * 2 + 1];
}

void spr(ll s, ll t, ll i) {
  if (lz[i]) {
    d[2 * i] = 0;
    d[2 * i + 1] = 0;
    lz[2 * i] = true;
    lz[2 * i + 1] = true;
    lz[i] = false;
  }
}

void upd(ll l, ll r, ll s, ll t, ll i) {
  if (l <= s && t <= r) {
    d[i] = 0;
    lz[i] = true;
    return;
  }
  mid;
  spr(s, t, i);
  if (l <= m) upd(l, r, s, m, 2 * i);
  if (r > m) upd(l, r, m + 1, t, 2 * i + 1);
  d[i] = d[i * 2] + d[i * 2 + 1];
}

ll get(ll l, ll r, ll s, ll t, ll i) {
  if (l > r) return 0;
  if (l <= s && t <= r) return d[i];
  mid;
  spr(s, t, i);
  ll sum;
  if (l <= m) sum = get(l, r, s, m, 2 * i);
  if (r > m) sum = sum + get(l, r, m + 1, t, 2 * i + 1);
  return sum;
}

signed main() {
  ll n, r, q, I, A;
  char c;
  f[0] = 1;
  scanf("%lld%lld\n", &n, &q);
  REP(i, 1, n - 1) f[i] = f[i - 1] * i;
  while (q--) {
    c = getchar();
    while (c != 'P' && c != 'Q') c = getchar();
    if (c == 'P') {
      scanf("%lld", &I);
      I--;
      REP(i, 1, n) vec.push_back(i);
      REP(i, 1, n) {
        A = I / f[n - i];
        I %= f[n - i];
        printf("%lld ", vec[A]);
        vec.erase(vec.begin() + A);
      }
      printf("\n");
    } else {
      r = 1;
      bd(1, n, 1);
      REP(i, 1, n) scanf("%lld", &p[i]);
      REP(k, 1, n) {
        upd(p[k], p[k], 1, n, 1);
        r += get(1, p[k] - 1, 1, n, 1) * f[n - k];
      }
      printf("%lld\n", r);
    }
  }
}

// https://www.luogu.com.cn/record/60986344

參考

相關文章