[APIO2022]T3 perm

考虑增量构造这个排列，设当前 LIS 个数为 \(cnt\)，一开始排列为空，\(cnt =1\)。

如果在序列末尾放上一个比之前所有数都大的数，那么新的 \(cnt'=2\times cnt\) 。

如果在序列末尾放上一个比之前所有数都小的数，那么新的 \(cnt'=cnt+1\) 。

这启发我们按 \(k\) 二进制的数位构造排列，除了最高位的 1 外，其它位上的所有 1 需要添加两个数，每一个 0 都需要添加一个数，总的排列长度就是 \(\lfloor \log_2 k \rfloor + \text{popcount}(k)-1\leq 118\) ，只能获得 91 分。

接下来的我们发现难以找到一个如刚才那个做法一样“简洁”的构造，于是思考如何乱搞。

考虑把排列 \(b\) 接到排列 \(a\) 的后面，再给排列 \(b\) 中每一个数都加上排列 \(a\) 的长度，不难发现新构造出的排列 LIS 个数为 \(cnt_a\times cnt_b\) 。

那么如果我们能将 \(k\) 分解为 \(\prod p_i\) 的形式，使得 \(\sum (\lfloor \log_2 p_i \rfloor + \text{popcount}(p_i)-1 )\leq 90\) ，就可以使用如上做法对每一个 \(p_i\) 构造一个排列。

事实上，对于一个足够大合数 \(k\) 尽管能构造出 \(k\) 的二进制数位上有很多 1 ，但对于它的因子们，数位上 1 的个数十分随机，满足题意的构造又特别多，因此对于许多 \(k\) 该做法总能找到一个合法的分解。

然而，对于有一些 \(k\) ，尤其是 \(k\) 是质数的情况，无法找到合法的分解。虽然如此，我们可以发现存在合法分解的 \(k\) 十分密集。

于是我们对于一个很小的 \(\delta\) ，找出 \(k-\delta\) 的合法分解并构造排列，再在序列末尾放上 \(\delta\) 个比前面所有数都小的数即可。

实现上使用了 Pollard Rho 算法对于 \(\delta\leq 30\) 随机地分解 \(k-\delta\) ，不用怎么调参，基本能够稳定迅速得出答案。

#include "perm.h"
#include <cstdio>
#include <random>
typedef long long ll;
typedef std::vector<int> vi;
std::mt19937 rng(std::random_device{}());
vi a;
int cur;
void con(ll k){
	int t;
	for(t=1;(1ll<<t)<=k;++t);
	int s=--t;
	k-=(1ll<<t);
	for(int i=0;i<=t;++i) s+=k>>i&1;
	int tmp=s;
	for(int i=0;i<t;++i){
		if(k>>i&1) a.emplace_back(--s+cur);
		a.emplace_back(i+cur);
	}
	if(k>>t&1) a.emplace_back(--s+cur);
	cur+=tmp;
}
int calc(ll k){
	return __builtin_popcountll(k)-__builtin_clzll(k)+62;
}
namespace PRs{
	ll gcd(ll a,ll b){
		if (!b) return a;
		return gcd(b,a%b);
	}
	ll abs(ll a){if (a<0) return -a;return a;}
	ll qt(ll a,ll b,ll p){return (a*b-(ll)((long double)a/p*b)*p+p)%p;}
	ll qp(ll a,ll b,ll p){
		ll res=1;
		while (b){
			if (b&1) res=qt(res,a,p);
			a=qt(a,a,p);
			b>>=1;
		}
		return res;
	}
	bool MR(ll o){
		if (o<3||o%2==0) return o==2;
		if(o%3==0) return o==3;
		if(o%5==0) return o==5;
		ll a=o-1; int t=0;
		while ((a&1)^1) a>>=1,++t;
		for (int ite=0; ite<4; ++ite){
			ll x=rng()%(o-2)+2,v=qp(x,a,o);int s;
			if (v==1) continue;
			for (s=0; s<t; ++s){
				if (v==o-1) break;
				v=qt(v,v,o);
			}
			if (s==t) return 0;
		}
		return 1;
	}
	ll f(ll a,ll c,ll p){return (qt(a,a,p)+c)%p;}
	ll PR(ll o){
		ll s=0,t=0;
		ll c=rng()%(o-1)+1,vl=1;
		for (int go=1; ; go<<=1,s=t,vl=1){
			for (int st=1; st<=go; ++st){
				t=f(t,c,o);
				vl=qt(vl,abs(t-s),o);
				if (st%63==0){
					ll d=gcd(vl,o);
					if (d>1) return d;
				}
			}
			ll d=gcd(vl,o);
			if (d>1) return d;
		}
	}
}
using PRs::PR;
using PRs::MR;
ll p[1000];int rk,tmp;
void factor(ll x){
	if(MR(x)){p[++rk]=x;return;}
	while(1){
		ll t=PR(x);
		if(t==x) continue;
		if(rng()&2){
			p[++rk]=t;p[++rk]=x/t;
			return;
		}
		else factor(t),factor(x/t);
	}
}
vi construct_permutation(ll k){
	a.clear();tmp=0;cur=0;
	if(calc(k)<=90){
		con(k);
		for(int i=tmp-1;~i;--i) a.emplace_back(i);
		return a;
	}
	while(1){
		while(1){
			int lim=90-tmp;
			rk=0;factor(k);
			for(int i=1;i<=rk;++i) lim-=calc(p[i]);
			if(lim>=0){
				cur=tmp;
				for(int i=1;i<=rk;++i) con(p[i]);
				for(int i=tmp-1;~i;--i) a.emplace_back(i);
				return a;
			}
			--k;++tmp;
			if(tmp>30) break;
		}
		k+=tmp;tmp=0;
	}
	return a;
}