Eigen使用demo

Eigen是一个高层次的C ++库，有效支持线性代数，矩阵和矢量运算，数值分析及其相关的算法

以下为Eigen库的具体实现功能说明

// Matrix-test.cpp : 此文件包含 "main" 函数。程序执行将在此处开始并结束。
//

#include "pch.h"
#include <iostream>
#include <Eigen/Dense>
#include <Eigen/Sparse>
#include<fstream>
#include <sstream>
#include<vector>
#include <string>


//这是基于eigen的矩阵使用demo，Eigen是一个高层次的C ++库，有效支持线性代数，矩阵和矢量运算，数值分析及其相关的算法。
//int main()中为常用规则   外部函数为部分具体功能（PCA、文本读写）服务




//以下4个函数为PCA算法功能函数
//样本均值化
void featurenormalize(Eigen::MatrixXd &X)  
{
	//样本均值化
   //计算每一维度均值
	Eigen::MatrixXd meanval = X.colwise().mean();
	Eigen::RowVectorXd meanvecRow = meanval;
	//样本均值化为0
	X.rowwise() -= meanvecRow;

}
//计算协方差
void computeCov(Eigen::MatrixXd &X, Eigen::MatrixXd &C)
{
	//计算协方差矩阵C = XTX / n-1;
	C = X.adjoint() * X;
	C = C.array() / X.rows() - 1;

}
//计算特征值和特征向量
void computeEig(Eigen::MatrixXd &C, Eigen::MatrixXd &vec, Eigen::MatrixXd &val)
{
	//计算特征值和特征向量，使用selfadjont按照对阵矩阵的算法去计算，可以让产生的vec和val按照有序排列
	Eigen::SelfAdjointEigenSolver<Eigen::MatrixXd> eig(C);
	vec = eig.eigenvectors();
	val = eig.eigenvalues();
}
//获取所需的特征向量数
int computeDim(Eigen::MatrixXd &val)
{
	//计算损失率，确定降低维数
	int dim;
	double sum = 0;
	for (int i = val.rows() - 1; i >= 0; --i)
	{
		sum += val(i, 0);
		dim = i;
		if (sum / val.sum() >= 0.95)
			break;
	}
	return val.rows() - dim;
}



//以下2个函数仅服务于eigen对txt文本读写   
//确定文本的行数
int getFileRows(const char *fileName)
{
	int count = 1;
	std::ifstream fileStream; 
	std::string tmp;// 行数计数器
	fileStream.open(fileName, std::ios::in);//ios::in 表示以只读的方式读取文件
	if(!fileStream)//文件打开失败:返回-1
	{        	
		return -1;
	}    
	else//文件存在    
	{
		char c;
		while (fileStream.get(c))
		{
			if (c == '\n')
			{
				count++;
			}
		}
		fileStream.close();
		return count;

	}

}
//确定文本的列数
int getFileColumns(const char * fileName) 
{
	int count = 0;//计数器
	std::ifstream fileStream;
	fileStream.open(fileName, std::ios::_Nocreate);
	double tmp = 0;		// 列数计数器
	char c;			//当前位置的字符
	c = fileStream.peek();	
	while (('\n' != c) && (! fileStream.eof()))	// 指针指向的当前字符，仅观测，不移动指针位置   
	{          
		fileStream >> tmp;
		count++;
		c = fileStream.peek();
	}  	    
	fileStream.close();
	return count;
}



int main()
{
	//定义
	Eigen::Matrix3d a;             //定义三维（3*3）double型 矩阵； 注意”d”表示double类型，”f”表示float类型，”i”表示整数，”c”表示复数；
	a <<2,2,1,3,1,5,3,2,3;       //赋值操作
	Eigen::Vector4f b(1, 2, 3, 4);           //定义4*1 float型 列向量,并赋值； 赋值也可以： b<<1, 2, 3, 4;    
	Eigen::RowVector2d c(1, 2);                    //   定义1*2 double型 行向量并赋值
	Eigen::MatrixXd m_1 = Eigen::MatrixXd::Random(3, 4);  //MatrixXd->动态   MatrixXd::Random(3,4)表示产生一个元素类型为double的3*4的临时矩阵对象
	Eigen::MatrixXd m_3(2, 2); m_3 << 1, 2, 3, 4;          //->动态  初始化（设置临时2*2矩阵）并赋值
	Eigen::Matrix3d m_2 = Eigen::Matrix3d::Random();      //Matrix3d->静态
    Eigen::MatrixXd m_4;
    m_4=Eigen::MatrixXd(100,90);//设置矩阵大小
    m_4=Eigen::MatrixXd::Ones(100,90);//矩阵元素都设置为1
    m_4=Eigen::MatrixXd::Zero(100,90);//矩阵元素都设置为0
    m_4=Eigen::MatrixXd::Identity(100,100);//矩阵设置为单位阵
    //计算
	a * m_2;//矩阵乘法
	a + m_2;//矩阵加法
	a - m_2;//矩阵减法
	a.MatrixBase::transpose();//转置矩阵transpose(),
	a.determinant();//矩阵的行列式
	a.determinant()*a.inverse();//伴随矩阵   发现伴随矩阵函数adjoint() 使用有问题！！！！！！！
	a.inverse();//矩阵的逆inverse()
	a.diagonal();//取对角元素 组成一个列向量
	a.conjugate();//共轭矩阵
	//这些函数返回操作后的结果，而不会对原矩阵的元素进行直接操作，如果要让原矩阵进行转换，则需要使用响应的InPlace函数，如transpoceInPlace()等
		
	//a.arry()：取出矩阵的所有元素，对每个元素进行操作
	a.array()+1.9;//元素都加上1.9  错误使用：a=a+1.9 !!!
	a.array() -= 2;//  即a.array()-2
	a.array().sin();//元素取sin
	a.array().exp();//指数e
	a.array().log();//取对数
	a.array().pow(2);//元素 ^2    	std::cout << a;
	a.array() / m_2.array();//两个矩阵对应点相除（对应点！！！）
	a.array()*m_2.array();//两个矩阵对应点相除（对应点！！！）
	
	
	//分块矩阵
	a.block<2,3>(1,0);  // 提出一个2*3分块矩阵，起始位置下标为原有矩阵(1,0)
	a.topRightCorner(2, 2);// 取右上角2*2矩阵
	a.bottomLeftCorner(2, 2);//取左下角角2*2矩阵  类似同理
	a.leftCols(2);//取左边2列
	a.rightCols(1);//取右边1列
	a.middleCols(0, 2);//中间从第0列向右数，取2列
	a.topRows(2);//从上向下取2行
	a.bottomRows(2);//从下向上取2行
	a.middleRows(0, 1);// 中间从第0行向下数，取1列
	
	//分块特例-vector 即向量的块操作
	b.head(2);     //获取向量的前n个元素
	b.tail(2);     //获取向量尾部的n个元素
	b.segment(1, 3);        //从第1元素开始 取3个元素   
	
	a.cwiseAbs();//取绝对值
	a.cwiseAbs2();  //元素平方（cwiseSqrt）
	Eigen::Vector3f v1(0, 2, 3), v2(1, 1, 4);
	v1.cwiseMin(v2);  //取两个中的min并赋给v1  同理cwiseMax取max	
	v1.cwiseProduct(v2);  //同位置元素相乘并赋给v1   1，2，12
	v1.cwiseQuotient(v2); //同位置元素相除并赋给v1   0，2，0.75 （ v1/v2 ）
		
	a.colwise().minCoeff();        
	//colwise()把每一列变成了一个块，然后这些块组成了一个行向量即分别对每列进行处理  colwise.minCoeff() 取每列最小值   思想:降维->行向量
	a.rowwise().maxCoeff();        //同上   取行最大值 降维->列向量
	a.rowwise().mean();//mean 矩阵元素平均数

	a.sum(); //所有系数的和
	a.prod();//所有系数的积
	a.all();  //且运算 
	a.any(); //或运算      
	a.col(0); a.row(2);     //取出第几列、行   differ cols() and rows()
	a.cast<int>();//强制转换int型 向下取整
    a.reverse();//元素倒置(i,j)->(j,i)

	//特征值与特征向量 计算
	Eigen::Matrix3d chara_eigen;
	chara_eigen << 4, 6, 0, -3, -5, 0, -3, -6, 1;
	Eigen::EigenSolver<Eigen::Matrix3d> eigen_solver(chara_eigen);
	eigen_solver.eigenvalues();//求特征值   形式为二维向量(1,0)(1,0)(2,0)。真实值为1,1,2。
	eigen_solver.eigenvectors() ;//求特征向量  形式同上  结果已单位化


	//PCA矩阵算法使用
	//PCA算法 PCA的目的：能够对一个多样本的多维特征向量构成的矩阵进行分析，分析主成分并去除维度之间的相关性，
	//使线性相关的向量组变成线性无关的向量组。  并且可以对样本进行降维，降高维向量映射到低维度空间，同时确保维度之间的信息损失尽可能小。
	Eigen::MatrixXd Xchara_eigen(3,3);//样本矩阵
	Xchara_eigen<< -1, 1, 0, -4, 3, 0, 1, 0, 2;
	Eigen::MatrixXd cov_eigen(3,3); //协方差
	Eigen::MatrixXd vector_eigen;  //特征向量
	Eigen::MatrixXd value_eigen;  //特征值
	featurenormalize(Xchara_eigen);	//样本均值化
	computeCov(Xchara_eigen, cov_eigen);	//计算协方差
	computeEig(cov_eigen, vector_eigen, value_eigen);	//计算特征值和特征向量
	int dim = computeDim(value_eigen);//生成所需特征向量
	Eigen::MatrixXd result_eigen = Xchara_eigen * vector_eigen.rightCols(dim);//生成目标矩阵


	//SVD矩阵分解
	Eigen::JacobiSVD<Eigen::MatrixXd> svd(Xchara_eigen, Eigen::ComputeThinU | Eigen::ComputeThinV);
	Eigen::Matrix3d V = svd.matrixV(), U = svd.matrixU();  //矩阵V
	Eigen::Matrix3d  S = U.inverse() * Xchara_eigen * V.transpose().inverse(); // S = U^-1 * A * VT * -1  
	

	
	//稀疏矩阵的相关操作 添加#include <Eigen/Sparse>
	Eigen::SparseMatrix<double> SpMat; // 声明一个列优先的双精度稀疏矩阵类型
	Eigen::SparseMatrix<double, Eigen::RowMajor> mat(1000, 2000);//声明一个 1000x2000 行优先的压缩的稀疏矩阵（元素类型：double）
	Eigen::Triplet<double> T; //三元组（行，列，值）
	mat.innerSize(); //内维度   行优先，所以为列数
	mat.outerSize(); //外维度
	mat.selfadjointView<Eigen::Upper>();//从上三角部分生成一个全的伴随矩阵
	mat.selfadjointView<Eigen::Lower>();//从下三角部分生成一个全的伴随矩阵
	mat.setZero();//去除所有的非零元素 即矩阵非零元素都变为零
	mat.nonZeros();     // 非零元素的个数   
	mat.outerSize();    // 列优先的稀疏矩阵的行数，行优先的稀疏矩阵的列数
	mat.innerSize();    // 列优先的稀疏矩阵的列数，行优先的稀疏矩阵的行数
	mat.norm();          // 矩阵的欧式范数
	mat.squaredNorm();  //矩阵的平方范数（平方和）
	mat.isVector();     // 检查是否是稀疏向量或稀疏矩阵
	mat.isCompressed(); // 检查是否是压缩格式
	mat.makeCompressed(); //压缩剩余的空间
	mat.coeffRef(0, 0); //插入值


//txt 矩阵文件读写
	int cols_num = 0;int rows_num = 0;
	std::ifstream infile("D:\\STUDY\\Computer\\C_training\\egine\\Matrix-test\\data.txt");  //读取txt路径
		if (infile.is_open())
	{
		std::cout << "Open file successfully." << std::endl;
		rows_num =getFileRows("D:\\STUDY\\Computer\\C_training\\egine\\Matrix-test\\data.txt");//得到矩阵行数
		cols_num = getFileColumns("D:\\STUDY\\Computer\\C_training\\egine\\Matrix-test\\data.txt");//得到矩阵列数
		
	}
	else
	{
		std::cout << "Unable to open file" << std::endl;
	}
	
	Eigen::MatrixXd set_mat = Eigen::MatrixXd(rows_num, cols_num);  //生成目标矩阵（确定行列）
	double mat_data;
	std::vector<double> data_save;
	while (infile >> mat_data)
	{
		data_save.push_back(mat_data);  //存入所有元素于vector
	}
	int count = 0;
	for (int row = 0; row < rows_num; row++)
	{
		for (int col = 0; col < cols_num; col++)
		{
			set_mat(row, col) = data_save[count]; //元素转存入目标矩阵
			count++;
		}
	}
	std::cout << "\nshow:\n";
	std::cout << set_mat; 
	
        return 0;

}