{ 
				pDataConst[j]=pDataConst[j]/pDataCoef[0];
				for (i=1; i<=n-1; i++)
				{ 
					u=i*n+i; 
					v=i*m+j;
					for (k=1; k<=i; k++)
						pDataConst[v]=pDataConst[v]-pDataCoef[(k-1)*n+i]*pDataConst[(k-1)*m+j];
					pDataConst[v]=pDataConst[v]/pDataCoef[u];
				}
			}
	    
			for (j=0; j<=m-1; j++)
			{ 
				u=(n-1)*m+j;
				pDataConst[u]=pDataConst[u]/pDataCoef[n*n-1];
				for (k=n-1; k>=1; k--)
				{ 
					u=(k-1)*m+j;
					for (i=k; i<=n-1; i++)
					{ 
						v=(k-1)*n+i;
						pDataConst[u]=pDataConst[u]-pDataCoef[v]*pDataConst[i*m+j];
					}
	            
					v=(k-1)*n+k-1;
					pDataConst[u]=pDataConst[u]/pDataCoef[v];
				}
			}
	    
			return true;
		}

		/**
		 * 求解大型稀疏方程組的全選主元高斯－約去消去法
		 * 
		 * @param mtxResult - CMatrix引用對象，返回方程組解矩陣
		 * @return bool 型，方程組求解是否成功
		 */
		public bool GetRootsetGgje(Matrix mtxResult)
		{ 
			int i,j,k,nIs=0,u,v;
			double d,t;
	    
			// 方程組屬性，將常數矩陣賦給解矩陣
			Matrix mtxCoef = new Matrix(mtxLECoef);
			mtxResult.SetValue(mtxLEConst);
			int n = mtxCoef.GetNumColumns();
			double[] pDataCoef = mtxCoef.GetData();
			double[] pDataConst = mtxResult.GetData();

			// 臨時緩沖區，存放變換的列數
			int[] pnJs = new int[n];

			// 消元
			for (k=0; k<=n-1; k++)
			{ 
				d=0.0;
				for (i=k; i<=n-1; i++)
				{
					for (j=k; j<=n-1; j++)
					{ 
						t=Math.Abs(pDataCoef[i*n+j]);
						if (t>d) 
						{
							d=t; 
							pnJs[k]=j; 
							nIs=i;
						}
					}
				}

				if (d == 0.0)
				{
					return false;
				}

				if (nIs!=k)
				{ 
					for (j=k; j<=n-1; j++)
					{ 
						u=k*n+j; 
						v=nIs*n+j;
						t=pDataCoef[u]; 
						pDataCoef[u]=pDataCoef[v]; 
						pDataCoef[v]=t;
					}
	            
					t=pDataConst[k]; 
					pDataConst[k]=pDataConst[nIs]; 
					pDataConst[nIs]=t;
				}
	        
				if (pnJs[k]!=k)
				{
					for (i=0; i<=n-1; i++)
					{ 
						u=i*n+k; 
						v=i*n+pnJs[k];
						t=pDataCoef[u]; 
						pDataCoef[u]=pDataCoef[v]; 
						pDataCoef[v]=t;
					}
				}

				t=pDataCoef[k*n+k];
				for (j=k+1; j<=n-1; j++)
				{ 
					u=k*n+j;
					if (pDataCoef[u]!=0.0) 
						pDataCoef[u]=pDataCoef[u]/t;
				}
	        
				pDataConst[k]=pDataConst[k]/t;
				for (j=k+1; j<=n-1; j++)
				{ 
					u=k*n+j;
					if (pDataCoef[u]!=0.0)
					{ 
						for (i=0; i<=n-1; i++)
						{ 
							v=i*n+k;
							if ((i!=k)&&(pDataCoef[v]!=0.0))
							{ 
								nIs=i*n+j;
								pDataCoef[nIs]=pDataCoef[nIs]-pDataCoef[v]*pDataCoef[u];
							}
						}
					}
				}
	        
				for (i=0; i<=n-1; i++)
				{ 
					u=i*n+k;
					if ((i!=k)&&(pDataCoef[u]!=0.0))
						pDataConst[i]=pDataConst[i]-pDataCoef[u]*pDataConst[k];
				}
			}
	    
			// 調整
			for (k=n-1; k>=0; k--)
			{
				if (k!=pnJs[k])
				{ 
					t=pDataConst[k]; 
					pDataConst[k]=pDataConst[pnJs[k]]; 
					pDataConst[pnJs[k]]=t;
				}
			}

			return true;
		}

		/**
		 * 求解托伯利茲方程組的列文遜方法
		 * 
		 * @param mtxResult - CMatrix引用對象，返回方程組解矩陣
		 * @return bool 型，方程組求解是否成功
		 */
		public bool GetRootsetTlvs(Matrix mtxResult)
		{ 
			int i,j,k;
			double a,beta,q,c,h;

			// 未知數個數
			int n = mtxLECoef.GetNumColumns();

			// 初始化解解向量
			mtxResult.Init(n, 1);
			double[] x = mtxResult.GetData();

			// 常數數組
			double[] pDataConst = mtxLEConst.GetData();

			// 建立T數組
			double[] t = new double[n];

			// 構造T數組
			for (i=0; i<n; ++i)
				t[i] = mtxLECoef.GetElement(0, i);

			// 臨時數組
			double[] s = new double[n];
			double[] y = new double[n];

			// 非托伯利茲方程組，不能用本方法求解
			a=t[0];
			if (a == 0.0)
			{ 
				return false;
			}

			// 列文遜方法求解
			y[0]=1.0; 
			x[0]=pDataConst[0]/a;
			for (k=1; k<=n-1; k++)
			{ 
				beta=0.0; 
				q=0.0;
				for (j=0; j<=k-1; j++)
				{ 
					beta=beta+y[j]*t[j+1];
					q=q+x[j]*t[k-j];
				}
	        
				if (a == 0.0)
				{ 
					return false;
				}

				c=-beta/a; 
				s[0]=c*y[k-1]; 
				y[k]=y[k-1];
				if (k!=1)
				{
					for (i=1; i<=k-1; i++)
						s[i]=y[i-1]+c*y[k-i-1];
				}

				a=a+c*beta;
				if (a == 0.0)
				{ 
					return false;
				}

				h=(pDataConst[k]-q)/a;
				for (i=0; i<=k-1; i++)
				{ 
					x[i]=x[i]+h*s[i]; 
					y[i]=s[i];
				}
	        
				x[k]=h*y[k];
			}
	    
			return true;
		}

		/**
		 * 高斯－賽德爾迭代法
		 *  
		 * @param mtxResult - CMatrix引用對象，返回方程組解矩陣
		 * @param eps - 控制精度
		 * @return bool 型，方程組求解是否成功
		 */
		public bool GetRootsetGaussSeidel(Matrix mtxResult, double eps)
		{ 
			int i,j,u,v;
			double p,t,s,q;

			// 未知數個數
			int n = mtxLECoef.GetNumColumns();

			// 初始化解向量
			mtxResult.Init(n, 1);
			double[] x = mtxResult.GetData();

			// 系數與常數
			double[] pDataCoef = mtxLECoef.GetData();
			double[] pDataConst = mtxLEConst.GetData();
	    
			// 求解
			for (i=0; i<=n-1; i++)
			{ 
				u=i*n+i; 
				p=0.0; 
				x[i]=0.0;
				for (j=0; j<=n-1; j++)
				{
					if (i!=j)
					{ 
						v=i*n+j; 
						p=p+Math.Abs(pDataCoef[v]);
					}
				}

				if (p>=Math.Abs(pDataCoef[u]))
					return false;
			}

			// 精度控制
			p=eps+1.0;
			while (p>=eps)
			{ 
				p=0.0;
				for (i=0; i<=n-1; i++)
				{ 
					t=x[i]; 
					s=0.0;
					for (j=0; j<=n-1; j++)
						if (j!=i) 
							s=s+pDataCoef[i*n+j]*x[j];

					x[i]=(pDataConst[i]-s)/pDataCoef[i*n+i];
					q=Math.Abs(x[i]-t)/(1.0+Math.Abs(x[i]));
					if (q>p) 
						p=q;
				}
			}
	    
			return true;
		}

		/**
		 * 求解對稱正定方程組的共軛梯度法
		 * 
		 * @param mtxResult - CMatrix引用對象，返回方程組解矩陣
		 * @param eps - 控制精度
		 */
		public void GetRootsetGrad(Matrix mtxResult, double eps)
		{ 
			int i,k;
			double alpha,beta,d,e;

			// 未知數個數
			int n = GetNumUnknowns();

			// 初始化解向量
			mtxResult.Init(n, 1);
			double[] x = mtxResult.GetData();

			// 構造臨時矩陣
			Matrix mtxP = new Matrix(n, 1);
			double[] p = mtxP.GetData();

			double[] pDataCoef = mtxLECoef.GetData();
			double[] pDataConst = mtxLEConst.GetData();

			double[] r = new double[n];

			for (i=0; i<=n-1; i++)
			{ 
				x[i]=0.0; 
				p[i]=pDataConst[i]; 
				r[i]=pDataConst[i]; 
			}
	    
			i=0;
			while (i<=n-1)
			{ 
				Matrix mtxS = mtxLECoef.Multiply(mtxP);
				double[] s = mtxS.GetData();
	        
				d=0.0; 
				e=0.0;
				for (k=0; k<=n-1; k++)
				{ 
					d=d+p[k]*pDataConst[k]; 
					e=e+p[k]*s[k]; 
				}
	        
				alpha=d/e;
				for (k=0; k<=n-1; k++)
					x[k]=x[k]+alpha*p[k];
	        
				Matrix mtxQ = mtxLECoef.Multiply(mtxResult);
				double[] q = mtxQ.GetData();
	        
				d=0.0;
				for (k=0; k<=n-1; k++)
				{ 
					r[k]=pDataConst[k]-q[k]; 
					d=d+r[k]*s[k]; 
				}
	        
				beta=d/e; d=0.0;
				for (k=0; k<=n-1; k++) 
					d=d+r[k]*r[k];
	        
				// 滿足精度，求解結束
				d=Math.Sqrt(d);
				if (d<eps)
					break;

				for (k=0; k<=n-1; k++)
					p[k]=r[k]-beta*p[k];
	        
				i=i+1;
			}
		}

		/**
		 * 求解線性最小二乘問題的豪斯荷爾德變換法
		 * 
		 * @param mtxResult - Matrix對象，返回方程組解矩陣
		 * @param mtxQ - Matrix對象，返回豪斯荷爾德變換的Q矩陣
		 * @param mtxR - Matrix對象，返回豪斯荷爾德變換的R矩陣
		 * @return bool 型，方程組求解是否成功
		 */
		public bool GetRootsetMqr(Matrix mtxResult, Matrix mtxQ, Matrix mtxR)
		{ 
			int i,j;
			double d;

			// 方程組的方程數和未知數個數
			int m = mtxLECoef.GetNumRows();
			int n = mtxLECoef.GetNumColumns();
			// 奇異方程組
			if (m < n)
				return false;

			// 將解向量初始化為常數向量
			mtxResult.SetValue(mtxLEConst);
			double[] pDataConst = mtxResult.GetData();

			// 構造臨時矩陣，用于QR分解
			mtxR.SetValue(mtxLECoef);
			double[] pDataCoef = mtxR.GetData();

			// QR分解
			if (! mtxR.SplitQR(mtxQ))
				return false;

			// 臨時緩沖區
			double[] c = new double[n];
			double[] q = mtxQ.GetData();

			// 求解
			for (i=0; i<=n-1; i++)
			{ 
				d=0.0;
				for (j=0; j<=m-1; j++)
					d=d+q[j*m+i]*pDataConst[j];
	    
				c[i]=d;
			}
	    
			pDataConst[n-1]=c[n-1]/pDataCoef[n*n-1];
			for (i=n-2; i>=0; i--)
			{ 
				d=0.0;
				for (j=i+1; j<=n-1; j++)
					d=d+pDataCoef[i*n+j]*pDataConst[j];
	        
				pDataConst[i]=(c[i]-d)/pDataCoef[i*n+i];
			}
	    
			return true;
		}

		/**
		 * 求解線性最小二乘問題的廣義逆法
		 * 
		 * @param mtxResult - Matrix對象，返回方程組解矩陣
		 * @param mtxAP - Matrix對象，返回系數矩陣的廣義逆矩陣
		 * @param mtxU - Matrix對象，返回U矩陣
		 * @param mtxV - Matrix對象，返回V矩陣
		 * @param eps - 控制精度
		 * @return bool 型，方程組求解是否成功
		 */
		public bool GetRootsetGinv(Matrix mtxResult, Matrix mtxAP, Matrix mtxU, Matrix mtxV, double eps)
		{ 
			int i,j;
	    
			// 方程個數和未知數個數
			int m = mtxLECoef.GetNumRows();
			int n = mtxLECoef.GetNumColumns();

			// 初始化解向量
			mtxResult.Init(n, 1);

			double[] pDataConst = mtxLEConst.GetData();
			double[] x = mtxResult.GetData();

			// 臨時矩陣
			Matrix mtxA = new Matrix(mtxLECoef);

			// 求廣義逆矩陣
			if (! mtxA.InvertUV(mtxAP, mtxU, mtxV, eps))
				return false;

			double[] pAPData = mtxAP.GetData();

			// 求解
			for (i=0; i<=n-1; i++)
			{ 
				x[i]=0.0;
				for (j=0; j<=m-1; j++)
					x[i]=x[i]+pAPData[i*m+j]*pDataConst[j];
			}
	    
			return true;
		}

		/**
		 * 
		 * @param mtxResult - Matrix對象，返回方程組解矩陣
		 * @param nMaxIt - 疊加次數
		 * @param eps - 控制精度
		 * @return bool 型，方程組求解是否成功
		 */
		public bool GetRootsetMorbid(Matrix mtxResult, int nMaxIt /*= 60*/, double eps)
		{ 
			int i, k;
			double q, qq;
	    
			// 方程的階數
			int n = GetNumUnknowns();

			// 設定迭代次數, 缺省為60
			i = nMaxIt;
	    
			// 用全選主元高斯消元法求解
			LEquations leqs = new LEquations(mtxLECoef, mtxLEConst);
			if (! leqs.GetRootsetGauss(mtxResult))
				return false;
			double[] x = mtxResult.GetData();
	    
			q=1.0+eps;
			while (q>=eps)
			{ 
				// 迭代次數已達最大值，仍為求得結果，求解失敗
				if (i==0)
					return false;
	        
				// 迭代次數減1
				i=i-1;
			
				// 矩陣運算
				Matrix mtxE = mtxLECoef.Multiply(mtxResult);
				Matrix mtxR = mtxLEConst.Subtract(mtxE);

				// 用全選主元高斯消元法求解
				leqs = new LEquations(mtxLECoef, mtxR);
				Matrix mtxRR = new Matrix();
				if (! leqs.GetRootsetGauss(mtxRR))
					return false;

				double[] r = mtxRR.GetData();
	        
				q=0.0;
				for ( k=0; k<=n-1; k++)
				{ 
					qq=Math.Abs(r[k])/(1.0+Math.Abs(x[k]+r[k]));
					if (qq>q) 
						q=qq;
				}
	        
				for ( k=0; k<=n-1; k++) 
					x[k]=x[k]+r[k];

			}
	    
			// 求解成功
			return true;
		}
	}
}