高斯消元模板！！！bzoj1013

/* 高斯消元模板题 n维球体确定圆心必须要用到n+1个点 设圆心坐标(x1,x2,x3,x4...xn),半径为C 设第i个点坐标为(ai1,ai2,ai3,,,ain)那么对应的方程为 (x1-ai1)^2+(x2-ai2)^2+...+(xn-ain)^2=C*C 如此可列出n+1个方程但是由于有 xi^2 在，无法高斯消元 所以将这n+1个方程上下相减，得 2(x[1]*a[i][1]-x[1]a[i+1][1])+(a[i][1]^2-a[i+1][1]^2)...=0 那么化简后就是 sum{2*x[j]*(a[i][j]-a[i+1][j])}=sum{a[i][j]^2-a[i+1][j]^2} 那么可以用高斯消元做了！ */ #include<bits/stdc++.h> using namespace std; double A[20][20],C[20][20],B[20];//系数矩阵，常数矩阵 int n; int main(){ cin>>n; for(int i=1;i<=n+1;i++) for(int j=1;j<=n;j++) cin>>A[i][j];//输入n+1个点的坐标 for(int i=1;i<=n;i++)//处理出增广矩阵 for(int j=1;j<=n;j++){ C[i][j]=2*(A[i][j] - A[i+1][j]); B[i]+=(A[i][j]*A[i][j] - A[i+1][j]*A[i+1][j]); } //高斯消元！ for(int i=1;i<=n;i++){ //找到xi系数不为0的第一个方程，并将其移到第i个方程处 for(int j=i;j<=n;j++){ if(C[j][i]>1e-8){ for(int k=1;k<=n;k++) swap(C[i][k],C[j][k]); swap(B[i],B[j]); } } //用xi的系数去消其余方程的系数 for(int j=1;j<=n;j++){ if(j==i)continue; double r=C[j][i]/C[i][i]; for(int k=1;k<=n;k++) C[j][k]-=r*C[i][k]; B[j]-=r*B[i]; } } for(int i=1;i<=n;i++) printf("%.3lf ",B[i]/C[i][i]); }

那么下面就是高斯消元的模板，其中C是系数矩阵，B是常数矩阵

//高斯消元！ for(int i=1;i<=n;i++){ //找到xi系数不为0的第一个方程，并将其移到第i个方程处 for(int j=i;j<=n;j++){ if(C[j][i]>1e-8){ for(int k=1;k<=n;k++) swap(C[i][k],C[j][k]); swap(B[i],B[j]); } } //用xi的系数去消其余方程的系数 for(int j=1;j<=n;j++){ if(j==i)continue; double r=C[j][i]/C[i][i]; for(int k=1;k<=n;k++) C[j][k]-=r*C[i][k]; B[j]-=r*B[i]; } }

 标准板子

void guess(int equ,int var){ //行，列 int k=1,col=1,max_r;//行 列 最大列 for(k=1,col=1;k<=equ&&col<var;k++,col++){ max_r=k; for(int i=k+1;i<=equ;i++){//寻找当前最大 if(fabs(mat[i][col])-fabs(mat[max_r][col])>eps) max_r=i; } if(max_r!=k){//如果不是，就换过来 for(int j=1;j<=var;j++) swap(mat[max_r][j],mat[k][j]); } if(fabs(mat[k][col])<eps){//如果已经为0，行不变，移到下一列 k--; continue; } for(int i=k+1;i<=equ;i++){ if(fabs(mat[i][col])>eps){ double t=mat[i][col]/mat[k][col]; mat[i][col]=0.0; for(int j=col+1;j<=var;j++) mat[i][j]-=mat[k][j]*t; } } } for(int i=equ;i>=1;i--){ if(fabs(mat[i][i])<eps) continue; double tmp=mat[i][var]; for(int j=i+1;j<var;j++) if(mat[i][j]!=0) tmp-=mat[i][j]*x[j]; x[i]=tmp/mat[i][i]; } }

 

转载于:https://www.cnblogs.com/zsben991126/p/10529098.html