线性代数
程序员文章站
2022-07-12 14:50:47
...
线性代数
矩阵
矩阵的逆
的逆矩阵是使得的矩阵。
逆矩阵可以用高斯消元的方式来求。
矩阵乘法
常用矩阵快速幂来解决线性递推式。
优化:(降常数)
// 以下文的参考代码为例
inline mat operator*(const mat& T) const {
mat res;
for (int i = 0; i < sz; ++i)
for (int j = 0; j < sz; ++j)
for (int k = 0; k < sz; ++k) {
res.a[i][j] += mul(a[i][k], T.a[k][j]);
res.a[i][j] %= MOD;
}
return res;
}
// 不如
inline mat operator*(const mat& T) const {
mat res;
int r;
for (int i = 0; i < sz; ++i)
for (int k = 0; k < sz; ++k) {
r = a[i][k];
for (int j = 0; j < sz; ++j) res.a[i][j] += T.a[k][j] * r;
res.a[i][j] %= MOD;
}
return res;
}
线性递推式
除了矩阵快速幂外介绍两个方法:
一、BM算法
#include <bits/stdc++.h>
using namespace std;
#define rep(i,a,n) for (long long i=a;i<n;i++)
#define per(i,a,n) for (long long i=n-1;i>=a;i--)
#define pb push_back
#define mp make_pair
#define all(x) (x).begin(),(x).end()
#define fi first
#define se second
#define SZ(x) ((long long)(x).size())
typedef vector<long long> VI;
typedef long long ll;
typedef pair<long long,long long> PII;
const ll mod=1e9+7;
ll powmod(ll a,ll b) {
ll res=1;a%=mod; assert(b>=0);
for(;b;b>>=1){
if(b&1)res=res*a%mod;a=a*a%mod;
}
return res;
}
// head
long long _,n;
namespace linear_seq
{
const long long N=10010;
ll res[N],base[N],_c[N],_md[N];
vector<long long> Md;
void mul(ll *a,ll *b,long long k)
{
rep(i,0,k+k) _c[i]=0;
rep(i,0,k) if (a[i]) rep(j,0,k)
_c[i+j]=(_c[i+j]+a[i]*b[j])%mod;
for (long long i=k+k-1;i>=k;i--)
if (_c[i])
rep(j,0,SZ(Md))
_c[i-k+Md[j]]=(_c[i-k+Md[j]]-_c[i]*_md[Md[j]])%mod;
rep(i,0,k) a[i]=_c[i];
}
long long solve(ll n,VI a,VI b)
{ // a 系数 b 初值 b[n+1]=a[0]*b[n]+...
// printf("%d\n",SZ(b));
ll ans=0,pnt=0;
long long k=SZ(a);
assert(SZ(a)==SZ(b));
rep(i,0,k) _md[k-1-i]=-a[i];_md[k]=1;
Md.clear();
rep(i,0,k) if (_md[i]!=0) Md.push_back(i);
rep(i,0,k) res[i]=base[i]=0;
res[0]=1;
while ((1ll<<pnt)<=n) pnt++;
for (long long p=pnt;p>=0;p--)
{
mul(res,res,k);
if ((n>>p)&1)
{
for (long long i=k-1;i>=0;i--)
res[i+1]=res[i];res[0]=0;
rep(j,0,SZ(Md))
res[Md[j]]=(res[Md[j]]-res[k]*_md[Md[j]])%mod;
}
}
rep(i,0,k) ans=(ans+res[i]*b[i])%mod;
if (ans<0) ans+=mod;
return ans;
}
VI BM(VI s)
{
VI C(1,1),B(1,1);
long long L=0,m=1,b=1;
rep(n,0,SZ(s))
{
ll d=0;
rep(i,0,L+1) d=(d+(ll)C[i]*s[n-i])%mod;
if (d==0) ++m;
else if (2*L<=n)
{
VI T=C;
ll c=mod-d*powmod(b,mod-2)%mod;
while (SZ(C)<SZ(B)+m) C.pb(0);
rep(i,0,SZ(B)) C[i+m]=(C[i+m]+c*B[i])%mod;
L=n+1-L; B=T; b=d; m=1;
}
else
{
ll c=mod-d*powmod(b,mod-2)%mod;
while (SZ(C)<SZ(B)+m) C.pb(0);
rep(i,0,SZ(B)) C[i+m]=(C[i+m]+c*B[i])%mod;
++m;
}
}
return C;
}
long long gao(VI a,ll n)
{
VI c=BM(a);
c.erase(c.begin());
rep(i,0,SZ(c)) c[i]=(mod-c[i])%mod;
return solve(n,c,VI(a.begin(),a.begin()+SZ(c)));
}
};
int main()
{
vector<int>v;
v.push_back(3);
v.push_back(9);
v.push_back(20);
v.push_back(46);
v.push_back(106);
v.push_back(244);
v.push_back(560);
v.push_back(1286);
v.push_back(2956);
v.push_back(6794);
int nCase;
scanf("%d", &nCase);
while(nCase--){
scanf("%lld", &n);
printf("%lld\n",1LL * linear_seq::gao(v,n-1) % mod);
}
return 0;
}
二、解特征方程
举例:
斐波那契数列:
其特征方程为:
解得其根为:
则有:
带入
解得
所以:
注:斐波那契数列取模的周期性:
对于模意义下的斐波那契数列,有皮萨诺周期 ( OEIS A001175 )可知,该周期长度不超过,且只有在满足的形式时才取到等号。
高斯消元
const double eps = 1e-9;//精度
const int maxn = 220;//方程个数
double a[maxn][maxn],ans[maxn];
//a[1][1]x1+a[1][2]x2+...+a[1][n]xn=ans[1];
int equ,var;//方程数、未知数个数
int gauss(){
int i,j,k,col,max_r;
for(k=0,col=0;k<equ&&col<var;k++,col++){
max_r=k;
for(i=k+1;i<equ;i++)
if(fabs(a[i][col])>fabs(a[max_r][col]))
max_r=i;
if(fabs(a[max_r][col])<eps)return 0;
if(k!=max_r){
for(j=col;j<var;j++)
swap(a[k][j],a[max_r][j]);
swap(ans[k],ans[max_r]);
}
ans[k]/=a[k][col];
for(j=col+1;j<var;j++)
a[k][j]/=a[k][col];
a[k][col]=1;
for(i=0;i<equ;i++)
if(i!=k){
ans[i]-=ans[k]*a[i][col];
for(j=col+1;j<var;j++)
a[i][j]-=a[k][j]*a[i][col];
a[i][col]=0;
}
}
return 1;
}
线性基
性质:
1、线性基的元素能相互异或得到原集合的元素的所有相互异或得到的值。
2、线性基是满足性质1的最小的集合。
3、线性基没有异或和为0的子集。
4、线性基中每个元素的异或方案唯一,也就是说,线性基中不同的异或组合异或出的数都是不一样的。
5、线性基中每个元素的二进制最高位互不相同。
构造线性基的方法如下:
对原集合的每个数转为二进制,从高位向低位扫,对于第位是1的,如果不存在,那么令并结束扫描,如果存在,令。
inline void insert(long long x) {
for (int i = 55; i + 1; i--) {
if (!(x >> i)) // x的第i位是0
continue;
if (!p[i]) {
p[i] = x;
break;
}
x ^= p[i];
}
}
上一篇: 8、telnet连接服务器
下一篇: SciPy 线性代数