【高斯消元-浮点数模板】高斯消元版本1(浮点数)

Hello!这里是“高斯消元”...

个人感觉就是高级一点的模拟...

先存代码,有时间再把内容补充上来

Ps.如果是要求整数解,我们就可以把1e-6换成0,在执行高斯消元的时候求最小公倍数

(发现好多人把"gauss"打成"guass",为此特地上网搜了一下233)

#include<cstdio>
#include<cstring>
#include<cmath>
#include<iostream>
#include<algorithm>
using namespace std;
const int MAXN=400;
double a[MAXN+5][MAXN+5],x[MAXN+5];
bool vis[MAXN+5];
int n,m;
void gauss()
{
	//一、正常求解部分——
	int now=1,k=1;//k枚举的是第几个方程,now是在求第now个未知数
	for(k=1;k<=n&&now<=m;k++,now++)
	{
		int j=k;
		for(int i=k+1;i<=n;i++)
			if(fabs(a[j][now])<fabs(a[i][now]))
				j=i;//找绝对值最大的系数
		if(fabs(a[j][now])<1e-6)//意思是"==0",因为是浮点数,所以有一点精度问题
		{//如果最大的都是0,则方程组里面第now个未知数被削完,则continue
			k--;//这里k--再continue后进入下一个循环,k++,即k不变
			continue;
		}
		if(j!=k)
			for(int i=1;i<=m+1;i++)//m+1:最后一列保存的是方程等号后面的的常数,也要交换
				swap(a[j][i],a[k][i]);
		for(int i=1;i<=n;i++)
		{
			//注意其实这里现在的第i行是不会变化的,所以我们第一步
			//找的时候是不能要前k个方程组
			if(fabs(a[i][now])>1e-6&&i!=k)
			{
				double t=a[i][now]/a[k][now];//乘的倍数,要让其它方程的第now个未知数为0
				for(int j=now;j<=m+1;j++)
					a[i][j]-=t*a[k][j];//乘了之后再用加减消元法
			}
		}
	}
	//二、处理无解的情况——eg."0·x=6"
	for(int i=k;i<=n;i++)
		if(fabs(a[i][m+1])>1e-6)//意思是"!=0"
		{
			printf("No solution");
			exit(0);
		}
	//三、处理自由元的情况——eg."0·x+0·y=0"
	//x[]数组是解 ,vis[i]表示第i个未知数是否有解
	if(k<=m)//只有k-1个方程组遍历过
	{
		for(int i=k-1;i>=1;i--)
		{
			int sum=0,id=0;//sum表示这个方程有多少未知数,如果sum=1,则id存的是未知数下标
			for(int j=1;j<=m;j++)
				if(fabs(a[i][j])>1e-6&&!vis[j])
					sum++,id=j;//不为0,且是未知数
			if(sum>1)//如果大于1,那么这个方程不能出解,则continue
				continue;
			double temp=a[i][m+1];
				for(int j=1;j<=m;j++)
					if(j!=id&&fabs(a[i][j])>1e-6)
						temp-=a[i][j]*x[j];//把已知元消掉
			x[id]=temp/a[i][id];//这个未知数是有解的,算出来
			vis[id]=1;//表示有解
		}
		return ;
	}
	//四、找普通情况的解——
	for(int i=m;i>=1;i--)
	{
		double temp=a[i][m+1];
		for(int j=i+1;j<=m;j++)
			temp-=a[i][j]*x[j];//消元
		x[i]=temp/a[i][i];
		vis[i]=1;
	}
}
int main()
{
	scanf("%d%d",&n,&m);
	for(int i=1;i<=n;i++)
		for(int j=1;j<=m+1;j++)
			scanf("%lf",&a[i][j]);
	gauss();
	for(int i=1;i<=m;i++)
	{
		if(vis[i]==0)
			printf("X[%d] not determined\n" , i);
		else
		{
			if(fabs(x[i])<1e-6)
				x[i]=0;
			printf("X[%d] = %.4lf\n",i,x[i]);
		}
	}
	return 0;
}

 

  • 0
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值