zoukankan      html  css  js  c++  java
  • UOJ#121. 【NOI2013】向量内积 随机化算法,矩阵

    原文链接www.cnblogs.com/zhouzhendong/UOJ121.html

    前言

    完蛋了我越来越菜了贺题都不会了。

    题解

    (O(n ^ 2 d)) 暴力送 60 分。

    Bitset 优化一下说不定更稳。可能有 85 分。

    来讲正解。

    注意下文中的 "p" 表示原题中的 "k"。

    首先我们来解决一个问题:

    如何在较低的复杂度下判定矩阵 A,B,C 是否满足 (A imes B = C)

    做法是:随机 O(1) 个行向量 (x),判定 (x imes A imes B = x imes C) 是否成立。本题中,行向量大概只需要 10 个左右。

    回到本题。

    考虑 p = 2 的情况。设矩阵 A 为 (n imes d) 的矩阵,且满足 (A_{i,j} = x_{i,j}) ,设矩阵 (B = A ^ {T})

    [self_ i = sum_{j = 1} ^ d A_{i,j} B_{j,i} ]

    那么,如果无解,则

    [A imes B = egin{bmatrix} self_1 & 1 & 1 & 1 & cdots & 1 \ 1 & self_2 & 1 & 1 & cdots & 1 \ 1 & 1 & self_3 & 1 & cdots & 1 \ 1 & 1 & 1 & self_4 & cdots & 1 \ vdots & vdots &vdots &vdots &ddots & vdots\ 1 & 1 & 1 & 1 & cdots & self_n end{bmatrix} ]

    所以,我们可以利用之前提到的方法在 (O(nd)) 的时间复杂度内判定是否有解。如果算出来的行向量的第 (i) 项前后不等,那么第 (i) 个向量和某一个向量的内积一定模p为0。

    当 p = 3 时,考虑到最终矩阵很难确定,因为模 3 意义下可能是 1 或 2。

    但是由于 (1 ^ 2 equiv 2 ^ 2 equiv 1 pmod 3),所以我们可以考虑取平方。

    可以发现

    [left ( sum _k A_{i,k} B_{k,j} ight ) ^ 2 = sum_{k_1} sum_{k_2} left ( A_{i,k_1} A_{i,k_2} ight) left ( B_{i,k_1} B_{i,k_2} ight ) ]

    所以我们可以改写矩阵的定义,设矩阵 (A') 为一个由 ([1cdots n]) 映射到 ((k_1,k_2) | k_1in [1cdots d] , k_2 in [1cdots d]) 的矩阵,矩阵 (B') 为一个由 ((k_1,k_2) | k_1in [1cdots d] , k_2 in [1cdots d]) 映射到 ([1cdots n]) 的矩阵。

    剩下的做法和原先一样。

    时间复杂度为 (O(n d ^ {p-1}))

    代码

    #include <bits/stdc++.h>
    #define clr(x) memset(x,0,sizeof x)
    #define For(i,a,b) for (int i=(a);i<=(b);i++)
    #define Fod(i,b,a) for (int i=(b);i>=(a);i--)
    #define fi first
    #define se second
    #define pb(x) push_back(x)
    #define mp(x,y) make_pair(x,y)
    #define outval(x) cerr<<#x" = "<<x<<endl
    #define outtag(x) cerr<<"---------------"#x"---------------"<<endl
    #define outarr(a,L,R) cerr<<#a"["<<L<<".."<<R<<"] = ";
    						For(_x,L,R)cerr<<a[_x]<<" ";cerr<<endl;
    using namespace std;
    typedef long long LL;
    LL read(){
    	LL x=0,f=0;
    	char ch=getchar();
    	while (!isdigit(ch))
    		f|=ch=='-',ch=getchar();
    	while (isdigit(ch))
    		x=(x<<1)+(x<<3)+(ch^48),ch=getchar();
    	return f?-x:x;
    }
    const int N=100005,D=105;
    int n,d,p;
    int x[N][D],y[N];
    int a[N][D],b[D][N];
    int s[N],t[N],c[D][D];
    int self[N];
    int Solve(int k){
    	For(i,1,n)
    		if (i!=k){
    			int s=0;
    			For(j,1,d)
    				s+=x[i][j]*x[k][j];
    			s%=p;
    			if (s==0){
    				printf("%d %d
    ",min(i,k),max(i,k));
    				return 1;
    			}
    		}
    	return 0;
    }
    int main(){
    	n=read(),d=read(),p=read();
    	For(i,1,n){
    		For(j,1,d)
    			x[i][j]=read()%p;
    		self[i]=0;
    		For(j,1,d)
    			self[i]+=x[i][j]*x[i][j];
    		self[i]%=p;
    		if (p==3)
    			self[i]=self[i]*self[i]%p;
    	}
    	For(i,1,n)
    		For(j,1,d)
    			a[i][j]=b[j][i]=x[i][j];
    	srand(233);
    	For(_,1,10){
    		int sum=0;
    		For(i,1,n)
    			y[i]=rand()%5?1:0,sum+=y[i];
    		if (p==2){
    			For(i,1,n)
    				s[i]=y[i];
    			clr(t);
    			For(i,1,n)
    				For(j,1,d)
    					t[j]+=s[i]*a[i][j];
    			For(i,1,d)
    				t[i]%=p;
    			clr(s);
    			For(i,1,n)
    				For(j,1,d)
    					s[i]+=t[j]*b[j][i];	
    			For(i,1,n)
    				s[i]%=p;
    		}
    		else {
    			For(i,1,n)
    				s[i]=y[i];
    			clr(c);
    			For(i,1,n)
    				For(j1,1,d)
    					For(j2,1,d)
    						c[j1][j2]+=s[i]*a[i][j1]*a[i][j2];
    			For(j1,1,d)
    				For(j2,1,d)
    					c[j1][j2]%=p;
    			clr(s);
    			For(i,1,n)
    				For(j1,1,d)
    					For(j2,1,d)
    						s[i]+=c[j1][j2]*b[j1][i]*b[j2][i];
    			For(i,1,n)
    				s[i]%=p;
    		}
    		For(i,1,n)
    			t[i]=(sum-(1-self[i])*y[i])%p;
    		For(i,1,n)
    			if (s[i]!=t[i])
    				return assert(Solve(i)),0;
    	}
    	puts("-1 -1");
    	return 0;
    }
    
  • 相关阅读:
    [YTU]_2436( C++ 习题 输出日期时间--友元类)
    [YTU]_2435 ( C++ 习题 输出日期时间--友元函数)
    病毒侵袭
    石子合并(区间DP经典例题)
    AC自动机模板2
    【模板】最近公共祖先(LCA)
    华华给月月出题
    线性筛素数
    华华开始学信息学
    华华和月月种树
  • 原文地址:https://www.cnblogs.com/zhouzhendong/p/UOJ121.html
Copyright © 2011-2022 走看看