【特徴多項式解線形繰返し】poj 2118
5491 ワード
叉姐論文:http://www.docin.com/p-724323397.html
a[i]=sigma(a[j]*b[k-j])第n項を求める
標準的な定数線形繰返しは,マトリクス乗算でo(k^3*logn)
しかし,特徴多項式で最適化すれば行列乗算を多項式乗算に置き換えることができ,o(k^2 logn)
まずhamilton-cayleyの定理によれば,一つの行列の特徴多項式はこの行列のゼロ多項式であるため,最高次がk次方項であると仮定し,特徴多項式のk次方項を0~k-1次方項で線形に表すことができ,同理k+1次方項を1~k次方項で線形に表すことができ,k次方項を取り外すことができ,k+1次数項は0~k-1次数項で線形に表出するようにまとめると、行列の任意の次数項は0~k-1次数項で線形に表出できるので、高速べき乗を用いると、o(k^2 logn)の時間で各項の係数を算出することができ、つまり任意の次数の行列はk個の行列がそれぞれベクトルに作用することに等価に加算することができ、そして,このk個の行列が関与する最高項がa[k−1+k−1]項であることに気づき,すなわちo(k^2)の時間で前の2 k−2項を算出できると,このk個のベクトルも決定され,要求されるベクトルも算出できる.
特徴多項式をどのように求めるかについては,この線形繰返し式に対して係数が存在する行ラプラス展開で,各余子式を対角行列の行列式(実際には対角上の値を動かさずに対角行列にすることができる)に容易にすることができ,最後に−1が2回乗算されることに気づく.
一般的な行列乗算に広がる問題について,最大のボトルネックは最後のk個の行列にベクトルをそれぞれ乗じることであり,一般的な繰返し式は前k^2項に及ぶため,初項を求めるだけでo(k^3)であり,一般的な対角行列に似たやり方では,ベース変換行列を求めるにはGauss消元もo(k^3)であることが容易に分かる.
fftで多項式乗算をすれば、最適化を続けることができるようですが、初項をどうするかはまだ分かりません.
hdu4471
線形繰返しにはn個の特殊点があるので,nセグメントに分けて繰返しを行い,tleを直接モーメント乗算するので,各セグメントは特徴多項式繰返しで最適化する.
a[i]=sigma(a[j]*b[k-j])第n項を求める
標準的な定数線形繰返しは,マトリクス乗算でo(k^3*logn)
しかし,特徴多項式で最適化すれば行列乗算を多項式乗算に置き換えることができ,o(k^2 logn)
まずhamilton-cayleyの定理によれば,一つの行列の特徴多項式はこの行列のゼロ多項式であるため,最高次がk次方項であると仮定し,特徴多項式のk次方項を0~k-1次方項で線形に表すことができ,同理k+1次方項を1~k次方項で線形に表すことができ,k次方項を取り外すことができ,k+1次数項は0~k-1次数項で線形に表出するようにまとめると、行列の任意の次数項は0~k-1次数項で線形に表出できるので、高速べき乗を用いると、o(k^2 logn)の時間で各項の係数を算出することができ、つまり任意の次数の行列はk個の行列がそれぞれベクトルに作用することに等価に加算することができ、そして,このk個の行列が関与する最高項がa[k−1+k−1]項であることに気づき,すなわちo(k^2)の時間で前の2 k−2項を算出できると,このk個のベクトルも決定され,要求されるベクトルも算出できる.
特徴多項式をどのように求めるかについては,この線形繰返し式に対して係数が存在する行ラプラス展開で,各余子式を対角行列の行列式(実際には対角上の値を動かさずに対角行列にすることができる)に容易にすることができ,最後に−1が2回乗算されることに気づく.
一般的な行列乗算に広がる問題について,最大のボトルネックは最後のk個の行列にベクトルをそれぞれ乗じることであり,一般的な繰返し式は前k^2項に及ぶため,初項を求めるだけでo(k^3)であり,一般的な対角行列に似たやり方では,ベース変換行列を求めるにはGauss消元もo(k^3)であることが容易に分かる.
fftで多項式乗算をすれば、最適化を続けることができるようですが、初項をどうするかはまだ分かりません.
#include
#include
#include
#include
#include
const int mo=10000;
using namespace std;
struct Coe{
int p[205];
Coe operator *(const Coe &b) const;
};
int ans[205],h[205],N,k,K[205];
int a[205],b[205];
Coe Coe::operator *(const Coe &b) const
{
Coe c;
for (int i=0;i<=k+k;i++) h[i]=0;
for (int i=0;i=k;i--)
for (int j=i-1;j>=i-k;j--)
(h[j]+=h[i]*K[k-(i-j)])%=mo;
for (int i=0;i1) b.p[1]=1;else b.p[0]=K[0];
for (;e;e>>=1) {
if (e&1) sum=sum*b;
b=b*b;
}
return sum;
}
int main()
{
for (;;) {
scanf("%d",&k);
if (!k) break;
memset(a,0,sizeof(a));
memset(b,0,sizeof(b));
for (int i=0;i=0;i--) scanf("%d",&b[i]);
scanf("%d",&N);
if (N=0;i--) {
K[i]=b[i];
// if ((k-1-i+1+1)&1) K[i]=-K[i];
// if ((k-1-i)&1) K[i]=-K[i];
}
Coe B=fgm(N-(k-1));
for (int i=k;i<=k+k;i++) {
a[i]=0;
for (int j=i-1;j>=i-k;j--) (a[i]+=a[j]*b[k-(i-j)])%=mo;
}
for (int i=k-1;i>=0;i--) ans[i]=0;
for (int i=0;i<=k-1;i++)
for (int j=i;j<=i+k-1;j++) (ans[j-i]+=a[j]*B.p[i]%mo)%=mo;
printf("%d
",(ans[k-1]%mo+mo)%mo);
}
return 0;
} hdu4471
線形繰返しにはn個の特殊点があるので,nセグメントに分けて繰返しを行い,tleを直接モーメント乗算するので,各セグメントは特徴多項式繰返しで最適化する.
#include
#include
#include
#include
#include
#define f64 "%lld
"
const int mo=1000000007;
const int lim=100;
using namespace std;
struct PP {
long long p[600];
}b,B;
int n,m,q,t,nk[600],tk[600];
int ck[600][600],c[600],u[600];
long long F[600],f[600],h[600],tmp[600];
int vis[600];
bool cmp(int i,int j)
{
if (nk[i]!=nk[j]) return nk[i]=lim+1;i--)
for (int j=i-1;j>=i-lim;j--)
(h[j]+=h[i]*c[i-j])%=mo;
for (int i=1;i<=lim;i++) a.p[i]=h[i];
}
void fgm(PP &sum,long long e)
{
for (int i=1;i<=lim;i++) sum.p[i]=B.p[i]=0;
sum.p[1]=1;
B.p[2]=1;
for (;e;e>>=1) {
if (e&1) mul(sum,B);
mul(B,B);
}
}
void doit(long long f[],long long e,long long F[])
{
//cout<=1;i--,j--) F[i]=f[j];
return ;
}
for (int i=2*lim+1;i<=2*lim+2*lim;i++) FF(f,i);
for (int i=2*lim,j=2*lim+2*lim; i>=1;i--,j--) f[i]=f[j];
//for (int i=1;i<=2*lim;i++) cout<=lim-i+1;j--,k++) {
tmp[i]=(tmp[i]+b.p[lim-k+1]*f[j])%mo;
}
}
for (int i=1;i<=lim;i++) F[i]=0;
for (int i=1;i<=lim;i++) F[2*lim-i+1]=tmp[i];
}
int main()
{
freopen("input.txt","r",stdin);
freopen("output.txt","w",stdout);
for (int test=1;scanf("%d%d%d",&n,&m,&q)==3;test++) {
for (int i=1;i<=2*lim;i++) vis[i]=0;
printf("Case %d: ",test);
for (int i=1;i<=m;i++) scanf(f64,&f[i]);
scanf("%d",&t);
for (int i=1;i<=t;i++) scanf("%d",&c[i]);
for (int i=t+1;i<=lim;i++) c[i]=0;
for (int i=1;i<=q;i++) {
scanf("%d%d",&nk[i],&tk[i]);
if (nk[i]<=2*lim) vis[nk[i]]=i;
for (int j=1;j<=tk[i];j++)
scanf("%d",&ck[i][j]);
}
++q;
nk[q]=n,tk[q]=t;
for (int i=1;i<=tk[q];i++) ck[q][i]=c[i];
for (int i=1;i<=q;i++) u[i]=i;
sort(u+1,u+q+1,cmp);
int op=m;
if (n<=m) {
printf(f64,f[n]);
continue;
}
for (int i=2*lim,j=m;i>=1;i--,j--)
if (j>0) f[i]=f[j];
else f[i]=0;
for (int i=1;i<=q;i++) {
if (nk[u[i]]<=m) continue;
doit(f,nk[u[i]]-op-1,F);
F[2*lim+1]=0;
for (int j=1;j<=tk[u[i]];j++) (F[2*lim+1]+=ck[u[i]][j]*F[2*lim+1-j])%=mo;
for (int j=1;j<=2*lim;j++) f[j]=F[j+1];
op=nk[u[i]];
//for (int j=2*lim;j>=1;j--) cout<