奇怪模板集合

Posted 2021-03-01 ldxcaicai

传说中 (O(n^{2.807})) 的矩阵乘法模板

实测大数据会比 (O(n^3)) 的矩乘快一点，卡常卡不过可以拿来用，并且博主写的比较丑，所以如果优化一下写法大概能跑得更快。

Code

#include<bits/stdc++.h>
#define ri register int
using namespace std;
typedef long long ll;
typedef vector<int> poly;
typedef pair<int,int> pii;
typedef pair<ll,int> pli;
typedef pair<int,ll> pil;
typedef pair<ll,ll> pll;
#define fi first
#define se second
#define pb push_back
#define rez resize
const ll Inf=1e18;
const int rlen=1<<20;
char buf[rlen],*ib=buf,*ob=buf;
#define gc() (((ib==ob)&&(ob=(ib=buf)+fread(buf,1,rlen,stdin)),ib==ob)?-1:*ib++)
inline int read() {
  static int ans;
  static char ch;
  ans=0,ch=gc();
  while(!isdigit(ch)) ch=gc();
  while(isdigit(ch)) ans=(ans<<3)+(ans<<1)+(ch^48),ch=gc();
  return ans;
}
inline ll readl() {
  static ll ans;
  static char ch;
  ans=0,ch=gc();
  while(!isdigit(ch)) ch=gc();
  while(isdigit(ch)) ans=(ans<<3)+(ans<<1)+(ch^48),ch=gc();
  return ans;
}
namespace modular {
  const int mod=1e9+7;
  inline int add(int a,int b) { return a+b<mod?a+b:a+b-mod; }
  inline int dec(int a,int b) { return a<b?a-b+mod:a-b; }
  inline int mul(int a,int b) { return (ll)a*b%mod; }
  inline void Add(int&a,int b) { a=a+b<mod?a+b:a+b-mod; }
  inline void Dec(int&a,int b) { a=a<b?a-b+mod:a-b; }
  inline void Mul(int&a,int b) { a=(ll)a*b%mod; }
  inline int ksm(int a,int p) { int ret=1;for(;p;p>>=1,Mul(a,a)) (p&1)&&(Mul(ret,a),1);return ret; }
  inline int Inv(int a) { return ksm(a,mod-2); }
  inline int sqr(int a) { return (ll)a*a%mod; }
  inline int cub(int a) { return (ll)a*a%mod*a%mod; }
} using namespace modular;
template<typename T> inline void ckmax(T&a,T b) { a<b?a=b:0; }
template<typename T> inline void ckmin(T&a,T b) { a>b?a=b:0; }
struct mat {
  int a[512][512],l;
  mat(int _l=0) { l=_l;for(ri i=0;i<l;++i) for(ri j=0;j<l;++j) a[i][j]=0; }
  inline int*operator[](const int&k) { return a[k]; }
  inline mat operator*(mat b) {
    mat c(l);
    if(l<=64) {
      for(ri i=0;i<l;++i) for(ri k=0;k<l;++k) if(a[i][k])
      for(ri j=0;j<l;++j) Add(c[i][j],mul(a[i][k],b[k][j]));
    }
    else {
      l>>=1;
      mat aa(l),bb(l),m1(l),m2(l),m3(l),m4(l),m5(l),m6(l),m7(l);
      //m1=(a11+a22)*(b11+b22);
      for(ri i=0;i<l;++i) for(ri j=0;j<l;++j) {
        aa[i][j]=add(a[i][j],a[i+l][j+l]);
        bb[i][j]=add(b[i][j],b[i+l][j+l]);
      } m1=aa*bb;
      //m2=(a21+a22)*b11;
      for(ri i=0;i<l;++i) for(ri j=0;j<l;++j) {
        aa[i][j]=add(a[i+l][j],a[i+l][j+l]);
        bb[i][j]=b[i][j];
      } m2=aa*bb;
      //m3=a11*(b12-b22);
      for(ri i=0;i<l;++i) for(ri j=0;j<l;++j) {
        aa[i][j]=a[i][j];
        bb[i][j]=dec(b[i][j+l],b[i+l][j+l]);
      } m3=aa*bb;
      //m4=a22*(b21-b11);
      for(ri i=0;i<l;++i) for(ri j=0;j<l;++j) {
        aa[i][j]=a[i+l][j+l];
        bb[i][j]=dec(b[i+l][j],b[i][j]);
      } m4=aa*bb;
      //m5=(a11+a12)*b22;
      for(ri i=0;i<l;++i) for(ri j=0;j<l;++j) {
        aa[i][j]=add(a[i][j],a[i][j+l]);
        bb[i][j]=b[i+l][j+l];
      } m5=aa*bb;
      //m6=(a21-a11)*(b11+b12);
      for(ri i=0;i<l;++i) for(ri j=0;j<l;++j) {
        aa[i][j]=dec(a[i+l][j],a[i][j]);
        bb[i][j]=add(b[i][j],b[i][j+l]);
      } m6=aa*bb;
      //m7=(a12-a22)*(b21+b22);
      for(ri i=0;i<l;++i) for(ri j=0;j<l;++j) {
        aa[i][j]=dec(a[i][j+l],a[i+l][j+l]);
        bb[i][j]=add(b[i+l][j],b[i+l][j+l]);
      } m7=aa*bb;
      // c11=m1+m4-m5+m7
      // c12=m3+m5
      // c21=m2+m4
      // c22=m1-m2+m3+m6
      for(ri i=0;i<l;++i) for(ri j=0;j<l;++j) {
        if((c[i][j]=((ll)m1[i][j]+m4[i][j]-m5[i][j]+m7[i][j])%mod)<0) c[i][j]+=mod;
        c[i+l][j]=add(m2[i][j],m4[i][j]);
        c[i][j+l]=add(m3[i][j],m5[i][j]);
        if((c[i+l][j+l]=((ll)m1[i][j]-m2[i][j]+m3[i][j]+m6[i][j])%mod)<0) c[i+l][j+l]+=mod;
      }
      l<<=1;
    }
    return c;
  }
  inline void print() { for(ri i=0;i<l;++i,puts("")) for(ri j=0;j<l;++j) cout<<a[i][j]<<‘ ‘; }
};
int main() {
  mat a(512);
  for(ri i=0;i<512;++i) for(ri j=0;j<512;++j) a[i][j]=i*256+j;
  a=a*a;
  cout<<a[0][0];
  return 0;
}