TopCoder SRM452 Div1Hard

問題

10進法で読んで、広義単調増加になる数を「increasing number」と呼ぶ。例えば1234、8899とかはincreasing numberだけど、314、2009とかは違う。
digits桁のincreasing number(以下、IS)で、divisorの倍数であるものの個数を、1,000,000,007で割った余りを求めよ。
digitsは1以上10^18以下、divisorは1以上500以下である。

解法

言うまでもなく、全探査では解けない。というかO(digits)でも無理。
とりあえず、ISを別な方法で表現することを考える。左から順に1,2,…,9がそれぞれ0個以上並んだものと考えると、右から9がv1個、8がv2-v1個、…、1がv9-v8個並んだものと考えることができ、これを計算すると、1/9*(10D - 9 + 10v1 + … + 10v8)となる。但し、D >= v8 >= … >= v1 >= 0 である。
これがdivisor(以下K)の倍数であるという条件がある。合同式評価において1/9とかいうのはうざいので取っ払うと、
(10D - 9 + 10v1 + … + 10v8) ≡ 0 (mod 9K)と書き換えられる。
右側の10v1 + … の部分に注目すると、DPで計算できそうに思えてくる。そのために、予め0<=N<9Kに対して、10X ≡ N (mod 9K)なる0<=X<=Dの数を求めておく(これは、10Xの周期を用いるとO(K)で実行できる)。
次に、10^…の部分において、v8 >= … >= v1 という大変うざい条件が付いているので、「v1, …, v8のセットを持ってきて、適切に並べる」ことにする。すると、この「セット」の数を数えればよい。
v1…の部分は、modで考えてしまう。各modごとに何個可能なものがあるかはさっき求めた。「どのmodまでやったか」「今和のmodはいくつか」「何個v1みたいのを使ったか」をキーにDPを行うと解ける。

コード

#include <cstdio>

using namespace std;

#define MOD 1000000007

class IncreasingNumber
{
	int modpow(int V, long long P, int M)
	{
		if(P==0) return 1;
		long long tmp = modpow(V, P/2, M);
		tmp = (tmp*tmp)%M;
		if(P%2==1) tmp = (tmp*V)%M;
		return tmp;
	}
	
	long long modc[4500];
	
	void make_cyc(long long D, int K)
	{
		//10^0, 10^1, ... , 10^D mod K
		int vp[4500], now;
		now = 1;
		for(int i=0;i<K;i++){
			vp[i] = -1;
			modc[i] = 0;
		}
		
		int pos = 0;
		for(;pos<=D;){
			if(vp[now]==-1){
				vp[now] = pos++;
				modc[now]++;
				now = (now*10)%K;
			}else{
				int cyc = pos - vp[now];
				D -= pos-1;
				for(int i=0;i<cyc;i++){
					modc[now] += D/cyc + ((D%cyc>i)?1:0);
					now = (now*10)%K;
				}
				break;
			}
		}
	}
	
public:
	int dp[2][4500][9];
	
	long long GCD(long long p, long long q)
	{
		if(q==0) return p;
		return GCD(q, p%q);
	}
	
	long long C(long long p, int q)
	{
		long long ret = 1;
		int used[8];
		for(int i=0;i<q;i++) used[i] = i+1;
		for(int i=0;i<q;i++){
			long long t = p-i;
			for(int j=0;j<q;j++){
				long long g = GCD(t, used[j]);
				t /= g;
				used[j] /= g;
			}
			t %= MOD;
			ret = (ret*t)%MOD;
		}
		return ret;
	}
	
	int countNumbers(long long D, int K)
	{
		make_cyc(D, K*9);
		
		int t = 0;
		for(int i=0;i<K*9;i++){
			for(int j=0;j<9;j++) dp[t][i][j] = 0;
		}
		dp[t][0][0] = 1;
		
		for(int i=0;i<K*9;i++) if(modc[i]){
			for(int j=0;j<K*9;j++){
				for(int k=0;k<9;k++) dp[1-t][j][k] = dp[t][j][k];
			}
			
			long long C_tbl[9];
			for(int j=1;j<9;j++) C_tbl[j] = C(modc[i]+j-1, j);
			for(int j=0;j<K*9;j++)
				for(int k=0;k<9;k++) if(dp[t][j][k]){
					for(int l=1;l<9-k;l++){
						dp[1-t][(j+i*l)%(9*K)][k+l] += (int)((long long)C_tbl[l] * (long long)dp[t][j][k] % MOD);
						dp[1-t][(j+i*l)%(9*K)][k+l] %= MOD;
					}
				}
				
			t = 1-t;
		}
		
		printf("%d\n", ((9*K)-modpow(10, D, 9*K)+9) % (9*K));
		return dp[t][((9*K)-modpow(10, D, 9*K)+9) % (9*K)][8];
	}
};