多核程序设计的相关基础知识----以误差扩散算法为例

    本文从基础入手，主要阐述基于桌面电脑的多核程序设计的基础知识，包括一些向量化运算，虚拟机算，多线程等的相关知识总结。

一.计算平台的分类

单指令单数据流机器（SISD）

传统的串行计算机，所有指令都是串行。

多指令单数据流机器（MISD）

多个指令流同时对一个数据流进行处理。但是大多数情况下，多个指令流处理多个数据才是更加有效的处理方式。

单指令多数据流机器（SIMD）

几乎所有的计算机都实现了SIMD功能，intel处理器中实现的MMX,SSE,SSE2,SSE3等扩展指令集

说到这里，我就多少说几句，最近在做这方面的优化，发现居然知网上面很多研究 ，SIMD编译优化的，其实debug和release下面的程序运行时间差别很大，visual studio默认开启了很多编译优化，如果对c语言的内部函数不是很熟悉，编译成release版本的程序已经是优化的不错了，但是针对SIMD指令优化最好的莫过于Intel自己家的编译器

windows下SIMD编译优化的例子：

http://bbs.csdn.net/topics/391894458

多指令多数据流机器（MIMD）

能够同时执行多个指令流，这些指令流分别对不同的数据流进行操作。MIMD是目前最流行的并行计算平台。例如 intel core duo双核处理器。

目前的计算机一般都属于SIMD机器或者MIMD机器，而这两种机器都提供了支持并行执行的硬件特性，因此软件开发人员能够非常方便的利用软件中存在的数据级和任务级并行性来提高程序性能。

如果要在应用程序中使用多线程技术，就必须对操作系统的限制有清楚的了解，也就是对系统的api有充分的了解，然而这对于开发通用高性能计算的程序确是一大障碍，我们不能换一个系统，就掌握一套api。

二. 虚拟环境：虚拟机和虚拟平台

在现在很多平台上运行的多线程环境其实是基于虚拟机的，并且目前计算的一个重要趋势是虚拟化。虚拟化技术主要有两种：

1.运行时虚拟化

典型的如：JAVA virtual machine，jvm，微软的通用语言运行时环境Common language runtime CLR

这些虚拟机都至少创建了三个线程：

执行线程

垃圾回收线程

编译线程（just-in-time 即时编译执行技术，将字节码编译成可执行的二进制代码）

一般来讲，这些虚拟机为任务创建的其他进程会以最优化的方式映射到其他可执行资源上。

2.系统虚拟化

VMM virtual machine monitor 虚拟机监视器对底层的平台进行必要的 虚拟化，从而使得每个虚拟机上的操作系统都感觉自己拥有整个硬件资源。

处理器虚拟化技术所带来的一个非常重要的好处就是能够剥离指令集结构（instruction set architecture，ISA）与处理器之间的必然联系。但是特权指令只能由操作系统进行执行。所以需要执行特权指令的时候需要向虚拟机监视器发送请求，得到响应后才能执行，这之中当然降低了效率。因此Intel开发了一系列ISA扩展来支持VMM有效的执行特权指令，这些指令就是大名鼎鼎的Intel虚拟化技术，用来提高VMM的性能（VMware，hamx等技术都依赖这个扩展进行实现，但是似乎和windows 的hyper-x有点冲突，实在是不太明白，还请过来人指点一二：

http://blog.csdn.net/wangyaninglm/article/details/50602482）

三.并行程序设计的基本概念

设计并行程序，程序员应该将应用程序看成是众多相互依赖的任务的集合。将应用程序划分成多个独立的任务，并确定这些任务之间的相互依赖关系的过程就称为分解（decomposition），分解问题 的方式主要有三种：

1.任务分解

2.数据分解

3.数据流分解

并行程序需要注意的几个问题：

1.同步（Synchronization）

2.通信

3.负载平衡

4.可扩展行（Scalability）

并行误差扩散算法程序实现如下：

程序主要按照论文和书中的代码实现，但是对于数据要进行一些特殊的处理

讨论帖子：讨论帖

// ERROR-diffusion.cpp : 定义控制台应用程序的入口点。
//

#include "stdafx.h"
#include "omp.h"
#include "windows.h"
#include <iostream>

#include <opencv2/core/core.hpp>    
#include <opencv2/highgui/highgui.hpp>  


#pragma comment(lib,"opencv_core2410d.lib")                  
#pragma comment(lib,"opencv_highgui2410d.lib")                  
#pragma comment(lib,"opencv_imgproc2410d.lib")     

using namespace std;  
using namespace cv;  

void error_diffusion(unsigned int width,
					 unsigned int height,
					 unsigned short **InputImage,
					 unsigned short **OutputImage
	)

{

	for (unsigned int i = 0;i < height-1; i++)
	{
		for(unsigned int j = 1;j < width-1; j++)
		{
			//计算输出像素的值
			if (InputImage[i][j]<128)
			{
				OutputImage[i][j] = 0;
			}
			else{ OutputImage[i][j] = 1;}

			//计算误差值
			int err = InputImage[i][j] - 255*OutputImage[i][j];
			//扩散误差
			int v;
			v=(int)InputImage[i  ][j+1];v+=err*7/16;if (v>255) v=255; if (v<0) v=0;InputImage[i  ][j+1]=(unsigned short)v;
			v=(int)InputImage[i+1][j-1];v+=err*3/16;if (v>255) v=255; if (v<0) v=0;InputImage[i+1][j-1]=(unsigned short)v;
			v=(int)InputImage[i+1][j  ];v+=err*5/16;if (v>255) v=255; if (v<0) v=0;InputImage[i+1][j  ]=(unsigned short)v;
			v=(int)InputImage[i+1][j+1];v+=err*1/16;if (v>255) v=255; if (v<0) v=0;InputImage[i+1][j+1]=(unsigned short)v;
		}
	}
}

int row = 288;
//int col = width;

void error_diffusion_omp(unsigned int width,
	unsigned int height,
	unsigned short **InputImage,
	unsigned short **OutputImage
	)

{
	int cpu_num = omp_get_num_procs();//cpu数
	 
	int col = width;
#pragma omp parallel private(row , col)//并行域
	{
		int thread_id = omp_get_num_threads();//每个线程的线程号
		Sleep(20*thread_id);//根据线程短延迟
	
#pragma omp for
		for (int i = 0; i<(height/cpu_num);i++)
		{
			 row  = row*cpu_num + thread_id;

			for ( col = 0;col<width;col++)
			{
				//计算输出像素的值
				if (InputImage[i][col]<128)
				{
					OutputImage[i][col] = 0;
				}
				else{ OutputImage[i][col] = 1;}

				//计算误差值
				int err = InputImage[i][col] - 255*OutputImage[i][col];
				//扩散误差
				InputImage[i][col+1] += err * 7/16;
				InputImage[i+1][col-1] += err * 7/16;
				InputImage[i+1][col] += err * 7/16;
				InputImage[i+1][col+1] += err * 7/16;

			}
		}

	

	}
}



int _tmain(int argc, _TCHAR* argv[])
{
	string str_name = "result.pgm";  
	

	Mat image_src = imread(str_name,0);
	imshow("original image",image_src);

	unsigned short **InputImage = new unsigned short *[image_src.rows];
	for (int i = 0;i<image_src.rows;i++)
	{
		InputImage[i] = new unsigned short [image_src.cols];
	}
	unsigned short **OutputImage = new unsigned short *[image_src.rows];
	for (int i = 0;i<image_src.rows;i++)
	{
		OutputImage[i] = new unsigned short [image_src.cols];
	}

	cout<<image_src.rows;
	cout<<image_src.cols;
	Mat image_dst(image_src.rows,image_src.cols,CV_8U);

	for(int y = 0;y < image_src.rows;y++)  
	{  
		uchar *ptr= image_src.ptr<uchar>(y); 
		
		for(int x = 0;x < image_src.cols;x++)  
		{  
			
			InputImage[y][x] = ptr[x];
			OutputImage[y][x] = 0;
			
		}  
	} 


	error_diffusion(image_src.cols,image_src.rows,InputImage,OutputImage);
	//error_diffusion_omp(image_src.cols,image_src.rows,InputImage,OutputImage);

	for(int y = 0;y < image_src.rows;y++)  
	{  
		uchar *ptr= image_dst.ptr<uchar>(y); 

		for(int x = 0;x < image_src.cols;x++)  
		{  
			if (OutputImage[y][x]==1)
			{
				ptr[x] = 255;
			}
			else
				ptr[x] = 0;

		}  
	} 
	imshow("error diffusion",image_dst);
	imwrite("result.bmp",image_dst);

	waitKey(0);


	return 0;
}

效果图：

参考文献：

SHAMEEM AKHTER(孟). 多核程序设计技术--通过软件多线程提升性能[M]. 电子工业, 2007.

张春柳, 李嘉, 熊琭. 基于OpenMP实现的误差扩散算法[J]. 软件产业与工程, 2015(1):44-48.