矩阵相乘只算出对角线的值

tho0o 2016-12-26 09:27:50
加精
有两个矩阵 A和B 他们都有几万行和列，我想只求出他们相乘对角线的值，但是我不想求出相乘后的矩阵再输出对角线（这样很浪费时间），我仅仅想做的是只算出对角线的值，然后输出。

用c++实现，有木有大神啊！
...全文
5876 30 打赏 收藏 转发到动态 举报
写回复
用AI写文章
30 条回复
切换为时间正序
请发表友善的回复…
发表回复
Junlys 2017-01-14
  • 打赏
  • 举报
 回复
积分怎么拿～积分怎么拿～积分怎么拿～积分怎么拿～积分怎么拿～
副组长 2017-01-02
  • 打赏
  • 举报
 回复
不错不错，是好东西，谢谢分享。
lkj2016 2016-12-31
  • 打赏
  • 举报
 回复
膜拜一下大神
NoEdUl 2016-12-30
  • 打赏
  • 举报
 回复
引用 26 楼 lunat 的回复:
题主问的是对角线计算，楼歪了。
无心插柳
ly_490380384 2016-12-30
  • 打赏
  • 举报
 回复
666666666666666666
ljheee 2016-12-29
  • 打赏
  • 举报
 回复
给下你完整代码
paschen 2016-12-29
  • 打赏
  • 举报
 回复
引用 20 楼 u012909435 的回复:
[quote=引用 18 楼 paschen 的回复:] [quote=引用 17 楼 u012947309 的回复:] [quote=引用 11 楼 paschen 的回复:] [quote=引用 10 楼 u012947309 的回复:] [quote=引用 5 楼 paschen 的回复:] 一个高效计算矩阵相乘的代码（看上去很复杂，但效率很高，利用了CPU缓存进行优化）： 
	
#define CACHE_LINE_SIZE 64
#define SM (CACHE_LINE_SIZE / sizeof(double))
	int i2, j2, k2;
	double* rres, *rm1, *rm2;
	for (int i = 0; i < N; i += SM)
		for (int j = 0; j < N; j += SM)
			for (int k = 0; k < N; k += SM)
				for (i2 = 0, rres = &res3[i][j], rm1 = &matrix1[i][k]; i2 < SM; ++i2, rres +=N, rm1 += N)
					for (k2 = 0, rm2 = &matrix2[k][j]; k2 < SM; ++k2, rm2 += N)
						for (j2 = 0; j2 < SM; ++j2)
							rres[j2] += rm1[k2] * rm2[j2];
计算结果在res3中，相乘的矩阵为matrix1与matrix2
然而并没有看懂.......至少跑了一下结果不正确诶 [/quote] 你应该没搞对，给下你完整代码[/quote] 

#include <vector>
#include <iostream>
using namespace std;

#define CACHE_LINE_SIZE 64
#define SM (CACHE_LINE_SIZE / sizeof(double))

int main()
{
	double matrix1[5][5];
	double matrix2[5][5];
	double res3[5][5];
	//initialize
	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			matrix1[i][j] = i*j;
			matrix2[i][j] = i*j;
		}
	}
	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			cout<<matrix1[i][j]<<" ";
		}
		cout<<endl;
	}
	cout<<endl;
	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			cout<<matrix2[i][j]<<" ";
		}
		cout<<endl;
	}
	cout<<endl;

	int N=5;
	int i2, j2, k2;
	double* rres, *rm1, *rm2;
	for (int i = 0; i < N; i += SM)
		for (int j = 0; j < N; j += SM)
			for (int k = 0; k < N; k += SM)
				for (i2 = 0, rres = &res3[i][j], rm1 = &matrix1[i][k]; i2 < SM; ++i2, rres +=N, rm1 += N)
					for (k2 = 0, rm2 = &matrix2[k][j]; k2 < SM; ++k2, rm2 += N)
						for (j2 = 0; j2 < SM; ++j2)
							rres[j2] += rm1[k2] * rm2[j2];

	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			cout<<res3[i][j]<<" ";
		}
		cout<<endl;
	}

	return 0;
}
[/quote] 

#include <vector>
#include <iostream>
using namespace std;

#define CACHE_LINE_SIZE 64
#define SM (CACHE_LINE_SIZE / sizeof(double))
#define N 8

int main()
{
	double matrix1[N][N];
	double matrix2[N][N];
	double res3[N][N]={0};
	//initialize
	for(int i=0;i<N;i++)
	{
		for(int j=0;j<N;j++)
		{
			matrix1[i][j] = i*j;
			matrix2[i][j] = i*j;
		}
	}
	for(int i=0;i<N;i++)
	{
		for(int j=0;j<N;j++)
		{
			cout<<matrix1[i][j]<<" ";
		}
		cout<<endl;
	}
	cout<<endl;
	for(int i=0;i<N;i++)
	{
		for(int j=0;j<N;j++)
		{
			cout<<matrix2[i][j]<<" ";
		}
		cout<<endl;
	}
	cout<<endl;

	int i2, j2, k2;
	double* rres, *rm1, *rm2;
	for (int i = 0; i < N; i += SM)
		for (int j = 0; j < N; j += SM)
			for (int k = 0; k < N; k += SM)
				for (i2 = 0, rres = &res3[i][j], rm1 = &matrix1[i][k]; i2 < SM; ++i2, rres +=N, rm1 += N)
					for (k2 = 0, rm2 = &matrix2[k][j]; k2 < SM; ++k2, rm2 += N)
						for (j2 = 0; j2 < SM; ++j2)
							rres[j2] += rm1[k2] * rm2[j2];

	for(int i=0;i<N;i++)
	{
		for(int j=0;j<N;j++)
		{
			cout<<res3[i][j]<<" ";
		}
		cout<<endl;
	}

	return 0;
}
矩阵维度要大于8，这种优化主要是减少cache miss，使用最大的cache line size 其次res3先初始化为0 矩阵较大时这种方式速度明显要快 在我机上测试1000*1000的矩阵乘法，用普通的方法用时13.17s，这种方法仅1.49s[/quote] 这个要求N能被SM整除吧 后面三个循环的跳出条件要改成 i2 < SM && i + i2 < N k2 < SM && k + k2 < N j2 < SM && j + j2 < N 这样是不是能处理不整除的情况了？[/quote] 是，一言中的
paschen 2016-12-29
  • 打赏
  • 举报
 回复
引用 22 楼 u012909435 的回复:
[quote=引用 5 楼 paschen 的回复:] 一个高效计算矩阵相乘的代码（看上去很复杂，但效率很高，利用了CPU缓存进行优化）： 
	
#define CACHE_LINE_SIZE 64
#define SM (CACHE_LINE_SIZE / sizeof(double))
	int i2, j2, k2;
	double* rres, *rm1, *rm2;
	for (int i = 0; i < N; i += SM)
		for (int j = 0; j < N; j += SM)
			for (int k = 0; k < N; k += SM)
				for (i2 = 0, rres = &res3[i][j], rm1 = &matrix1[i][k]; i2 < SM; ++i2, rres +=N, rm1 += N)
					for (k2 = 0, rm2 = &matrix2[k][j]; k2 < SM; ++k2, rm2 += N)
						for (j2 = 0; j2 < SM; ++j2)
							rres[j2] += rm1[k2] * rm2[j2];
计算结果在res3中，相乘的矩阵为matrix1与matrix2
我试了一下，按照您这个方法 OpenMP+cache blocking+loop interchange+vectorize 在Xeon X5675上测试 4096*4096矩阵乘法 trival的实现大概是1171s 并行优化之后是4s 在12个thread的时候，加速了200多倍 4096*4096的矩阵乘法应该有64G次乘法+64G次加法。。。。 不过，这个性能算是达到了16Gflop/s吗，还是32Gflop/s，还是咋算。。。[/quote] 具体是不是这种算我也说不上来，你可以用LinX工具来测试FLOPS
cattpon 2016-12-29
  • 打赏
  • 举报
 回复
看看是什么~
xxiaoccen 2016-12-29
  • 打赏
  • 举报
 回复
引用 5 楼 paschen 的回复:
一个高效计算矩阵相乘的代码（看上去很复杂，但效率很高，利用了CPU缓存进行优化）： 
	
#define CACHE_LINE_SIZE 64
#define SM (CACHE_LINE_SIZE / sizeof(double))
	int i2, j2, k2;
	double* rres, *rm1, *rm2;
	for (int i = 0; i < N; i += SM)
		for (int j = 0; j < N; j += SM)
			for (int k = 0; k < N; k += SM)
				for (i2 = 0, rres = &res3[i][j], rm1 = &matrix1[i][k]; i2 < SM; ++i2, rres +=N, rm1 += N)
					for (k2 = 0, rm2 = &matrix2[k][j]; k2 < SM; ++k2, rm2 += N)
						for (j2 = 0; j2 < SM; ++j2)
							rres[j2] += rm1[k2] * rm2[j2];
计算结果在res3中，相乘的矩阵为matrix1与matrix2
我试了一下，按照您这个方法 OpenMP+cache blocking+loop interchange+vectorize 在Xeon X5675上测试 4096*4096矩阵乘法 trival的实现大概是1171s 并行优化之后是4s 在12个thread的时候，加速了200多倍 4096*4096的矩阵乘法应该有64G次乘法+64G次加法。。。。 不过，这个性能算是达到了16Gflop/s吗，还是32Gflop/s，还是咋算。。。
lunat 2016-12-29
  • 打赏
  • 举报
 回复
题主问的是对角线计算，楼歪了。
GaryCV 2016-12-28
  • 打赏
  • 举报
 回复
厉害厉害，参考
xxiaoccen 2016-12-28
  • 打赏
  • 举报
 回复
引用 18 楼 paschen 的回复:
[quote=引用 17 楼 u012947309 的回复:] [quote=引用 11 楼 paschen 的回复:] [quote=引用 10 楼 u012947309 的回复:] [quote=引用 5 楼 paschen 的回复:] 一个高效计算矩阵相乘的代码（看上去很复杂，但效率很高，利用了CPU缓存进行优化）： 
	
#define CACHE_LINE_SIZE 64
#define SM (CACHE_LINE_SIZE / sizeof(double))
	int i2, j2, k2;
	double* rres, *rm1, *rm2;
	for (int i = 0; i < N; i += SM)
		for (int j = 0; j < N; j += SM)
			for (int k = 0; k < N; k += SM)
				for (i2 = 0, rres = &res3[i][j], rm1 = &matrix1[i][k]; i2 < SM; ++i2, rres +=N, rm1 += N)
					for (k2 = 0, rm2 = &matrix2[k][j]; k2 < SM; ++k2, rm2 += N)
						for (j2 = 0; j2 < SM; ++j2)
							rres[j2] += rm1[k2] * rm2[j2];
计算结果在res3中，相乘的矩阵为matrix1与matrix2
然而并没有看懂.......至少跑了一下结果不正确诶 [/quote] 你应该没搞对，给下你完整代码[/quote] 

#include <vector>
#include <iostream>
using namespace std;

#define CACHE_LINE_SIZE 64
#define SM (CACHE_LINE_SIZE / sizeof(double))

int main()
{
	double matrix1[5][5];
	double matrix2[5][5];
	double res3[5][5];
	//initialize
	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			matrix1[i][j] = i*j;
			matrix2[i][j] = i*j;
		}
	}
	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			cout<<matrix1[i][j]<<" ";
		}
		cout<<endl;
	}
	cout<<endl;
	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			cout<<matrix2[i][j]<<" ";
		}
		cout<<endl;
	}
	cout<<endl;

	int N=5;
	int i2, j2, k2;
	double* rres, *rm1, *rm2;
	for (int i = 0; i < N; i += SM)
		for (int j = 0; j < N; j += SM)
			for (int k = 0; k < N; k += SM)
				for (i2 = 0, rres = &res3[i][j], rm1 = &matrix1[i][k]; i2 < SM; ++i2, rres +=N, rm1 += N)
					for (k2 = 0, rm2 = &matrix2[k][j]; k2 < SM; ++k2, rm2 += N)
						for (j2 = 0; j2 < SM; ++j2)
							rres[j2] += rm1[k2] * rm2[j2];

	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			cout<<res3[i][j]<<" ";
		}
		cout<<endl;
	}

	return 0;
}
[/quote] 

#include <vector>
#include <iostream>
using namespace std;

#define CACHE_LINE_SIZE 64
#define SM (CACHE_LINE_SIZE / sizeof(double))
#define N 8

int main()
{
	double matrix1[N][N];
	double matrix2[N][N];
	double res3[N][N]={0};
	//initialize
	for(int i=0;i<N;i++)
	{
		for(int j=0;j<N;j++)
		{
			matrix1[i][j] = i*j;
			matrix2[i][j] = i*j;
		}
	}
	for(int i=0;i<N;i++)
	{
		for(int j=0;j<N;j++)
		{
			cout<<matrix1[i][j]<<" ";
		}
		cout<<endl;
	}
	cout<<endl;
	for(int i=0;i<N;i++)
	{
		for(int j=0;j<N;j++)
		{
			cout<<matrix2[i][j]<<" ";
		}
		cout<<endl;
	}
	cout<<endl;

	int i2, j2, k2;
	double* rres, *rm1, *rm2;
	for (int i = 0; i < N; i += SM)
		for (int j = 0; j < N; j += SM)
			for (int k = 0; k < N; k += SM)
				for (i2 = 0, rres = &res3[i][j], rm1 = &matrix1[i][k]; i2 < SM; ++i2, rres +=N, rm1 += N)
					for (k2 = 0, rm2 = &matrix2[k][j]; k2 < SM; ++k2, rm2 += N)
						for (j2 = 0; j2 < SM; ++j2)
							rres[j2] += rm1[k2] * rm2[j2];

	for(int i=0;i<N;i++)
	{
		for(int j=0;j<N;j++)
		{
			cout<<res3[i][j]<<" ";
		}
		cout<<endl;
	}

	return 0;
}
矩阵维度要大于8，这种优化主要是减少cache miss，使用最大的cache line size 其次res3先初始化为0 矩阵较大时这种方式速度明显要快 在我机上测试1000*1000的矩阵乘法，用普通的方法用时13.17s，这种方法仅1.49s[/quote] 这个要求N能被SM整除吧 后面三个循环的跳出条件要改成 i2 < SM && i + i2 < N k2 < SM && k + k2 < N j2 < SM && j + j2 < N 这样是不是能处理不整除的情况了？
tho0o 2016-12-28
  • 打赏
  • 举报
 回复
引用 17 楼 u012947309 的回复:
[quote=引用 11 楼 paschen 的回复:] [quote=引用 10 楼 u012947309 的回复:] [quote=引用 5 楼 paschen 的回复:] 一个高效计算矩阵相乘的代码（看上去很复杂，但效率很高，利用了CPU缓存进行优化）： 
	
#define CACHE_LINE_SIZE 64
#define SM (CACHE_LINE_SIZE / sizeof(double))
	int i2, j2, k2;
	double* rres, *rm1, *rm2;
	for (int i = 0; i < N; i += SM)
		for (int j = 0; j < N; j += SM)
			for (int k = 0; k < N; k += SM)
				for (i2 = 0, rres = &res3[i][j], rm1 = &matrix1[i][k]; i2 < SM; ++i2, rres +=N, rm1 += N)
					for (k2 = 0, rm2 = &matrix2[k][j]; k2 < SM; ++k2, rm2 += N)
						for (j2 = 0; j2 < SM; ++j2)
							rres[j2] += rm1[k2] * rm2[j2];
计算结果在res3中，相乘的矩阵为matrix1与matrix2
然而并没有看懂.......至少跑了一下结果不正确诶 [/quote] 你应该没搞对，给下你完整代码[/quote] 

#include <vector>
#include <iostream>
using namespace std;

#define CACHE_LINE_SIZE 64
#define SM (CACHE_LINE_SIZE / sizeof(double))

int main()
{
	double matrix1[5][5];
	double matrix2[5][5];
	double res3[5][5];
	//initialize
	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			matrix1[i][j] = i*j;
			matrix2[i][j] = i*j;
		}
	}
	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			cout<<matrix1[i][j]<<" ";
		}
		cout<<endl;
	}
	cout<<endl;
	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			cout<<matrix2[i][j]<<" ";
		}
		cout<<endl;
	}
	cout<<endl;

	int N=5;
	int i2, j2, k2;
	double* rres, *rm1, *rm2;
	for (int i = 0; i < N; i += SM)
		for (int j = 0; j < N; j += SM)
			for (int k = 0; k < N; k += SM)
				for (i2 = 0, rres = &res3[i][j], rm1 = &matrix1[i][k]; i2 < SM; ++i2, rres +=N, rm1 += N)
					for (k2 = 0, rm2 = &matrix2[k][j]; k2 < SM; ++k2, rm2 += N)
						for (j2 = 0; j2 < SM; ++j2)
							rres[j2] += rm1[k2] * rm2[j2];

	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			cout<<res3[i][j]<<" ";
		}
		cout<<endl;
	}

	return 0;
}
[/quote] 大神受教了，是不是还有SSE加速的
paschen 2016-12-28
  • 打赏
  • 举报
 回复
引用 17 楼 u012947309 的回复:
[quote=引用 11 楼 paschen 的回复:] [quote=引用 10 楼 u012947309 的回复:] [quote=引用 5 楼 paschen 的回复:] 一个高效计算矩阵相乘的代码（看上去很复杂，但效率很高，利用了CPU缓存进行优化）： 
	
#define CACHE_LINE_SIZE 64
#define SM (CACHE_LINE_SIZE / sizeof(double))
	int i2, j2, k2;
	double* rres, *rm1, *rm2;
	for (int i = 0; i < N; i += SM)
		for (int j = 0; j < N; j += SM)
			for (int k = 0; k < N; k += SM)
				for (i2 = 0, rres = &res3[i][j], rm1 = &matrix1[i][k]; i2 < SM; ++i2, rres +=N, rm1 += N)
					for (k2 = 0, rm2 = &matrix2[k][j]; k2 < SM; ++k2, rm2 += N)
						for (j2 = 0; j2 < SM; ++j2)
							rres[j2] += rm1[k2] * rm2[j2];
计算结果在res3中，相乘的矩阵为matrix1与matrix2
然而并没有看懂.......至少跑了一下结果不正确诶 [/quote] 你应该没搞对，给下你完整代码[/quote] 

#include <vector>
#include <iostream>
using namespace std;

#define CACHE_LINE_SIZE 64
#define SM (CACHE_LINE_SIZE / sizeof(double))

int main()
{
	double matrix1[5][5];
	double matrix2[5][5];
	double res3[5][5];
	//initialize
	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			matrix1[i][j] = i*j;
			matrix2[i][j] = i*j;
		}
	}
	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			cout<<matrix1[i][j]<<" ";
		}
		cout<<endl;
	}
	cout<<endl;
	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			cout<<matrix2[i][j]<<" ";
		}
		cout<<endl;
	}
	cout<<endl;

	int N=5;
	int i2, j2, k2;
	double* rres, *rm1, *rm2;
	for (int i = 0; i < N; i += SM)
		for (int j = 0; j < N; j += SM)
			for (int k = 0; k < N; k += SM)
				for (i2 = 0, rres = &res3[i][j], rm1 = &matrix1[i][k]; i2 < SM; ++i2, rres +=N, rm1 += N)
					for (k2 = 0, rm2 = &matrix2[k][j]; k2 < SM; ++k2, rm2 += N)
						for (j2 = 0; j2 < SM; ++j2)
							rres[j2] += rm1[k2] * rm2[j2];

	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			cout<<res3[i][j]<<" ";
		}
		cout<<endl;
	}

	return 0;
}
[/quote] 

#include <vector>
#include <iostream>
using namespace std;

#define CACHE_LINE_SIZE 64
#define SM (CACHE_LINE_SIZE / sizeof(double))
#define N 8

int main()
{
	double matrix1[N][N];
	double matrix2[N][N];
	double res3[N][N]={0};
	//initialize
	for(int i=0;i<N;i++)
	{
		for(int j=0;j<N;j++)
		{
			matrix1[i][j] = i*j;
			matrix2[i][j] = i*j;
		}
	}
	for(int i=0;i<N;i++)
	{
		for(int j=0;j<N;j++)
		{
			cout<<matrix1[i][j]<<" ";
		}
		cout<<endl;
	}
	cout<<endl;
	for(int i=0;i<N;i++)
	{
		for(int j=0;j<N;j++)
		{
			cout<<matrix2[i][j]<<" ";
		}
		cout<<endl;
	}
	cout<<endl;

	int i2, j2, k2;
	double* rres, *rm1, *rm2;
	for (int i = 0; i < N; i += SM)
		for (int j = 0; j < N; j += SM)
			for (int k = 0; k < N; k += SM)
				for (i2 = 0, rres = &res3[i][j], rm1 = &matrix1[i][k]; i2 < SM; ++i2, rres +=N, rm1 += N)
					for (k2 = 0, rm2 = &matrix2[k][j]; k2 < SM; ++k2, rm2 += N)
						for (j2 = 0; j2 < SM; ++j2)
							rres[j2] += rm1[k2] * rm2[j2];

	for(int i=0;i<N;i++)
	{
		for(int j=0;j<N;j++)
		{
			cout<<res3[i][j]<<" ";
		}
		cout<<endl;
	}

	return 0;
}
矩阵维度要大于8，这种优化主要是减少cache miss，使用最大的cache line size 其次res3先初始化为0 矩阵较大时这种方式速度明显要快 在我机上测试1000*1000的矩阵乘法，用普通的方法用时13.17s，这种方法仅1.49s
NoEdUl 2016-12-28
  • 打赏
  • 举报
 回复
引用 11 楼 paschen 的回复:
[quote=引用 10 楼 u012947309 的回复:] [quote=引用 5 楼 paschen 的回复:] 一个高效计算矩阵相乘的代码（看上去很复杂，但效率很高，利用了CPU缓存进行优化）： 
	
#define CACHE_LINE_SIZE 64
#define SM (CACHE_LINE_SIZE / sizeof(double))
	int i2, j2, k2;
	double* rres, *rm1, *rm2;
	for (int i = 0; i < N; i += SM)
		for (int j = 0; j < N; j += SM)
			for (int k = 0; k < N; k += SM)
				for (i2 = 0, rres = &res3[i][j], rm1 = &matrix1[i][k]; i2 < SM; ++i2, rres +=N, rm1 += N)
					for (k2 = 0, rm2 = &matrix2[k][j]; k2 < SM; ++k2, rm2 += N)
						for (j2 = 0; j2 < SM; ++j2)
							rres[j2] += rm1[k2] * rm2[j2];
计算结果在res3中，相乘的矩阵为matrix1与matrix2
然而并没有看懂.......至少跑了一下结果不正确诶 [/quote] 你应该没搞对，给下你完整代码[/quote] 

#include <vector>
#include <iostream>
using namespace std;

#define CACHE_LINE_SIZE 64
#define SM (CACHE_LINE_SIZE / sizeof(double))

int main()
{
	double matrix1[5][5];
	double matrix2[5][5];
	double res3[5][5];
	//initialize
	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			matrix1[i][j] = i*j;
			matrix2[i][j] = i*j;
		}
	}
	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			cout<<matrix1[i][j]<<" ";
		}
		cout<<endl;
	}
	cout<<endl;
	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			cout<<matrix2[i][j]<<" ";
		}
		cout<<endl;
	}
	cout<<endl;

	int N=5;
	int i2, j2, k2;
	double* rres, *rm1, *rm2;
	for (int i = 0; i < N; i += SM)
		for (int j = 0; j < N; j += SM)
			for (int k = 0; k < N; k += SM)
				for (i2 = 0, rres = &res3[i][j], rm1 = &matrix1[i][k]; i2 < SM; ++i2, rres +=N, rm1 += N)
					for (k2 = 0, rm2 = &matrix2[k][j]; k2 < SM; ++k2, rm2 += N)
						for (j2 = 0; j2 < SM; ++j2)
							rres[j2] += rm1[k2] * rm2[j2];

	for(int i=0;i<5;i++)
	{
		for(int j=0;j<5;j++)
		{
			cout<<res3[i][j]<<" ";
		}
		cout<<endl;
	}

	return 0;
}
xxiaoccen 2016-12-28
  • 打赏
  • 举报
 回复 

for(int i = 0; i < nrows; i += SM) {
    for(int j = 0; j < ncols; j += SM) {
        for(int ii = 0; ii < SM; ii++) {
            for(int jj = 0; jj < SM; jj++) {
                diag[i + ii] += a[i + ii][j + jj] * b[j + jj][i + ii];
            }
       }
    }
}
只要对角线
xxiaoccen 2016-12-28
  • 打赏
  • 举报
 回复
引用 11 楼 paschen 的回复:
[quote=引用 10 楼 u012947309 的回复:] [quote=引用 5 楼 paschen 的回复:] 一个高效计算矩阵相乘的代码（看上去很复杂，但效率很高，利用了CPU缓存进行优化）： 
	
#define CACHE_LINE_SIZE 64
#define SM (CACHE_LINE_SIZE / sizeof(double))
	int i2, j2, k2;
	double* rres, *rm1, *rm2;
	for (int i = 0; i < N; i += SM)
		for (int j = 0; j < N; j += SM)
			for (int k = 0; k < N; k += SM)
				for (i2 = 0, rres = &res3[i][j], rm1 = &matrix1[i][k]; i2 < SM; ++i2, rres +=N, rm1 += N)
					for (k2 = 0, rm2 = &matrix2[k][j]; k2 < SM; ++k2, rm2 += N)
						for (j2 = 0; j2 < SM; ++j2)
							rres[j2] += rm1[k2] * rm2[j2];
计算结果在res3中，相乘的矩阵为matrix1与matrix2
然而并没有看懂.......至少跑了一下结果不正确诶 [/quote] 你应该没搞对，给下你完整代码[/quote] 这是cache blocking吗？
NoEdUl 2016-12-27
  • 打赏
  • 举报
 回复
引用 5 楼 paschen 的回复:
一个高效计算矩阵相乘的代码（看上去很复杂，但效率很高，利用了CPU缓存进行优化）：

	

#define CACHE_LINE_SIZE 64

#define SM (CACHE_LINE_SIZE / sizeof(double))

	int i2, j2, k2;

	double* rres, *rm1, *rm2;

	for (int i = 0; i < N; i += SM)

		for (int j = 0; j < N; j += SM)

			for (int k = 0; k < N; k += SM)

				for (i2 = 0, rres = &res3[i][j], rm1 = &matrix1[i][k]; i2 < SM; ++i2, rres +=N, rm1 += N)

					for (k2 = 0, rm2 = &matrix2[k][j]; k2 < SM; ++k2, rm2 += N)

						for (j2 = 0; j2 < SM; ++j2)

							rres[j2] += rm1[k2] * rm2[j2];


计算结果在res3中，相乘的矩阵为matrix1与matrix2



然而并没有看懂.......至少跑了一下结果不正确诶
NoEdUl 2016-12-27
  • 打赏
  • 举报
 回复
引用 5 楼 paschen 的回复:
一个高效计算矩阵相乘的代码（看上去很复杂，但效率很高，利用了CPU缓存进行优化）： 
	
#define CACHE_LINE_SIZE 64
#define SM (CACHE_LINE_SIZE / sizeof(double))
	int i2, j2, k2;
	double* rres, *rm1, *rm2;
	for (int i = 0; i < N; i += SM)
		for (int j = 0; j < N; j += SM)
			for (int k = 0; k < N; k += SM)
				for (i2 = 0, rres = &res3[i][j], rm1 = &matrix1[i][k]; i2 < SM; ++i2, rres +=N, rm1 += N)
					for (k2 = 0, rm2 = &matrix2[k][j]; k2 < SM; ++k2, rm2 += N)
						for (j2 = 0; j2 < SM; ++j2)
							rres[j2] += rm1[k2] * rm2[j2];
计算结果在res3中，相乘的矩阵为matrix1与matrix2
最近刚好要做大矩阵分解-。-参考一下
加载更多回复（10）
矩阵的生成及对角线运算
在二维数组的使用过程中，我们知道二位数组可表示着矩阵，当然也可表示多个字符串存储在二维数组里。本文主要针对矩阵的赋值进行讨论，从用户输入和随机赋值等方面着手。
矩阵的转置，逆矩阵，行列式的计算，伴随矩阵等
行列式的操作： 逆矩阵：就是两个矩阵相乘是单位矩阵 对角矩阵相乘：就是对角线元素相乘 当两个矩阵相乘不是单位矩阵： 伴随矩阵：是有代数余子式拼成的 为什么伴随矩阵会出现？为什么伴随矩阵的形式是这样的？ 因为行列式的乘法： 根据矩阵的乘法可以看到： 行列式是一个数字，当改行元素跟本身的代数余子式相乘积的情况下，才能非0。绝大情况下是0，只有少数部分情况对上的时候才非0。 所以伴随矩阵出现之...
对角矩阵，投影矩阵，三角矩阵的特征值性质
特征值为0的矩阵是奇异的，也就不可求逆，行列式为0.换句话说如果矩阵所有的特征值都非0的话，那么该矩阵则可以求逆。需要注意的是，矩阵的特征值，矩阵对角线元素的值，矩阵的主元（pivot），这三者是完全不同的概念。利用特征值的乘积也可以验证投影矩阵为奇异矩阵，行列式为0.反过来，矩阵行列式为0，说明至少存在一个特征值为0。根据经验，三角矩阵（triangular）的特征值就是对角线上的元素。很明显矩阵U的特征值就是其对角线的元素，但是和矩阵A的特征值大小是不一样的。需要注意的是，特征值是可以重复的。
【矩阵乘法】逆矩阵
这东西很吼啊，简介有时遇到一些题目，要用到矩阵的前缀积+前缀逆矩阵，来计算区间的矩阵积， 逆矩阵和逆元差不多，两两相乘即可达到除效果，待定系数法直接设逆矩阵为a,b,c,d….，它与原矩阵相乘以后得出的矩阵是 除了对角线为1外，其余为0 的矩阵E，所以直接代入高斯消元解即可， 复杂度：O(n6)O(n^6)初等变换法我们发现，最标准的高斯消元，其实是可以表示为两矩阵相乘的样子的， 那么我们可以大
漫步线性代数二十七——矩阵对角化
现在我们开始实质性的计算，它非常简单并且在随后的几篇文章里都会用到。特征向量对角化一个矩阵：3、假设n×nn\times n矩阵有nn个线性无关的特征向量，如果这些向量是矩阵SS的列，那么S−1ASS^{-1}AS是一个对角矩阵Λ\Lambda，AA的特征值在Λ\Lambda的对角线上： S−1AS=Λ=⎡⎣⎢⎢⎢⎢⎢λ1λ2⋱λn⎤⎦⎥⎥⎥⎥⎥(1)\begin{equation} S^{-1}
模式及实现

5,530

社区成员

4,169

社区内容

发帖
与我相关
我的任务
社区描述
C/C++ 模式及实现
社区管理员
  • 模式及实现社区
加入社区
  • 近7日
  • 近30日
  • 至今

加载中

查看更多榜单
社区公告
暂无公告

试试用AI创作助手写篇文章吧

+ 用AI写文章