LU分解，Javascript代码

Posted 2020-08-12

///A 为矩阵，这里写成一维数组，如 【1】，【1，2，3，4】
function GetLU(a)
{
	var n = a.length;//矩阵的总数据数目
	var s = Math.sqrt(n);//矩阵的阶数
	var L = new Array(n);
	var U = new Array(n);

	if (GetDet(a) != 0)
{

	var allOrderNotEqulesZero = true;
	for (var i = 0; i < s; i++)
	{
		if (GetDet(GetOrderMatrix(a, i)) == 0)
		{
			allOrderNotEqulesZero  = false;
			break;
		}
	}	
	
        if (allOrderNotEqulesZero)



{
	for (var i = 0; i < s; i++)
	{
		for (var j = 0; j < s; j++)
		{
			if (i == j)
			{
				L[i*s + j] = 1;
			}
			else if (i < j)
			{
			     L[i*s+j] = 0;
			}	
            		if (i > j)
                		U[i*s+j] = 0;

		            U[0*s+j] = a[0*s+j];
            		L[i*s+0] = a[i*s+0] / U[0*s+0];
		}
	}

	 for (var k = 1; k < s; k++)
    	{

        for (var j = k; j < s; j++)
        {
            tmp = 0;
            for (var m = 0; m < k; m++)
            {
                tmp += L[k*s+m] * U[m*s+j];
            }

            U[k*s+j] = a[k*s+j] - tmp;
        }

        for (var i = k+1; i < s; i++)
        {
            tmp = 0;
            for (var m = 0; m < k; m++)
            {
                tmp += L[i*s+m] * U[m*s+k];
            }

            L[i*s+k] = ( a[i*s+k] - tmp ) / U[k*s+k];
        }

    }
}
else 
{
	alert(‘该矩阵的顺序主子式中有零！‘)
	console.log(GetLastMatrix(a));
	for (var i = 0; i < s; i++)
	{
		for (var j = 0; j < s; j++)
		{
			if (j > i)
			{
				U[i * s + j] = a[i * s + j];
				L[i * s + j] = 0;
			}
else if (j == i)
{
	U[i * s + j] = a[i * s + j];
				L[i * s + j] = 1;
}
else
{
	L[i * s + j] = a[i * s + j];
				U[i * s + j] = 0;
}
		}
	}
}

	console.log(L);
	console.log(U);
}
else
{
	alert(‘该矩阵为奇异矩阵，无法进行LU分解！‘);
}
		
}

function GetDet(a)
{
	var n = a.length;//矩阵的总数据数目
	var s = Math.sqrt(n);//矩阵的阶数
	var sum = 0;
	if (n == 1)
	{
		return a[0];
	}

	for (var i = 0; i < s; i++)
	{
		if (i % 2 == 0)
		{
			sum += a[i] * GetDet(GetSubMatrix(a, i))
		}
		else
		{
			sum -= a[i] * GetDet(GetSubMatrix(a, i))
		}
	}

	return sum;
}

///获取比较第n列之后的矩阵
function GetLastMatrix(a)
{
	var n = a.length;//矩阵的总数据数目
	var s = Math.sqrt(n);//矩阵的阶数
	var p = [];
	for (var i = 0; i < s; i++)
	{
		for (var j = 0; j < s; j++)
		{
			if (i == j)
			{
				p.push(1);

			}
			else {
				p.push(0);
			}
		}
	}

	for (var i = 0; i < s -1; i++)
	{
		a = GetNewMatrix(a, i, p); 
	}

	console.log(p);
	return a;
}

function GetNewMatrix(a, k, p)
{
	var n = a.length;//矩阵的总数据数目
	var s = Math.sqrt(n);//矩阵的阶数
	var result = [];
	var maxRowNum = 0;
	var MAXValue = Math.abs(a[ k * s + k]);
 	for (var i = 1; i < s; i++)
	{
		if ( i >= k && Math.abs(a[i * s + k]) > MAXValue )
		{
			maxRowNum = i;
			MAXValue =  Math.abs(a[i * s + k]);
		}
	}

	//console.log(‘maxRowNum:‘ + maxRowNum);
	//console.log(‘MAXValue:‘ + MAXValue);

	if (maxRowNum != k)
	{
		for (var i = 0; i < s; i++)
		{
			
			var temp = a[k * s + i];
			a[k * s + i] = a[maxRowNum * s + i]
			a[maxRowNum * s + i] = temp;

			var tmp1 = p[k * s + i];
			p[k * s + i] = p[maxRowNum * s + i]
			p[maxRowNum * s + i] = tmp1;
		}

		//console.log(a);
		//console.log(p);

		if (MAXValue != 0)
		{
			for (var i = k + 1; i < s; i++)
			{
				
				if (a[i * s + k] != 0)
				{
					a[i * s + k] = parseFloat(a[i * s + k] / a[k * s + k]) 
					for (var j = k + 1; j< s; j++)
					{
						a[i * s + j] = parseFloat(a[i * s + j]- a[i * s + k] * a[k * s + j])
					}		
				}
			}
		}
	}else
	{
		if (MAXValue != 0)
		{
			for (var i = k + 1; i < s; i++)
			{
				
				if (a[i * s + k] != 0)
				{
					a[i * s + k] = parseFloat(a[i * s + k] / a[k * s + k]) 
					for (var j = k + 1; j< s; j++)
					{
						a[i * s + j] = parseFloat(a[k * s + j] - a[i * s + k] * a[i * s + j])
					}		
				}
			}
		}
	}

	//console.log(a);	
	return a;
}

///a为原始矩阵，k为所在的列数
function GetSubMatrix(a, k)
{
	var n = a.length;//矩阵的总数据数目
	var s = Math.sqrt(n);//矩阵的阶数
	var result = []; 
	for (var i = 1; i < s; i++)
	{
		for (var j = 0; j < s; j++)
		{
			if (k != j)
			{
				result.push(a[i*s +j]);
			}	
		}
	}

	return result;
}

///获取顺序主子式的矩阵, k为阶数
function GetOrderMatrix(a, k)
{
	var n = a.length;//矩阵的总数据数目
	var s = Math.sqrt(n);//矩阵的阶数
	var result = []; 
	for (var i = 0; i < s; i++)
	{
		for (var j = 0; j < s; j++)
		{
			if (i <= k && j <= k)
			{
				result.push(a[i*s +j]);
			}	
		}
	}

	return result;	
}