挺好的數位dp……
先說一下我個人的做法:
經過觀察,發現這題按照以往的思路從後往前遞增,不怎麽好推,然後我就大膽猜想,從前往後推,發現很好推啊,維護四個變量,從開始位置到現在有了i個數
f[i]:所有數的所有未包含最後一位的子串的和
s[i]:所有數的所有後綴子串的和
c[i]:所有數的所有後綴子串的個數
n[i]:所有數共有多少個
他們的轉移依次是(k為進制數)
f[i]=f[i-1]*k+s[i-1]*k
s[i]=s[i-1]*k*k+c[i-1]*k*(k-1)/2+n[i-1]*k*(k-1)/2
c[i]=c[i-1]*k+n[i-1]*k
n[i]=n[i-1]*k
我們發現對於最高位低於上界的數,我們可以在確定最高位上是1~9之後用上面的轉移一遍O(n)dp算出來.如果最高位等於上界的話,我們的轉移不太一樣,但是也只不過是把某些k改為了這一位的上屆,而且如果本位未達到上屆,往後轉移還是老樣子,然而每次都要從前往後走一遍,會T,不過,這很明顯是個可以用矩陣乘法優化的dp,因為他的轉移方式每次都一樣,所以我們就可以加速了,然而這是4*4的矩陣再加上一個log,吃不消啊,但是我們可以預處理轉移i(1<=i<=max(n,m))次的矩陣,這樣就可以做到O(4^3*n)了,又因為這個矩陣是個上三角矩陣,所以我們加一些矩陣乘法時的優化就可以有有著一個10左右常數的O(n)的做法了,我們解決了這道題!!!

#include <cstdio>
#include <cstring>
#include <algorithm>
char xB[(1<<15
 
)+10],*xS,*xT;
#define gtc (xS==xT&&(xT=(xS=xB)+fread(xB,1,1<<15,stdin),xS==xT)?0:*xS++)
template <typename _t>
inline void read(_t &x){
  register char ch=gtc;bool ud=false;
  for(x=0;ch<‘0‘||ch>‘9‘;ch=gtc)if(ch==‘-‘)ud=true;
  for(;ch>=‘0‘&&ch<=‘9‘;x=(x<<1
 
)+(x<<3)+ch-‘0‘,ch=gtc);
  if(ud)x=-x;
}
typedef long long LL;
const int P=20130427;
const int N=100010;
int a[4][4],b[4],s[N][4][4],temp_a[4][4],temp_b[4],c[4],d[4];
inline void get(int x[][4],int y){
  memset(temp_a,0,sizeof(a));
  register int i,j,k;
  for(i=0;i<4;++i)
    for(j=0;j<4;++j)
      if(a[i][j])
        for(k=0;k<4;++k)
          if(x[j][k])
            temp_a[i][k]=(temp_a[i][k]+(LL)x[j][k]*a[i][j])%P;
  memcpy(s[y],temp_a,sizeof(s[y]));
}
inline void run(int x[][4]){
  memset(temp_b,0,sizeof(temp_b));
  register int i,j;
  for(i=0;i<4;++i)
    for(j=0;j<4;++j)
      if(x[i][j])
        temp_b[i]=(temp_b[i]+(LL)x[i][j]*d[j])%P;
  memcpy(c,temp_b,sizeof(c));
}
int bit,digit[N],k,n,m,len;
inline int calc(){
  int ans=0,i;
  d[3]=0,d[2]=(LL)k*(k-1)/2%P,d[1]=k-1,d[0]=k-1;
  for(i=1;i<bit;++i)
    run(s[i-1]),ans=(ans+c[3]+c[2])%P;
  memset(b,0,sizeof(b)),b[0]=1;
  for(i=bit;i>0;--i){
    d[0]=((LL)b[0]*(digit[i]-(i==bit))%P);
    d[1]=((LL)b[1]*(digit[i]-(i==bit))+d[0])%P;
    d[2]=((LL)k*b[2]%P*(digit[i]-(i==bit))+(LL)b[1]*((LL)digit[i]*(digit[i]-1)/2%P)+(LL)b[0]*((LL)digit[i]*(digit[i]-1)/2%P))%P;
    d[3]=((LL)b[3]*(digit[i]-(i==bit))+(LL)b[2]*(digit[i]-(i==bit)))%P;
    run(s[i-1]);
    ans=(ans+c[3]+c[2])%P;
    b[3]=(b[3]+b[2])%P;
    b[2]=((LL)k*b[2]+(LL)(b[1]+b[0])*digit[i])%P;
    ++b[1];
  }
  return (ans+b[3]+b[2])%P;
}
int main(){
  read(k);int i,j,ans=0;
  a[3][3]=k,a[3][2]=k;
  a[2][2]=(LL)k*k%P,a[2][1]=((LL)k*(k-1)/2)%P,a[2][0]=((LL)k*(k-1)/2)%P;
  a[1][1]=k,a[1][0]=k;
  a[0][0]=k;
  s[0][0][0]=s[0][1][1]=s[0][2][2]=s[0][3][3]=1;
  for(read(n),i=n;i>0;--i)read(digit[i]);
  read(m),len=std::max(n,m);
  for(i=1;i<=len;++i)get(s[i-1],i);
  if(n==1)ans=(ans-(LL)digit[1]*(digit[1]-1)/2%P+P)%P;
  else{
    for(--digit[1],i=1;i<=n;++i)
      if(digit[i]<0)digit[i]+=k,--digit[i+1];
      else break;
    while(digit[n]==0)--n;
    bit=n,ans=(ans-calc()+P)%P;
  }
  for(i=m;i>0;--i)read(digit[i]);
  if(m==1)ans=(ans+(LL)digit[1]*(digit[1]+1)/2%P)%P;
  else bit=m,ans=(ans+calc())%P;
  printf("%d\n",ans);
  return 0;
}

【BZOJ 3326】[Scoi2013]數數 數位dp+矩陣乘法優化