[其他] 矩阵表达式计算器

slimay

神奇的矩阵表达式运算工具
下载   http://cmd1152.ys168.com/  文件区  矩阵表达式计算器.zip
MAT.EXE
摘要：
=============================================================================
矩阵表达式计算器，支持矩阵四则运算,乘方运算计算。
=============================================================================

用法：
-----------------------------------------------------------------------------
mat "[expression1] [expression2] [expression3] [.] ..."
-----------------------------------------------------------------------------

示例:
-----------------------------------------------------------------------------
mat "[A=1,0;1,1] [B=A^(2)] [C=A^(-1)+B] [(A^(-1))-C] [det(A)]"
-----------------------------------------------------------------------------

输入格式：
    {E = en(2)}      <==>   [1   0]
                            [0   1]
    {A = 1,,2;3,8,7} <==>   [1   0   2]
                            [3   8   7]

~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
矩阵函数:
    +, -, *, ^, A^(-1), A^(n), cp(A), rot(A), en(A)
    gs(A), diag(A), lad(A), tri(A)
    tr(A), r(A), det(A)
    lf(n): 控制精度
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
常数类
        pi    3.1415926535897932
        e     2.7182818284590452       

通用类
        +     加
        -     减
        *     乘
        /     除
        %     取余数
        ^     次方
        rand  随机数
        round 四舍五入
        int   取整
        ceil  向上舍入
        floor 向下舍入
        sqrt  开方
        lg    常用对数，以10为底
        ln    自然对数
        exp   e的次幂
        deg   度转弧度

三角函数类
        sin、cos、tan   
        arcsin、arccos、arctan

双曲函数类
        sinh、cosh、tanh
        arcsinh、arccosh、arctanh

~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

英译：
The matrix calculator, Copyright@2021~2023 Slimay
usage:
    mat "[expression1] [expression2] [expression3] [.] ..."
example:
    mat "[A=1,0;1,1] [B=A^(2)] [C=A^(-1)+B] [(A^(-1))-C] [det(A
    You can use the letters A to O to name the matrix
    {E = en(2)}      <==>   [1   0]
                            [0   1]
    {A = 1,,2;3,8,7} <==>   [1   0   2]
                            [3   8   7]
matrix functions:
    +, -, *, ^, A^(-1), A^(n), cp(A), rot(A), en(A)
    gs(A), diag(A), lad(A), tri(A)
    tr(A), r(A), det(A)
    lf(n): control print precision

mathematical function:
    pi = 3.1415, e = 2.7182, rand, abs, sqrt, lg, ln, exp
    sin, cos, tan, arcsin, arccos, arctan
version 1.0.5
COPYRIGHT@2021~2023 BY SLIMAY, VERSION 1.0
MAT.EXE

核心代码

//高斯若当主函数 
Matrix* Gauss_Jordan_Elimination(Matrix* mat, int eliminateMode)
{
	int eliminateModeBak =  eliminateMode;
	
	if(eliminateMode == ELIMINATE_INVERSE2)
	{
		eliminateMode = ELIMINATE_INVERSE;
	}
 
	if(eliminateMode == ELIMINATE_GET_DET)
	{
		eliminateMode = ELIMINATE_TRIANGLE;
	}
	
	//获取第一个矩阵列数 
	int n1 = ( *(mat + 0) );
	//获取第一个矩阵行数 
	int m1 = ( *(mat + 1) );
 
 
	//低于二阶的矩阵没法阶梯化 
	if(n1 <2 || m1 < 2)
	{
		PRINTF_ERROR("The 1-order matrix can't run this.\n", "");
		exit(1);
		return real_matrix(1.0f / (*(mat + 2)));
	}
 
	//用传入矩阵的备份做 高斯若当处理 
	Matrix* mat1 = (Matrix *)calloc(n1 * m1 + MAT_PRELEN, sizeof(double));
	memcpy(mat1, mat, sizeof(double) * (n1 * m1 + MAT_PRELEN) ); 
 
	//行交换存储器 
	double exChangeLineTmp[MAX_LEN];
	
	Matrix* matInverse = NULL; 
	
	if(eliminateMode == ELIMINATE_INVERSE)
	{
		//动态分配一个存放逆矩阵的空间 
		matInverse = (Matrix *)calloc(n1 * m1 + MAT_PRELEN, sizeof(double));
		//写入列行数 
		*(matInverse + 0) =  (double)n1;
		*(matInverse + 1) =  (double)m1;
		// 先置一同阶单位矩阵 
		for(int i = 0; i < GET_MIN(n1, m1); i++) 
		{
			*(matInverse + MAT_PRELEN + i + i * n1) = 1.0f;
		}
	}
	
	//行交换的次数 
	int exChangeTimes = 0;
	//高斯若当消去法置换成上三角型矩阵 
	int i = 0, j = 0, j_start = 0;
	for(i = 0; i < GET_MIN(n1, m1); i++) 
	{	
		for(j = j_start; (j < m1) && (i < n1); j++) 
		{
			//找到主元直接中断内循环 
			if( (*(mat1 + MAT_PRELEN  + i + j_start * n1)) != 0)
			{
				break;
			}
			//找到主元为0，则通过行交换使得主元不为0 
			if( (*(mat1 + MAT_PRELEN  + i + j*n1) != 0) && (j != j_start))
			{
				
				memcpy(exChangeLineTmp, mat1 + MAT_PRELEN + j*n1, n1 * sizeof(double));			
				// i行与j行互换 
				memcpy(mat1 + MAT_PRELEN + j * n1,    mat1 + MAT_PRELEN + i * n1,            n1 * sizeof(double));
				memcpy(mat1 + MAT_PRELEN + i * n1,    exChangeLineTmp,                       n1 * sizeof(double));
				
				// 如果是求逆矩阵模式 
				if(eliminateMode == ELIMINATE_INVERSE)
				{ 
					memcpy(exChangeLineTmp, matInverse + MAT_PRELEN + j*n1, n1 * sizeof(double));			
					// i行与j行互换 
					memcpy(matInverse + MAT_PRELEN + j * n1,    matInverse + MAT_PRELEN + i * n1,            n1 * sizeof(double));
					memcpy(matInverse + MAT_PRELEN + i * n1,    exChangeLineTmp,                             n1 * sizeof(double));
				}
						
				exChangeTimes ++;
				break;
			}
		}
 
		//如果该列主元都为0，直接跳过，处理下一列 
		if((j == m1) && (*(mat1 + MAT_PRELEN  + i + j_start * n1) == 0))
		{
			// 这里根据开关选择是否采用阶梯消去法，阶梯消去用于非方阵的m*n型矩阵秩计算 
			if((eliminateModeBak != ELIMINATE_LADDER))
			{
				j_start ++;
			}
			continue;
		}
 
		//如果该列主元不为0，则本行除以主元 
		if( (*(mat1 + MAT_PRELEN  + i + j_start * n1)) != 0)
		{
			// 这里根据开关选择是否采用主元归一，
			if((eliminateMode == ELIMINATE_INVERSE))
			{
				// 主元归一, 避免误差 
				double mainPara = (*(mat1 + MAT_PRELEN  + i + j_start*n1));
				for(int k = i; k < n1; k++) 
				{
					(*(mat1 + MAT_PRELEN  + k + j_start*n1)) /= mainPara;
					
					// 精度处理
					_ZERO_( *(mat1 + MAT_PRELEN  + k + j_start*n1) );
				}
				
				if(eliminateMode == ELIMINATE_INVERSE)
				{
					for(int k = 0; k < n1; k++) 
					{
						(*(matInverse + MAT_PRELEN  + k + j_start*n1)) /= mainPara;
 
						// 精度处理
						_ZERO_( *(matInverse + MAT_PRELEN  + k + j_start*n1) );	
					}
				}
			}
		}
 
		//开始消元 
		for(j = ( (eliminateMode==ELIMINATE_DIAG) || (eliminateMode==ELIMINATE_INVERSE) )?0 :(j_start + 1); (j < m1) && (i < n1); j++) 
		{
			// 主元不消去自己 
			if(j == j_start)
			{
				continue;
			}
			
			//计算消去系数
			double disLinePara =  - (*(mat1 + MAT_PRELEN  + i + j*n1)) / (*(mat1 + MAT_PRELEN  + i + j_start*n1));
				
			//用主元行消去其他行 
			for(int k = i; k < n1; k++) 
			{
				(*(mat1 + MAT_PRELEN  + k + j*n1))  +=  disLinePara * (*(mat1 + MAT_PRELEN  + k + j_start*n1));
					
				// 精度处理
				_ZERO_(*(mat1 + MAT_PRELEN  + k + j*n1));
			}
			
			// 如果是求逆矩阵模式 
			if(eliminateMode == ELIMINATE_INVERSE)
			{
				//计算消去系数 (这里沿用左边的消去系数，因为 AE -> EA^(-1),故两边乘以的消去系数一致
				double disLineParaInverse = disLinePara ;
					
				//用主元行消去其他行 
				for(int k = ((eliminateMode==ELIMINATE_DIAG) || (eliminateMode==ELIMINATE_INVERSE) ?0 :i); k < n1; k++) 
				{
					(*(matInverse + MAT_PRELEN  + k + j*n1))  += disLineParaInverse * (*(matInverse + MAT_PRELEN  + k + j_start*n1));
					// 精度控制 
					_ZERO_( *(matInverse + MAT_PRELEN  + k + j*n1) );
				}
			} 
		}
 
		//下沉一行
		j_start ++;	
	}
	
	// 如果是det求值模式 
	if(eliminateModeBak == ELIMINATE_GET_DET)
	{	
		double detValue = 1.0f;	
		for(int j = 0; j < m1; j++) 
		{
			 detValue  *=  (*(mat1 + MAT_PRELEN  + j + j*m1));
		}
		
		// 精度控制 
		_ZERO_(detValue);
				
		//释放备份矩阵内存 
		free(mat1);
			
		//返回det值 1阶矩阵 
		return  real_matrix(detValue);
	}
	
	//如果是求逆矩阵模式 
	if(eliminateModeBak == ELIMINATE_INVERSE)
	{
		//获取首非零行个数 
		for(int j = 0; j < n1; j++) 
		{
			double mPara = (*(mat1 + MAT_PRELEN  + j + j*n1));
			if( mPara == 0.0f )
			{
				PRINTF_ERROR("The matrix has no inverse mat.\n", "");
				exit(1);
			}
		}
		
		//释放备份矩阵内存 
		free(mat1);	
		return  matInverse;
	}
 
	//如果是解方程组模式 
	if(eliminateModeBak == ELIMINATE_INVERSE2)
	{
		free(matInverse);
		return mat1;
	}
 
	//返回行交换次数 
	return mat1;
} COPY