gauss_jordan.cpp

来自「数值计算的例子和参考书籍」· C++ 代码 · 共 40 行

bool gauss_jordan(void)
{
	int i, j, k, ik;
	for(k = 1; k <= n; k++)
	{
		double mkk, max = -inf;
		ik = -1;
		for(i = k; i <= n; i++)
		{
			if(a[i][k] > max)
			{
				max = a[i][k];
				ik = i;
			}
		}
		if(fabs(max) < eps) return false;
		if(ik != k)
		{
			for(i = 1; i <= n; i++) swap(a[k][i], a[ik][i]);
			swap(b[k], b[ik]);
		}
		for(i = 1; i <= n; i++)
		{
			if(i != k) a[i][k] /= -a[k][k];
			else mkk = 1 / a[k][k];
		}
		for(i = 1; i <= n; i++)
		{
			if(i != k)
			{
				for(j = k + 1; j <= n; j++) a[i][j] = a[i][j] + a[i][k] * a[k][j];
				b[i] = b[i] + a[i][k] * b[k];
			}
		}
		for(j = k; j <= n; j++) a[k][j] = a[k][j] * mkk;
		b[k] = b[k] * mkk;
	}
	return true;
}