高斯-约旦消元
解线性方程组
例题:线性方程组
步骤:
-
选出未被更新的行中第 \(k\) 列绝对值最大的值,令为主元
-
把主元所在行移到当前行,加减消元消去主元
-
重复12
-
结束后若存在找不到主元(找到是0) 的情况,那就遍历没处理的行,如果有常数项非 \(0\) 则无解,否则无数解
点击查看代码
#include<bits/stdc++.h>
using namespace std;
const double eps=1e-9;
double A[55][55];
int n;
bool eq(double x,double y){
return fabs(x-y)<eps;
}
int main(){
cin>>n;
for(int i=0;i<n;i++)
for(int j=0;j<=n;j++)
cin>>A[i][j];
int nl=0;
for(int k=0;k<n;k++){
int I=nl;
for (int i=nl+1;i<n;i++){
if(fabs(A[i][k])>fabs(A[I][k])) I=i;
}
if(eq(A[I][k],0)) {continue;}
for(int j=0;j<=n;j++) swap(A[nl][j],A[I][j]);
for(int i=0;i<n;i++){
if(i==nl) continue;
double mul=A[i][k]/A[nl][k];
for(int j=k;j<=n;j++){
A[i][j]-=mul*A[nl][j];
}
}
nl++;
}
if(nl<n){
while(nl<n)
if(!eq(A[nl++][n],0))
{cout<<-1;return 0;}// 无解-1
cout<<0;
}
else {
for(int i=0;i<n;i++){
printf("x%d=%.2lf\n",i+1,A[i][n]/A[i][i]);
}
}
}
计算矩阵逆
把一个单位矩阵拼在给出矩阵后面,高斯约旦消元后再额外将结果的对角线化1,则此时原先单位矩阵位置就是逆矩阵
点击查看代码
#include<bits/stdc++.h>
using namespace std;
#define int long long
const int mod=1e9+7;
const double eps=1e-9;
int A[405][805];
int n;
bool eq(int x,int y){
return abs(x-y)<eps;
}
int qpow(int x,int k){
int ans=1;
while(k){
if(k&1) ans=ans*x%mod;
x=x*x%mod;
k>>=1;
}
return ans%mod;
}
signed main(){
cin>>n;
for(int i=0;i<n;i++)
for(int j=0;j<n;j++)
cin>>A[i][j];
for(int i=0;i<n;i++){
A[i][i+n]=1;
}
int nl=0;
for(int k=0;k<n;k++){
int I=nl;
for (int i=nl+1;i<n;i++){
if(abs(A[i][k])>abs(A[I][k])) I=i;
}
if(eq(A[I][k],0)) {cout<<"No Solution\n";return 0;}
swap(A[nl],A[I]);
int kk=qpow(A[nl][k],mod-2);
for(int i=0;i<n;i++){
if(i==nl) continue;
int mul=A[i][k]*kk%mod;
for(int j=k;j<n*2;j++){
A[i][j]=(mod+A[i][j]-mul*A[nl][j]%mod)%mod;
}
}
for(int j=0;j<2*n;j++){
A[nl][j]=(A[nl][j]*kk)%mod;
}
nl++;
}
for(int i=0;i<n;i++){
for(int j=n;j<n*2;j++)cout<<(A[i][j]+mod)%mod<<" ";
cout<<"\n";
}
}
计算行列式
这玩意真的还是高斯消元吗...
点击查看代码
#include<bits/stdc++.h>
#define int long long
using namespace std;
const int N=605;
int mod,OK=1;;
int A[N][N];
int n;
int det_calc(int a[N][N],int n){
for(int i=0;i<n;i++){
for(int j=i+1;j<n;j++){
while(a[i][i]){
int K=a[j][i]/a[i][i];
for(int k=i;k<n;k++){
((a[j][k]-=K*a[i][k]%mod)+=mod)%=mod;
}
swap(a[i],a[j]);OK*=-1;
}swap(a[i] ,a[j]),OK*=-1;
}
}
int ans=1;
for(int i=0;i<n;i++) ans*=a[i][i],ans%=mod;
ans=ans*OK+mod;ans%=mod;
return ans ;
}
signed main(){
cin>>n>>mod;
for(int i=0;i<n;i++)
for(int j=0;j<n;j++)
cin>>A[i][j];
cout<<det_calc(A,n);
}