bzoj3583 杰杰的女性朋友 || bzoj4362 Graph

http://210.33.19.103/problem/2174

很显然是矩阵快速幂的题，设有in和ou矩阵，设in矩阵的转置为in'

显然可以直接暴力求出任意两点间走一步路径条数，然后求其d次幂，但是这样子复杂度不对

注意到设任意两点间走一步路径条数的矩阵为A，那么A=ou*in'，A^d=(ou*in')^d=ou*(in'*ou)^(d-1)*in'（当然d==1时直接特判掉）

in'*ou就是一个K*K的矩阵了，复杂度很对的样子

好像还有点不对。。。题意要求的是前缀和，问题也不大，

设P=in'*ou

题目要求的是ou*P^0*in'+ou*P^1*in'+..+ou*P^(d-1)*in'=ou*(P^0+P^1+..+P^(d-1))*in'

这个P^0+P^1+..+P^(d-1)，搞一个矩阵套矩阵可以算（的确是可以用的23333）

搞两个矩阵X,Y，X=(I,I)，Y是一个2*2的矩阵，第一行为(P,P)，第二行为(0,I)（I表示K阶单位矩阵，0表示K阶零矩阵）

那么求Z=X*Y^(d-1)，Z的第一行第二列的值即为P^0+P^1+..+P^(d-1)

（以上与求值的前缀和的方法是一样的）

明明是思路很清晰的题啊，然而我又是一天被续掉了。。。

鬼畜卡常题啊。。。

此题最好，最正确的方法是不用任何结构体，手写所有矩阵相乘

绝对不能用vector来存不同大小的矩阵!不然会T飞

绝对不要尝试在这类地方用模板类来实现嵌套矩阵！搞不出来的

嵌套矩阵很难写，常数又大，后来我干脆把这个嵌套的矩阵拆开（2*2的矩阵,其中每个元素是K*K的矩阵，就拆成(2*K)*(2*K)的矩阵），目前可以确认是对的

O3卡过代码：

  1 #pragma GCC optimize(3)
  2 #include<cstdio>
  3 #include<algorithm>
  4 #include<cstring>
  5 #include<vector>
  6 #include<cassert>
  7 using namespace std;
  8 #define fi first
  9 #define se second
 10 #define mp make_pair
 11 #define pb push_back
 12 typedef long long ll;
 13 typedef unsigned long long ull;
 14 typedef pair<int,int> pii;
 15 //#define assert(x)
 16 const int md=1000000007;
 17 struct M_NN
 18 {
 19     int d[1010][1010];int x,y;
 20     void resize(int a,int b)
 21     {
 22         x=a;y=b;
 23         memset(d,0,sizeof(d));
 24         //d.clear();
 25         //d.resize(x+1);
 26         //for(int i=1;i<=x;i++)  d[i].resize(y+1);
 27     }
 28 };
 29 struct M_NK
 30 {
 31     int d[1010][42];int x,y;
 32     void resize(int a,int b)
 33     {
 34         x=a;y=b;
 35         memset(d,0,sizeof(d));
 36         //d.clear();
 37         //d.resize(x+1);
 38         //for(int i=1;i<=x;i++)  d[i].resize(y+1);
 39     }
 40 };
 41 struct M_KN
 42 {
 43     int d[42][1010];int x,y;
 44     void resize(int a,int b)
 45     {
 46         x=a;y=b;
 47         memset(d,0,sizeof(d));
 48         //d.clear();
 49         //d.resize(x+1);
 50         //for(int i=1;i<=x;i++)  d[i].resize(y+1);
 51     }
 52 };
 53 struct M_KK
 54 {
 55     int d[42][42];int x,y;
 56     void resize(int a,int b)
 57     {
 58         x=a;y=b;
 59         memset(d,0,sizeof(d));
 60         //d.clear();
 61         //d.resize(x+1);
 62         //for(int i=1;i<=x;i++)  d[i].resize(y+1);
 63     }
 64 };
 65 /*
 66 void init_0(M &p)
 67 {
 68     p.resize(p.x,p.y);
 69 }
 70 */
 71 M_NK operator*(const M_NK &a,const M_KK &b)
 72 {
 73     assert(a.y==b.x);
 74     M_NK c;c.resize(a.x,b.y);
 75     int i,j,k;
 76     for(i=1;i<=a.x;i++)
 77         for(j=1;j<=b.x;j++)
 78             for(k=1;k<=b.y;k++)
 79                 c.d[i][k]=(c.d[i][k]+ll(a.d[i][j])*b.d[j][k]%md)%md;
 80     return c;
 81 }
 82  
 83 M_NN operator*(const M_NK &a,const M_KN &b)
 84 {
 85     assert(a.y==b.x);
 86     M_NN c;c.resize(a.x,b.y);
 87     int i,j,k;
 88     for(i=1;i<=a.x;i++)
 89         for(j=1;j<=b.x;j++)
 90             for(k=1;k<=b.y;k++)
 91                 c.d[i][k]=(c.d[i][k]+ll(a.d[i][j])*b.d[j][k]%md)%md;
 92     return c;
 93 }
 94 M_KK operator*(const M_KN &a,const M_NK &b)
 95 {
 96     assert(a.y==b.x);
 97     M_KK c;c.resize(a.x,b.y);
 98     int i,j,k;
 99     for(i=1;i<=a.x;i++)
100         for(j=1;j<=b.x;j++)
101             for(k=1;k<=b.y;k++)
102                 c.d[i][k]=(c.d[i][k]+ll(a.d[i][j])*b.d[j][k]%md)%md;
103     return c;
104 }
105 M_KK operator*(const M_KK &a,const M_KK &b)
106 {
107     assert(a.y==b.x);
108     M_KK c;c.resize(a.x,b.y);
109     int i,j,k;
110     for(i=1;i<=a.x;i++)
111         for(j=1;j<=b.x;j++)
112             for(k=1;k<=b.y;k++)
113                 c.d[i][k]=(c.d[i][k]+ll(a.d[i][j])*b.d[j][k]%md)%md;
114     return c;
115 }
116 M_KK o,tt2,t1,t2;
117 M_NK ou1;M_KN in1;
118 M_KK p;
119 int in[1010][22],ou[1010][22];
120 M_KK poww(const M_KK &a,int b)
121 {
122     M_KK ans=o,base=a;
123     assert(a.x==a.y);
124     for(;b;base=base*base,b>>=1)
125         if(b&1)
126             ans=ans*base;
127     return ans;
128 }
129  
130 int n,K,m;
131 int solve(int u,int v,int d)
132 {
133     if(d==0)    return u==v;
134     M_KK tt=t1*poww(t2,d-1);
135     int i,j;
136     for(i=1;i<=K;i++)
137     {
138         for(j=1;j<=K;j++)
139         {
140             tt2.d[i][j]=tt.d[i][j+K];
141         }
142     }
143     M_NK t2=ou1*tt2;int an=u==v;
144     for(i=1;i<=t2.y;i++)
145         an=(an+ll(t2.d[u][i])*in1.d[i][v]%md)%md;
146     return an;
147 }
148 int main()
149 {
150     //freopen("/tmp/3583/6.in","r",stdin);
151     //freopen("/tmp/3583/6.ans","w",stdout);
152     int i,j,u,v,d;
153     scanf("%d%d",&n,&K);
154     for(i=1;i<=n;i++)
155     {
156         for(j=1;j<=K;j++)
157         {
158             scanf("%d",&ou[i][j]);
159         }
160         for(j=1;j<=K;j++)
161         {
162             scanf("%d",&in[i][j]);
163         }
164     }
165     ou1.resize(n,K);
166     for(i=1;i<=n;i++)
167     {
168         for(j=1;j<=K;j++)
169         {
170             ou1.d[i][j]=ou[i][j];
171         }
172     }
173     in1.resize(K,n);
174     for(i=1;i<=K;i++)
175     {
176         for(j=1;j<=n;j++)
177         {
178             in1.d[i][j]=in[j][i];
179         }
180     }
181     p=in1*ou1;
182     t1.resize(K,2*K);
183     for(i=1;i<=K;i++)    t1.d[i][i]=t1.d[i][i+K]=1;
184     t2.resize(2*K,2*K);
185     for(i=1;i<=K;i++)
186     {
187         for(j=1;j<=K;j++)
188         {
189             t2.d[i][j]=t2.d[i][j+K]=p.d[i][j];
190         }
191     }
192     for(i=1;i<=K;i++)    t2.d[i+K][i+K]=1;
193     tt2.resize(K,K);
194     o.resize(2*K,2*K);
195     for(i=1;i<=2*K;i++)  o.d[i][i]=1;
196     scanf("%d",&m);
197     while(m--)
198     {
199         scanf("%d%d%d",&u,&v,&d);
200         printf("%d\n",solve(u,v,d));
201         //return 0;
202     }
203     return 0;
204 }
View Code

一份有分，可以跑的代码：

  1 #pragma GCC optimize(3)
  2 #include<cstdio>
  3 #include<algorithm>
  4 #include<cstring>
  5 #include<vector>
  6 //#include<cassert>
  7 using namespace std;
  8 #define fi first
  9 #define se second
 10 #define mp make_pair
 11 #define pb push_back
 12 typedef long long ll;
 13 typedef unsigned long long ull;
 14 typedef pair<int,int> pii;
 15 #define vector myvec
 16 template<typename T>
 17 struct myvec
 18 {
 19     T *p;int sz;
 20     myvec():p(0),sz(0){}
 21     ~myvec(){delete[] p;}
 22     myvec &operator=(const myvec &b)
 23     {
 24         delete[] p;
 25         sz=b.sz;
 26         p=new T[sz];
 27         for(T *i1=p,*i2=b.p;i1!=p+sz;++i1,++i2) *i1=*i2;
 28         return *this;
 29     }
 30     myvec(const myvec &b):p(0),sz(0){*this=b;}
 31     void resize(int d){sz=d;delete[] p;p=new T[d]();}
 32     T &operator[](int d){return p[d];}
 33     const T &operator[](int d)const{return p[d];}
 34 };
 35 const ll md=1000000007;
 36 struct M1
 37 {
 38     vector<vector<ll> > d;int x,y;
 39     M1();
 40     void resize(int xx,int yy)
 41     {
 42         x=xx;y=yy;
 43         d.resize(x+1);
 44         for(int i=1;i<=x;i++)    d[i].resize(y+1);
 45     }
 46 };
 47 struct M2
 48 {
 49     M1 d[3][3];int x,y;
 50     M2();
 51     void resize_sub(int xx,int yy)
 52     {
 53         int i,j;
 54         for(i=1;i<=x;i++)
 55             for(j=1;j<=y;j++)
 56                 d[i][j].resize(xx,yy);
 57     }
 58 };
 59 void init_0(M1 &p)
 60 {
 61     //p.d.clear();
 62     p.resize(p.x,p.y);
 63 }
 64 void init_0(M2 &p)
 65 {
 66     int i,j;
 67     for(i=1;i<=2;i++)
 68         for(j=1;j<=2;j++)
 69             init_0(p.d[i][j]);
 70 }
 71 void init_1(M1 &p)
 72 {
 73     init_0(p);//assert(p.x==p.y);
 74     for(int i=1;i<=p.x;i++)  p.d[i][i]=1;
 75 }
 76 void init_1(M2 &p)
 77 {
 78     init_0(p);//assert(p.x==p.y);
 79     for(int i=1;i<=p.x;i++)  init_1(p.d[i][i]);
 80 }
 81 M1::M1(){x=y=0;init_0(*this);}
 82 M2::M2(){x=y=0;init_0(*this);}
 83 M1 operator+(const M1 &a,const M1 &b)
 84 {
 85     //assert(a.x==b.x&&a.y==b.y);
 86     M1 c;c.resize(a.x,a.y);
 87     int i,j;
 88     for(i=1;i<=a.x;i++)
 89         for(j=1;j<=a.y;j++)
 90             c.d[i][j]=(a.d[i][j]+b.d[i][j])%md;
 91     return c;
 92 }
 93 M1 operator*(const M1 &a,const M1 &b)
 94 {
 95     //assert(a.y==b.x);
 96     M1 c;c.resize(a.x,b.y);
 97     int i,j,k;
 98     for(i=1;i<=a.x;i++)
 99         for(j=1;j<=b.x;j++)
100             for(k=1;k<=b.y;k++)
101                 c.d[i][k]=(c.d[i][k]+a.d[i][j]*b.d[j][k])%md;
102     return c;
103 }
104 M2 operator*(const M2 &a,const M2 &b)
105 {
106     //assert(a.y==b.x);
107     M2 c;c.x=a.x;c.y=b.y;c.resize_sub(a.d[1][1].x,a.d[1][1].y);
108     int i,j,k;
109     for(i=1;i<=a.x;i++)
110         for(j=1;j<=b.x;j++)
111             for(k=1;k<=b.y;k++)
112                 c.d[i][k]=(c.d[i][k]+a.d[i][j]*b.d[j][k]);
113     return c;
114 }
115 M1 ou1,in1,p;
116 int ou[1010][22],in[1010][22];
117 M2 poww(const M2 &a,ll b)
118 {
119     //printf("tt%d %d\n",a.x,a.y);
120     M2 ans,base=a;/*assert(a.x==a.y);*/ans.x=ans.y=a.x;
121     ans.resize_sub(a.d[1][1].x,a.d[1][1].y);init_1(ans);
122     for(;b;base=base*base,b>>=1)
123         if(b&1)
124             ans=ans*base;
125     return ans;
126 }
127 int n,K,m;
128  
129 /*
130 void out(const M1 &a)
131 {
132             puts("start");
133                         int i,j;
134                                         for(i=1;i<=a.x;i++)
135                                                                     {
136                                                                                                             for(j=1;j<=a.y;j++)
137                                                                                                                                                                 {
138                                                                                                                                                                                                                                     printf("%lld ",a.d[i][j]);
139                                                                                                                                                                                                                                                                                                             }
140                                                                                                                                                         puts("");
141                                                                                                                                                                                                         }
142                                                             puts("end");
143                                                                                 fflush(stdout);
144 }
145  
146  
147 void out(const M2 &p)
148 {
149     puts("st");
150     printf("xy%d %d\n",p.x,p.y);
151     for(int i=1;i<=p.x;i++)
152         for(int j=1;j<=p.y;j++)
153             out(p.d[i][j]),puts("");
154     puts("ed");
155 }
156 */
157 M2 t1,t2;
158 int solve(int u,int v,int d)
159 {
160     if(d==0)    return u==v;
161     return ((u==v)+(ou1*((t1*poww(t2,d-1)).d[1][2])*in1).d[u][v])%md;
162 }
163 int main()
164 {
165     int i,j,u,v,d;
166     scanf("%d%d",&n,&K);
167     for(i=1;i<=n;i++)
168     {
169         for(j=1;j<=K;j++)
170         {
171             scanf("%d",&ou[i][j]);
172         }
173         for(j=1;j<=K;j++)
174         {
175             scanf("%d",&in[i][j]);
176         }
177     }
178     ou1.resize(n,K);
179     for(i=1;i<=n;i++)
180     {
181         for(j=1;j<=K;j++)
182         {
183             ou1.d[i][j]=ou[i][j];
184         }
185     }
186     in1.resize(K,n);
187     for(i=1;i<=K;i++)
188     {
189         for(j=1;j<=n;j++)
190         {
191             in1.d[i][j]=in[j][i];
192         }
193     }
194     p=in1*ou1;
195     t1.x=1;t1.y=2;t1.resize_sub(K,K);
196     //t1.d[1][1]=M1(1,K,K);
197     init_1(t1.d[1][1]);init_1(t1.d[1][2]);
198     t2.x=t2.y=2;t2.resize_sub(K,K);
199     t2.d[1][1]=t2.d[1][2]=p;
200     //t2.d[2][2]=M1(1,K,K);
201     init_1(t2.d[2][2]);
202     //out((t1*poww(t2,0)).d[1][2]);
203     scanf("%d",&m);
204     while(m--)
205     {
206         scanf("%d%d%d",&u,&v,&d);
207         printf("%d\n",solve(u,v,d));
208     }
209     return 0;
210 }
View Code

 

posted @ 2018-09-07 20:39  hehe_54321  阅读(380)  评论(0)  编辑  收藏  举报
AmazingCounters.com