预处理矩阵（共7篇）
预处理矩阵 篇1
1 图象的矩阵表示
在电脑显示器中, 平面图象其实就是由矩阵来表示的, 像素与分辨率是平面图像处理中两个最为重要的基本概念。像素是组成图像的最基本单位, 也就是矩阵的元素。像素是我们描述图像大小的依据, 也就是矩阵的大小。分辨率是指单位面积里所包含的像素的数目, 基本上分辨率会随着单位面积含有像素值的提高而有所提高, 输出的画面质量也会相应地提高, 相反, 在分辨率极低的情况下, 画面就会呈现锯齿状。位图使用我们称为像素的一格一格的小点来描述图像。计算机屏幕其实就是一张包含大量像素点的网格矩阵。
2 基于矩阵运算的图象变换
当我们进行位图编辑时, 其实是在一点一点的定义图像矩阵中的所有像素点元素的数值, 图象的变化对应于矩阵的运算, 比如说图象的相减对应于矩阵的相减, 图象的亮度变化对应于矩阵去乘一个系数, 将位图对应的图象矩阵去乘常数0.5, 图象的形状没有变化, 只是亮度减弱了。
3 图象的小波变换
图象的小波变换也是图象的一种变换, 也通过图象矩阵的运算来实现。
假设f x, 乙y乙是表示一个二维图像, ψx, 乙y乙是一个小波函数 (其中, x, y分别表示行和列, f x, 乙y乙的值表示矩阵元素的值) 。它满足小波变换的许可条件, 具有某种程度的正则性并对于特定的图像和特定的应用要求具有合适的波形。如前所述, f x, 乙y乙的小波变换为:
我们通过对矩阵的上述运算对图象进行了三层小波分解, 进而提取其各层的低频系数的图象。一幅图象的信息主要集中于低频部分, 这一特性在图象压缩中将起到重要作用。
4 基于奇异值分解的人脸识别技术
近年来人脸识别技术经过不断地研究获得了成功, 并在实际应用领域中得到尝试。科学技术的飞速发展, 理论不断地完善, 并在实验成功后, 将进入到实际运用阶段。从原理的角度而言, 人脸识别技术主要是从两个方向进行研究, 即从宏观理论的角度, 基于整体的研究方法。其在研究过程中, 考虑更多的是模式的整体属性, 其中包括特征脸、弹性图匹配 (Elastic Graph Matching) 以及奇异值分解的方法。此外, 还将与模型和网络相对应的方法建立起来, 比如隐马尔可夫模型 (Hidden Markov Method) 、神经网络等等。与整体研究相对应的是基于特征分析的方法。这种方法是将人脸作为基准点, 并计算出相对比率。为了制定出识别特征向量, 还要将描述人脸的各项参数, 如脸部特征的形状参数和脸部的类别参数等等的识别特征进行处理。人脸识别中, 将奇异值分解 (Singular Value Decomposition) 纳入其中, 并以此为依据发展出了多种不同的识别方法。然而, 研究表明, 即便是图像不同, 也会形成相同的奇异值。奇异值不足以作为识别的依据, 在于正交矩阵由特征向量素构成, 而图像的主要信息正包含在其中。只有对于图像进行重构, 才能够形成识别效应, 那么, 要实现奇异值在人脸识别中的组合价值, 就要将奇异值与特征向量综合利用。从实验的结果表明, 要获得很高的识别率, 运用这种方法不需要将过多的奇异值和训练样本运用于其中。明显地降低了计算复杂度。在进行图像奇异值分解的时候, 运用人脸数字, 可以在重构图像的时候, 运用奇异值和与其对应的特征向量, 将原有图像中的面部表情、姿势以及光照等等所形成的噪音剔除, 以避免高频信息出现。在进行重构图像识别的时候, 将其作为模板来执行。保留奇异值以及特征向量的数目越高, 其识别率就会越高。但是随着数目的越来越多, 其对于识别的贡献度将会相应地降低。为了提高识别率, 可以增加训练样本。当训练样本多, 识别率就会降低, 相反, 训练样本较少的时候, 识别率就会相应地升高。当然也存在奇异值和训练样本较少而识别率相当高的情况。采用高于PCA的方法, 提高了识别率的同时, 还可以将工作量减低。在进行识别的时候, 要将第一个奇异值和第二个奇异值以及与其相对应向量值保留下来, 此时, 单幅图像的识别率甚至于可以超过60%。这种方法之有效性甚至于可以超越任何的方法。当然了, 在识别奇异值和特征向量的时候, 针对于具体事项还有待于进一步研究。
5 基于矩阵分解的图像融合技术
图像融合是近年来图像工程中一个新兴的研究领域, 其主要思想就是采用一定的算法, 把两个或多个不同图像融合成一个新的图像, 从而使融合的图像具有更高的可信度、较少的模糊、更好的可理解性, 它是多传感器信息融合的一个重要分支, 作为一种有效的信息融合技术, 已经广泛应用于机器视觉、医疗诊断、军事、遥感等领域。
图像融合方法中具有代表性的主要是加权平均法和小波变换法。加权平均法是最简单的图像融合方法, 即将参与融合的源图像对应像素的灰度值进行加权平均, 生成一幅新的图像, 这种方法最为简单, 但其同时效果也最差, 适用面比较窄。小波变换法具有良好的空域和频域的局域性, 可以较好地保留多聚焦图像的高频信息, 但仍在一定程度上会丢失原始图像中的一些特征信息 (如边缘信息) 。利用非负矩阵分解 (Non-negative Matrix Factorization, NMF) 的思想来达到对图像中特征信息的保留的目的, 与小波变换的图像融合方法相比较, 图像的特征信息损失较少, 其融合结果也更接近于理想的实际图像。将NMF算法思想用于图像融合, 该算法能够自动地进行融合, 基本不需要人工干预, 而且适用面较广, 能够较多地保留图像中的边缘特征等信息, 通过可见光图像与毫米波图像融合实验以及曝光过度图像和曝光不足图像融合实验的结果表明, 基于NMF算法能够得到比小波变换方法更好的融合效果。
6 结论
矩阵理论是数学理论中的一个重要分支, 它在许多方面有着诸多应用, 本文通过分析举例浅显的谈论了矩阵运算在图象处理中的一些应用, 希望能对有志于深入研究图像处理的读者有所启示。
参考文献
[1]Aggarwal JK Multisensor Fusionfor Computer Vision[M].Berlin Heidellberg:Springer Verlag, 1993.
[2]Varshney PK1 Multisensor data fusion[J].Electronics&Communication Engineering Journal, 1997.
预处理矩阵 篇2
在这个信息爆炸的时代,对于海量数据处理的解决方案以及其他诸多相关问题,在云计算概念出现之后,人们有了新的思路和希望。
1 处理方案
1.1 数据输入
由于Map过程需要将整数对(air,brj)作为输入,而矩阵以文件的形式存储在HDFS中,所以需要对输入文件进行处理,为Map过程提供输入参数。这个预处理的过程(如图一):首先从HDFS中读取文件M0和M1,并把数据分别读入左矩阵m[0]和右矩阵m[1]中。然后将左矩阵的行和右矩阵的列分割成多个Split块。矩阵分割完毕后,再将每个Split块形成多个键值对<(row,col),(air,brj)>并传送给map任务。
如图二中所示的两个矩阵,经分割后形成的某组键值,如图三所示。
1.2 Map阶段
map接收输入的键值对<(row,col),(air,brj)>之后,将value中的两个数相乘,相乘结果如图四所示。然后对map的操作结果进行合并处理,即将中间结果中有相同key的
合并处理结束后,我们需要根据map的输出
1.3 Reduce阶段
此阶段将处理map输出的数据对。我们根据key值,把具有相同的key的数据对的value值放入一个列表中,然后形成新的数据对
1.4 数据输出
在这个阶段,我们需要把reduce方法处理并输出的value值按照积矩阵的格式写入文件中。
2 程序实现
2.1 程序中的基本结构类
(1)matrix类用于存储矩阵,其类图如图八所示。当Input文件中的矩阵数据被读出以后,便经过转换成为整型,然后按照行列再依次存入到为整型的二维数组matrix中,供下面的初始化过程使用。
(2)Int Pair类实现Writable Comparable接口用于存储整数对的,其类图如图九所示。这些整数对主要是为表示后面对两个matrix进行分割后形成的各个块中的数据矩阵元素的位置和生成一个积矩阵元素的多个数据对。
2.2 初始化类
数据输入主要通过martix Input Format、matrixInput Split、martix Record Reader这些类来实现。在这里,我们通过构造matrix Input Format类取出文件中的数据,并把数据按照格式要求进行类型转化,然后存到数组matrix[0]和matirx[1]中。之后对两个矩阵按照matrix[0]的行和matirx[1]的列分割成多块,其中每一块为matrix Input Split类型。该类型继承Input Split,每个matrix Input Split包括用来生成一个积矩阵元素的多个数据对。这些数据对存储了矩阵元素的位置和生成一个积矩阵元素的多个数据对中的一个,它们都为Int Pair类型。而matrix Record Reader负责将每个matrix Input Split块分解成多个
2.3 map阶段的类
map阶段主要通过Matrix Mapper、Matrix Combiner、First Partioner这些类来实现。Matrix Mapper类继承Mapper,它主要是实现其map方法,map主要用来处理输入
2.4 reduce阶段的类
reduce阶段主要通过Matrix Reducer类实现。Matrix Reducer类主要是实现reduce方法,框架将map输出的中间结果根据相同的key2组合成
3 结束语
云计算的关键技术Map Reduce是一种主要应用在大型计算机集群上,并行处理海量数据的框架模型。用户通过指定一个map函数将输入数据转变为一系列中间键值对,再通过自定义的reduce函数将具有相同键的值聚集起来,输出结果。在很多处理海量数据的领域,都用这种模型来表示。
参考文献
[1]Anthony T.Velte著.周庆辉,等译.云计算实践指南[M].北京:机械工业出版社,2010.
预处理矩阵 篇3
快速傅立叶变换(FFT)/快速反傅立叶变换(IFFT)是时域和频域转换的基本运算,广泛应用于现代数字信号处理的各个领域。自从Winograd在1976年提出的Winograd快速傅立叶变换算法(WFTA)以来,FFT/IFFT的运算量已得到了很大幅度的减少[1],各种算法也层出不穷。1984年,杜哈梅尔(PDohamel)和霍尔曼(H.Hollmann)提出了更有效的分裂基的快速算法。这些算法经过改进形成一整套算法,即现在的FFT算法。为了实现高速FFT,通常在算法、系统架构和运算单元这三个层面上进行优化。总结目前的FFT的算法,大致可以分为基-2、基-4、基-8、高基算法及混合算法等[2];根据数据通道数与架构不同,大致可以分为RβMDC[3]、RβSDF[4]、R2i SDF[5]、RβSDC[6]、MRMDC[7]等结构。
基β单通道延迟反馈(RβSDF)结构优点是对寄存器的需求是最少的,但它主要应用于单输入单输出OFDM系统,处理数据周期比较长,在多通道应用中没有优势。R2i SDF结构与RβSDF结构基本上相同,但它的加法器的数量会更少。基-β单通道延迟切换(RβSDC)使用一个改良的基-β运算结构,用可编程的1/β个基-β蝶形运算达到较高的乘法器利用率[6]。其对寄存器的需求介于RβSDF与RβMDC之间,并行性上与RβSDF存在相同的问题。基-β多通道延迟切换(RβMDC)结构将输入序列分解成β个并行数据流向前流动,处理数据的并行度高,吞吐率高,常应用于需实时处理的领域。但对寄存器与复数乘法器等需求大,硬件成本高。混合基多通道延迟切换(MRMDC)结构可以看成是RβMDC与RβSDF结构的组合,结构多变且复杂,但可以处理N不是β幂次方的情况,应用比较灵活。因此,在高速应用中,经常会用MDC结构,UWB系统就是一个很好的例子。本文针对MDC结构中需大量延迟单元问题,提出了用存储器矩阵转置方式实现数据的交换,减少了MDC结构中延迟单元的使用,减少了硬件资源的开支,在高基应用中有明显的优势。
文章第1节介绍算法推导与RβMDC结构流水线FFT处理器。第2节介绍存储器矩阵的结构,组成与数据交换过程,以及用存储器矩阵实现流水线FFT处理器的RβMATC结构。第3节介绍ATC结构在MRMDC结构中的应用并给出仿真结果。第4节得出结论。
1 RβMDC结构流水线FFT处理器
RβMDC结构是流水线FFT处理器最常用的结构,这一节对它的算法推导与结构做个介绍。
1.1 算法推导
设输入序列x(n)点数N=βL,L为整数。则N点DFT定义为:
其中, WNn*k是旋转因子,0≤n、k<N-1,n、k用β进制数表示为:
将(2)、(3)式入(1)式可得:
按频率抽取法则进行抽取,且因N=βL,则
将(5)式代入(4)式并整理可得递推关系:
从(6)式中可以看出,Xm( k0, k1, .. .km-1,nL-m -1,...,n1,n0)
Xm-1是( k0, k1,. . .km-2, nL-m,. . ., n1, n0) 乘上一个旋转因子后的β点DFT,其输入变量为nL-m,输出变量为km-1。所以结构设计时,可以用级连的方式实现。
1.2 结构
图1为N=64点的R2MDC、R4MDC与R8MDC结构图,图中略去了旋转因子存储ROM与控制逻辑。RβMDC结构将输入序列分解成β个并行数据流向前流动,数据元素间需保持正确的“间距”进入蝶形运算。正确的“间距”是通过延迟单元实现的,图中已经标出了各级所需延迟单元的个数。因延迟单元个数差异很大,无法用存储器的宏单元实现,只能用寄存器实现。RβMDC结构要使用(β-1)(logβN-1)个乘法器、logβN个基β蝶形运算器与(β+1) N /2-β个寄存器。从中可以看出,β值越大,并行度越大,所需要的基β蝶形运算器会越少,所需时钟数越小,但需要的延迟单元数越多。
2 ATC在基-β流水线FFT处理器中的应用
基β多通道延迟切换(RβMDC)结构中要用很多延迟单元,且换向器的结构比较复杂,用ATC直接代替RβMDC结构中的延迟单元可以很好地解决这些问题,这种结构我们称之为基-β矩阵转置交换(Multi-Path Array Transposed Commutator--RβMATC)。下面以64点基-4 FFT处理器为例先说明一下存储器矩阵转置结构,这一结构我们称之为R4MATC结构。
2.1 存储器矩阵结构
R4MATC需用到4*4存储器矩阵,如图2所示。用Bankij表示i行与j列交叉点对应的存储器模块。Bank内单元数比较多时候可以由随机存储器(RAM)实现,较少时用寄存器实现。
Memory array(MA)的4行或4列可以同时写入数据或输出数据。数据可以从行方向写入,然后从列方向读出;也可以从列方向写入,再从行方向性读出。通过这行列的切换,读出数据序列与FFT处理器抽取顺序是相吻合的。所以存储器矩阵的读与写之间相当于矩阵行列之间的转置关系,所以称之为矩阵转置交换(ATC)。
存储器矩阵的读写控制很简单,输入64点数据以输入先后标记为x0~x63,第一帧数据先从行方向写入MA,数据在MA中的存储位置如图3所示。行写入时,用一个计数器(命名为WR_count)产生Bank与Bank内单元地址,按地址的自然顺序写入。列读出时,用另一个计数器(命名为RD_count)产生读地址,按地址的倒位序读出。列读操作开始后,就允许第二帧数据的输入。第二帧数据从列方向写入,写入顺序与第一帧数据的读出顺序相同。第二帧数据全部到位后,可以启动从行方向读出过程,读写操作过程相同,只是方向不同。第三帧数据的操作过程与第一帧相同,行列的读写顺序按此切换即可。
2.2 RβMATC结构
64点采用R4MATC结构的FFT处理器框图如图4所示。图中省去了旋转因子存储ROM与读写控制等部分。用到三个MA都采用4*4矩阵,MA1与MA2矩阵的每个Bank需有4个存储单元,而MA3的每个B a n k只有一个存储单元, 所以6 4 点F F T采用R4MATC结构共需64+64+16=144个单元,比R4MDC结构少了12个单元。
对MA1的操作过程如图3所示,对MA2与MA3的操作过程与此类似,不同之处在于MA2与MA3的写入也是并行的,且MA3的每个Bank只有一个存储单元。从以上过程可以看出,R4MATC结构与R4MDC结构并无太大差别,数据在蝶形中运算过程相同,数据的流动顺序相同。
图4 用于流水线处理器的R4MATC结构(N=64)
推广到一般情况,采用RβMATC结构需采用β*β的存储器矩阵。N=βL,则需要L级基-β蝶形运算单元与L个存储器矩阵。前2个存储器矩阵每个Bank需要βL-2个存储单元。如L>2,则从第三级开始,后一级B a n k内单元个数是前一级的1/β ,直到最后一级Bank内单元个数只有一个。可以看出,N一定时,β值越大,所需蝶形运算级数越少,所需存储器单元数越少。表1、表2分别给出了基-2、基-4、基-8在不同点数下RβMATC与RβMDC所需存储单元的对比数据。从中可以看出,基-2因通道数少,MATC反而使用便多的存储单元。但除此之外,β值越大,MATC结构优势越大。另外N值越大,MATC结构优势越大。
表3给出了R8MDC与R8MATC两种结构在不同点数下需要存储器的个数, 从中可以看出,R8MATC结构对存储器的需求较少,相当于R8MDC结构的45.7%左右。
3 ATC结构在MRMDC结构中的应用
从以上的介绍中可以了解到,RβMATC结构应用存在局限性,它并不适用于N不是β的幂次方的情况。在这一节中,我们介绍在MRMDC结构中如何利用ATC结构简化处理器的结构与减少延迟单元的使用量,并以UWB-OFDM系统中128点FFT处理器为例进行说明。
3.1 算法推导
UWB-OFDM系统中使用128点FFT处理器,除基-2外无法用RβMDC结构,而采用基-2算法,处理器需要7级蝶形运算,需要较多的复数乘法器,也会有较大的延时。因此,大家都采用MRMDC结构实现,基-(8+16)算法也是其中的一种结构,算法推导过程如下。
对输入序列X(n),128点离散傅立叶变换(DFT)的公式为:
其中,Wn*k128是旋转因子,n=0…127,k=0…127。
算法采用8*16方式分解,即令n=16*n1+ n0,k=8*k1+k0,其中n1=0…7,n0=0…15,k1=0…15,k0=0…7。
则
第一级是一个n1=0基-8运算,可以用一个并行基-8蝶形运算单元实现,第二级是一个16点FFT,可以用n0=0一个串行基-16蝶形运算单元实现。
3.2 处理器结构
FFT处理器由一个存储器矩阵、一级并行基-8蝶形运算单元与一级R16SDF等组成,如图5所示。存储器矩阵是一个8*8的存储器矩阵,矩阵的每一个点为一个Bank模块,每个Bank由两个存储单元组成,共128个单元。基-8蝶形运算单元采用并行流水结构,其输出数据乘以旋转因子后可以直接传送给基-16蝶算单元。基-16蝶算单元采用8个并行的R16SDF结构。而R16SDF采用了4级的R2SDF实现,其结构如图7所示。因输入数据的顺序问题,SDF结构的延迟单元分布与传统的SDF结构不同。整个处理器可以实现输入数据是顺序的,输出也是顺序的,输出数据可以不用存储器存储。从结构框图中可以计算出延迟单元的使用量为128+15*8=248,完整的复数乘法器为7个,另加8个较简单一点的乘法器,总体硬件开支较少。
3.3 存储器矩阵结构
存储器矩阵的结构要根据N的值与数据抽取方式而改变,在所设计的结构中,第一级基8是按16抽取的,且为了FFT处理器是按自然序输出,对存储矩阵的接法进行了调整,如图7所示。
存储器矩阵采用8*8矩阵,i行与j列交叉点对应的存储模块记为Bankij,每个Bank有二个存储单元,都由双端口寄存器组成。存储器矩阵(Memory array)可以从行方向写入数据,然后从列方向读出数据;也可以从列方向写数据,再从行方向性读出数据。
输入数据按时间顺序记为a0~a127,可以单通道或8个通道并行输入。8个通道并行输入时,输入数据线的连接关系要调整,如图8所示。两个方向都是1与4交叉,3与6相交叉。写入时按Bank内地址的自然顺序写入数据。因此,从行方向输入为例时,数据在存储矩阵中的存储位置如图8所示。从列方向读出时,也按Bank内地址的自然顺序输出,第一个时钟并行读出的8个数据刚好是a0、a16、a32、a48、a64、a80、a96、a112,与FFT的抽取关系相同。以j=0列为例,16个时钟读取的数据顺序是a0、a8、a4、a12、a2、a10、a6、a14、a1、a9、a5、a13、a3、a11、a7、a15。对于所设计的处理器结构,可以通过调整R16SDF延迟单元的个数来适应这种顺序,各级延迟单元的个数在图6中已经标明。
从以上分析可以看出,所设计的处理器结构可以采用矩阵转置方式进行数据的切换,并采用并串混合结构实现,硬件开支较少。在数据延迟方面,如并行基-8运算单元内部按加一级流水计算,则所设计的结构延迟了16个时钟周期。
3.4 Modelsim仿真与DC验证
我们用Verilog HDL语言进行RTL实现,用寄存器实现存储器矩阵。RTL代码仿真结果与Matlab运算结果进行了对比,验证了算法的正确性。并将RTL代码基于TSMC 90 nm标准工艺库,采用Design Compiler进行综合,得到门级评估面积与功耗,其结果如表4所示。芯片的内核面积为0.6056mm2。在52MHz频率下工作就可以满足UWB-OFDM系统的要求,内核功耗为5.9m W。在1.2V,25℃条件下,最大工作时钟可达200MHz。表4同时也列出了二个有关UWB-OFDM FFT处理器文献的数据。与文献[10]对比,本文推荐处理器在面积与功耗上都有较大优势。与文献[8]对比,本文推荐处理器由于在计算过程采用12位数据宽度,在面积略大于它,但在功耗方面会略优于它。
4 结论
本文介绍了用存储器矩阵转置切换实现FFT处理器的数据存储与抽取过程,存储器矩阵转置输出的数据顺序恰好符合FFT处理器抽取数据要求,并能简化控制逻辑与减少延时单元的使用量,适用于各种RβMDC处理器结构。本文还以UWB-OFDM系统的128点FFT为例,介绍了一种N非β幂次方的情况下存储器矩阵转置切换的实现方法。介绍的128点FFT处理器采用基-(8+16)算法,硬件上采用存储器矩阵+并行蝶8处理器+R16SDF结构实现。同时对该结构用Verilog HDL语言编程实现,存储器矩阵用寄存器实现,并基于TSMC 90 nm标准工艺库,采用Design Compiler进行综合,综合结果表明所设计的FFT处理器在面积和功耗上具有一定优势,数据吞吐率大,能够满足MB-OFDM UWB标准中FFT处理器的要求。
参考文献
[1]R.G.Lyons,Understanding Digital Signal Processing,Publishing House of Electronics Industry,China,pp.87-122,2010.
[2]PQ.Cheng.Digital Signal Processing.Publishing House of Tsinghua University,China,pp.143-194,2001.
[3]K.-J.Yang,S.-H.Tsai,G.C.H.Chuang,MDC FFT/IFFT Processor With Variable Length for MIMO-OFDM Systems.Very Large Scale Integration(VLSI)Systems,IEEE Transactions on.21(4):720-731,2013.
[4]E.H.Wold,A.M.Despain.Pipeline and parallelpipeline FFT processors for VLSI implementation.IEEE Trans.Comput.,C-33(5):414-426,May 1984.
[5]C.Vennila,G.Lakshminarayanan,Seok-Bum Ko.Dynamic Partial Reconfigurable FFT for OFDM Based Communication Systems.Circuits Syst Signal Process 31:1049-1066,2012.
[6]G.Bi and E.V.Jones.A pipelined FFT processor for word-sequential data.IEEE Trans.Acoust.,Speech,Signal Processing,37(12):1982-1985,Dec.1989.
[7]Lee,Jea Hack;Kim,Eun Ji;Myung Hoon Sunwoo.Low complexity FFT-IFFT processor for highspeed OFDM system using efficient multiplier scheduling.Circuits and Systems(ISCAS),2012IEEE International Symposium on.pp.1520-1523.May 2012.
[8]Song Nien Tang,Jni Wei Tsai,Tsin Yuan Chang.A2.4-GS/s FFT Processor for OFDM-Based WPAN Applications[J].IEEE Journal of Circuits and Systems,2010,57(6):451-455.
[9]高亚军,主编.基于FPGA的数学信号处理.北京:电子工业出版社,2012.02.
预处理矩阵 篇4
矩阵乘法是很多应用问题的核心计算方法, 是线性代数中最常用的算法之一, 被广泛地应用到物理、生物、化学等应用领域。伴随着计算机硬件以及软件的快速发展, 多核心的CPU正在快速普及到各个部门, 当前硬件设计和软件开发的热门问题是计算机上的并行计算。并行计算的出现对传统的串行计算是一个巨大的改变, 现在multi-core (多核) 以及multi-threading (多线程) 处理器的发展非常迅速, 迫切需要解决矩阵乘法在多核多线程处理器上并行计算问题。
当前, 国内外对多核多线程的矩阵乘法的研究主要集中在以下几个方面:
矩阵乘法的分层。在国内外的文献中, 对于矩阵乘法的优化主要从矩阵乘法的分层来进行优化。文献[1]~[4]都提到了BLAS分层的方法。文献[1]总结了用于建立高性能BLAS代码的理论和实际的方法, 将矩阵乘运算划分成3层。第一层是向量运算, 第二层是矩阵-向量乘法, 第三层是矩阵-矩阵乘法。文献[2]提出了一种简单并且高效的方法, 这种方法可以把基于Cache这种结构的矩阵乘的高性能运算转变成为采用矩阵乘运算的通用的结构。通过把GEMM划分成GEPP这种形式来实现对矩阵乘法的优化。文献[3]和文献[4]提出了Packing的概念, 分别将矩阵A、B存放进连续的内存空间中, 这样可以减少TLB MISS和额外的开销。
在文献[5]中, 对矩阵乘法划分方法进行了总结, 对GEMM作了6中方案的划分, 通过分析得出:对矩阵A和B同时进行划分的时候矩阵乘法得到的性能是最优的。文献[6]对通用的矩阵乘法C=AB+C (A为m*k矩阵, B为k*n矩阵, C为m*n矩阵) 作了分析, 考虑一种特殊的情况:k很小, m, n相对很大, 通过分析得出这是一种重要的矩阵乘运算的特殊情况。
本文基于Open SPARC T2处理器平台, 采用分块算法、循环展开等优化途径对其进行矩阵乘法的优化, 使其性能得到了极大地提高, 本文所采用的优化技术也同样适用于其他的一些高性能处理器结构。
2、Open SPARC T2处理器结构
Open SPARC T2处理器包含8个核, 每个核有8个线程。图1是该处理器的结构, 针对软件优化, Open SPARC T2体系结构主要具有以下特点:
(1) 拥有独立的功能部件, 包括指令读取部件 (IFU) 、执行部件 (EXU) 、存储部件 (LSU) 、浮点功能部件 (FGU) 、陷阱逻辑部件 (TLU) 和内存管理部件 (MMU) 。其中, IFU包括取指令单元、指令选择单元和指令译码单元;EXU执行除了整数的乘法和除法外所有的整数运算和逻辑运算;FGU包括三条流水线, 分别是:浮点执行流水线 (FPX) 、图形运算流水线 (FGX) 和执行浮点除法和开方的流水线 (FPD) 。 (2) 每个核包含2条整数流水线、1条浮点流水线和1条存储流水线。8个线程共享浮点流水线和存储流水线; (3) 每个核内有独立的一级数据cache和一级指令cache。L1I Cache大小为16KB, 8路组相联块大小32字节。L1D Cache大小为8KB, 4路组相联, 块大小32字节; (4) 指令TLB为64项全相联, 数据TLB为128项全相联。8个线程共享L1I、L1D和TLB, 通过TLB中的自动释放机制使多个线程更新TLB。 (如图1)
3、Open SPARC T2处理器上的矩阵乘法优化
下面针对Open SPARC T2处理器来讨论如何对矩阵乘法进行优化。
3.1 分块算法
前面介绍了Open SPSRCA T2处理器的体系结构, 8个核共64个线程共享L2 Cache, 为了能够匹配运算器的计算速度, 需要首先把内存中的A, B矩阵的数据放入Cache中, 寄存器直接从Cache中取数进行计算, 在实际应用中, 矩阵A, B都是数据量非常庞大的矩阵, 而处理器整个L2的容量只有4MB, 不可能直接将整个矩阵全放入Cache, 所以必须对矩阵进行分块以便让它们能够放入Cache中。分块的大小也是有讲究的, 如果过大, 会造成Cache容量不足;如果分块过小, 处理器将内存的矩阵块调入Cache的次数会增多, 这无疑会加大处理器的负担。
按照BLAS的分层模型, 矩阵乘法可以划分为3层, 第一层是GEMM, 第二层是GEPP、GEMP和GEPM, 第三层是GEBP、GEPB和GEPDOT, 而在第二层的划分中, GEPP优于另外两种, 因为GEPP是对矩阵A、B都进行了分块, 因此, 本文所采用的分块也是基于GEPP这种划分。这三层划分方法分别如图2 a、b、c所示。
在GEMM中, 对矩阵A按列划分, 矩阵B按行划分, 得到图b的GEPP形式, 在GEPP划分中, 设A、B、C三个分块矩阵的大小分别为mc*kc, kc*n, mc*n, 把B和C按列划分, A按行划分, 每nr列B和C分别记作Bj、Ck, 如图3所示, 每次取出mc*kc的A分块 (Ail) , kc*nr的B分块 (Blj) 放入一个进程来计算Cij=Ail*Blj, 划分方法如图3所示。
3.2 优化技术
对于矩阵乘的优化, 循环展开也是一种很重要的优化方法, 通过减少循环分支次数来达到减少循环开销的目的;而且可以尽可能的重用数据来减少处理器访存的次数;可以拉长循环代码以便进行指令调度;使得更多的指令在不同的功能部件上并行执行来提高指令的并行性。在前面已经介绍过了矩阵乘的串行优化算法, 矩阵乘法算法实际上是做了3层循环, 具体的实现过程是每次读取矩阵A的一行和矩阵B的一列来求出对应的C矩阵中的元素。每次读取A中的一个数据和B中的一个数据, 这两个数据来进行乘法和加法运算, 如果矩阵规模很大的话, 要多次从内存读取数据, 这样做的效率很低。因此需要对外层循环进行展开, 让每次读取的数据尽可能多的参与到计算中。
3.3 并行算法
在进行矩阵乘法的运算时候, 考虑到在实际的工作中, 矩阵都是相当大的, 这就需要我们对矩阵进行分块, 每个线程执行块A和块B的乘加运算。Open SPSRCA T2处理器一共有64个线程, 这样就可以同时在64个线程中进行运算, 可以大大的提高处理器的运算效率。因此, 在实际编程过程中, 可以采用OPENMP多线程技术来实现矩阵乘法的并行算法。从理论上来讲, 线程数越多, 这个时候处理器可以进行更多的并行计算, 从而得到更大的运算效率。
在2.1中, 我们通过矩阵的分块算法分析, 研究了在一个线程中矩阵乘法是如何进行计算的, 那么下面看看在64个线程中矩阵乘法是如何进行运算的。
在矩阵乘法的执行过程中进行Packing, 可以减少额外的开销和TLB Miss, Packing就是把3.1中提到的分块Ail和Blj存放在连续的空间中, 如图4所示, 假设A、B、C均为4阶矩阵, 一般情况下, C21=A21B12+A22B22+A23B32+A24B42, 在这里, Ail就是mc*kc的A分块, Blj就是kc*nr的B分块。在并行算法中, 把分块Ail放入Packing过的中, 把分块Blj放入Packing过的中, 对于进程P (r, c) , 在它里面完成局部计算, 这样就可以实现在64个线程中的矩阵乘法并行算法。比如说对于进程P (1, 1) , 它所完成的是C=C+A1lBl1 (l=0, 1, 2, 3) 。具体算法见图4:
4、实验结果和分析
基于Open MP并行环境, 用C语言实现了前面部分所描述的算法, 并分别在SPARC-T2处理器上, 对算法进行了评估。所做测试的矩阵乘法阶数为8000阶, 也即M=N=K=8000, 分别测试它们在不同的线程中, CPU处理浮点运算的能力。图6是实验中SPARC-T2处理器, 分块大小采用256×256×128, 分块在不同数量的线程中CPU处理浮点运算的性能对比。 (如图6)
矩阵分块大小对矩阵乘法的性能会产生影响, 对比单线程, 并行矩阵乘算法在8核64线程的加速比达到21.9%。图7是不同的矩阵分块大小对性能影响的对比, 可以得出, 当分块大小为256×256×128时, 性能最好。 (如图7)
5、结论及未来工作
本文针对OPENSPARC T2处理器做了优化, 从分块和循环展开对矩阵乘法做了分析, 还分别从矩阵乘法的分块大小对性能的影响和矩阵分块位于不同线程数量时候对性能的影响这两个方面做了测试, 可以得到的结果是, 处理器的线程越多, 处理器的运算处理能力更强。对比单线程, 并行矩阵乘算法在8核64线程的加速比达到21.9%。测试的矩阵规模为8000阶, 当矩阵分块大小为256×256×128时, 矩阵乘法的性能最优。经过上述优化后, Open SPARC T2处理器的性能发挥峰值性能的53.9%。
在未来的工作中, 用汇编语言来实现核心模块的功能可以使矩阵乘法获得更好的性能。
摘要：矩阵乘法是很多应用问题的核心计算模块, 在OpenSPARC T2处理器平台上, 对矩阵乘法算法进行了设计优化, 针对矩阵乘法访存特点, 利用处理器本身8核64线程的特征, 基于Open MP并行编程模型设计了矩阵乘多线程并行算法, 并对访存和块大小进行了优化, 采用C语言编程, 对比单核单线程, 并行矩阵乘算法在8核64线程的加速比达到21.9%, 发挥峰值性能的53.9%。
关键词：矩阵乘法,多线程处理器,Cache
参考文献
[1]Goto K.van de Geijn R A.Anatomy of high-performance ma-trix multiplication[J].ACM Transactions on Mathematical Software, 2008, 34 (3) :Article 12 (1-25) .
[2]Goto K.van de Geijn R A.High-Performance Implementation ofthe Level-3 BLAS[J].ACM Transactions on Mathematical Software, Vol.V (3) :Article 12 (1-18) .
[3]Gunnels J A, Henry G M, van de Geijn R A.A family of high-performance matrix multiplication algorithms[C]//Proceedingsof the International Conference on Computational Science一Part I.London, UK:Springer, 2001:51-60.
[4]Dongarra, J.J., Du Croz, J., Hammarling, S., and Duff, I.1990.A set oflevel 3 basic linear algebra subprograms.ACM Trans.Math.Soft.16, 1 (March) , 1-17.
[5]Gunnels J, Lin C, Morrow G, et al.A flexible class of parallelmatrix multiplication algorithms[C]//First Merged Interna-tional Parallel Processing Symposium and Symposium on Paralleland Distributed Processing.Washington, USA:IEEE Computer Society, 1998, 12:110-116.
[6]Marker B, van Zee F G, Goto K, et al.Toward scalable matrixmultiply on multithreaded architectures[C]i Proceedings ofthe 13th International European Conference on Parallel andDistributed Computing.Rennes, France:ACM Press, 2007:748-757.
[7]Chandra R, Menon R, Dagum L, et al.Parallel Programming in Open-MP, Morgan Kaufman Publishers, Oct.2000.
[8]MalyshkinV.Parallelcomputing technologies[C]//8th inter-na-tional conference, PaCT 2005.Krasnoyarsk, Russia, 2005.Berlin;NewYork:Springer, 2005.
预处理矩阵 篇5
奈奎斯特采样定理要求对信号进行采样时采样频率要大于或者等于信号带宽的两倍。而这样得到的数据具有较大的冗余。压缩感知技术的出现，打破了采样定理的限制。压缩感知技术表明：只要信号是稀疏的或可压缩的，则数据采样可以低于奈奎斯特采样频率，而且原始信号可以由采样值精确的恢复。从而大大减少了数据的处理量。通过挖掘信号的部分先验信息，可以减少观测值方面的工作。本文中利用一部分的大系数的稀疏系数的位置设置了一个阈值矩阵，之后这个矩阵与高斯随机矩阵构成自适应观测矩阵，达到改进重构性能的目的。
2 压缩感知的原理：
标准的压缩感知的数学模型为：
其中θ为未知信号x在稀疏域Φ中的系数，y是观测矩阵Ψ对x进行先行观测后的得到的观测向量。若θ中的非零的向量为K(K<
在此选择阈值迭代算法(IHT)解决此问题。
通过挖掘信号的的部分先验信息，可以减少观测值的长度，当观测值的维数及重构的算法相同时，信号稀疏度越高，重构的精度越高。显然信号的未知支撑集合比原始信号更稀疏。利用自适应的观测矩阵，可以将稀疏域的小系数更接近于零，从而优化重构建效果。
3 观测矩阵
观测矩阵的设计性能直接影响压缩感知的质量，观测通过观测矩阵Ψ对信号x进行映射，或者压缩感知矩阵T对θ进行线性投影，若T满足受限等距特性(RIP)，则此时L0可以由θ重构，从而的到信号x的重构。自适应观测是借助信号x或者稀疏系数θ的某些先验特性，生成相应的观测矩阵，这类矩阵性能上显著由于随机观测矩阵和确定观测矩阵。文中以高斯随机观测矩阵为原始矩阵由θ的稀疏特性对矩阵的特性进行改进，使得改进后的矩阵满足受限等距特性(RIP)的同时，使得θ的小系数更接近于0。
3.1 观测矩阵的构造
设W=supp(θ);K=|W|，W1是信号的已知支撑集部分，W2是信号的未知支撑集部分，K1=|W1|，W=W1+W2，观测矩阵的自适应过程如下：
1)生成M×N的高斯随机矩阵，作为原始矩阵Ψ;
2)将Ψ最大值归一化并且按照已知和未知支撑部分分为两块;
3)选取θ中的大系数的K0个并把他们的位置记录在矩阵A中，记A为{a1,a2，…ak0};得到域阀矩阵，Φ是信号的稀疏域变换矩阵。
∆中的元素的选取规则是：若元素的的行列数不等则为零，若元素行列数相等且列数位于大系数矩阵A中，则此元素等于1，除此之外的元素赋值为极小数10-q其中q<-10，新的自适应观测矩阵是不难看出我们使用θ的部分大系数的位置信息W1，从而使优化过的自适应观测矩阵的行数缩减了K0行，从而减少了K0个观测值。在保证才数据传输量变化不大的前提下，应用后面相应的重构算法，使得应用自适应矩阵的重构效果明显好于应用高斯随机矩阵的重构效果。
3.2 压缩感知矩阵的受限等距特性（RIP)
应用自适应观测矩阵代替原有矩阵于(1)中，可知压缩感知模型变化为：
其中，,通过对压缩感知系数θ进行线性观测后可以得到观测矩阵。且当受限等距常数δk满足：若C∈(cj)j∈Γ,|Γ|≤K时有
则称压缩感知矩阵T满足RIP。
由于T、都是高斯随机矩阵，故可的RIP常数和T的RIP常数关系有有。的RIP常数的RIP常数大小关系为
4 自适应观测矩阵对压缩感知重构性能的改进
下面我们介绍使用的阈值迭代算法。IHT算法源于梯度下降算法，其采用高斯随机矩阵Ψ最为观测矩阵，传感矩阵为T=，从而观测向量值y=Tθ+e(e为观测噪声)。在运用自适应矩阵代后，传感矩阵变为，观测值向量为,IHT的迭代算法如下：
其中分别为第n+1和第n次迭代后得到的重构稀疏系数。因为可得：。IHT算法具有实现速度快，重构效果好的特点。下面对于自适应观测矩阵对IHT算法重构性能方面的改进作具体的分析。
4.1 自适应观测矩阵对IHT重构误差的改进
对于IHT重构的性能判断，其中一条是其重构误差的大小。IHT的重构误差的定义式：
其中β2k与T2k的RIP常数δ2k的关系为，应用自适应观测矩阵后IHT的重构误差为，根据参考资料[4]可知，对于由于所以的到以下结论：
可以看出，使用自适应观测矩阵后的IHT重构误差会明显小于用高斯随机观测矩阵的IHT重构误差。
4.2 仿真结果
对于信号经稀疏分解后的系数的组成往往是有限个的大系数和其余的小系数，而往往这些小系数并不为零，只是相对较小，有些甚至和零相差较大，而这种近视稀疏必然会降低压缩感知的性能。当采用自适应观测矩阵对稀疏系数θ进行观测时，可以看作是压缩感知矩阵进行线性观测。从而的稀疏性要优于θ的稀疏性。在仿真中，稀疏分解采用离散余弦变换(DCT)，重构算法采用IHT算法，信噪比取多帧的平均值。仿真软件为Matlab 2011b，分别对灰度图像rice,peppers和football进行实验。结果表明应用自适应观测矩阵的到的重构图像效果要明显优于使用高斯随机观测矩阵重构的图像。
5 小结
本文利用自适应观测矩阵，对信号进行重构。自适应观测矩阵以高斯随机观测矩阵为原始矩阵，结合信号的稀疏稀疏的部分先验信息，利用小部分的大稀疏系数的位置信息设计了一个域阀矩阵。在运算的同时减少了观测矩阵的行数，是的数据的传输代价增加较小。最后利用DCT变换和IHT算法进行仿真，证明了自适应观测矩阵的重构性能要优于高斯随机矩阵的重构性能。
关键词：自适应观测矩阵,压缩感知,阈值迭代算法
参考文献
[1]John G.Proakis,Dimitris G.Manolakis.数字信号处理[M].电子工业出版社,第四版方艳梅等译,2007年
[2]Katsuhiko Ogata.现代控制工程[M].电子工业出版社,第五版,卢伯英等,2011年
[3]U.Meyer-Baese.数字信号处理的FPGA实现[M].清华大学出版社,第三版,刘凌译,2011年版
预处理矩阵 篇6
关键词：欧式距离矩阵,矩阵填充,奇异值分解,低秩
近年来EDM重建问题得到许多学者的关注和研究, 它主要应用于机器学习, 多维尺度分析, 核磁共振分子构象等方面。根据给定的几个成对节点间的距离如何有效地重建低维几何结构的节点?这就是欧氏距离矩阵填充所要解决的问题。本文利用低维空间节点距离矩阵的低秩性, 将缺少的数据元素进行有效重建, 得到准确、结构性良好的欧式距离矩阵 (EDM) 。
1相关理论
1.1欧式距离矩阵
下面考虑EDM填充问题。由于欧式距离矩阵具有低秩性, 可通过秩最小化的凸优化问题实现矩阵重建, 其数学表达式为
1.2固定点迭代算法
本文采用固定点迭代算法 (FPCA) 进行求解, 该算法的核心思想是算子分裂技术。
结合式 (13) 和式 (18) 给出用FPC算法解决式 (9) 问题迭代式:
2数值结果
本文提出固定点迭代 (FPC) 算法重建目标矩阵DM, 为解决欧式距离矩阵填充问题提供了一个有效方案。该算法主要用来实现秩最小化矩阵填充问题, 通过选择适当参数, 运行程序用Matlab语言编写, 重建效果显著。
评估EDM重建准确率的一个重要参数为相对误差, 而采用FPC算法实现重建的目标欧式距离矩阵, 其相对误差量级均在10-4, 已然是非常准确的重构结果。如图1、2所示, 选取DM是一组秩为5, 维数不同的欧式距离矩阵。从图1知, 矩阵维数越大, 最小采样率值越小, 即仅需采集非常少量的数据便可准确的重建EDM;图2中随着程序运行时间随着维数的增大而增长。这说明采样FPC算法重建欧式距离矩阵, 在采样率、运行时间及精准度方面均具有优越性, 尤其在重建低秩大型矩阵上采样数目极少且运行速率较快。
3结论
本文利用FPC算法将程序收敛到秩最小化来解决部分欧式距离矩阵重建问题。但是由于待重建的目标矩阵是未知的, 那么要想得到准确的重建效果, 就需要分析观测元素数量, 质点坐标分布结构等问题。大多数EDM在许多方面均有广泛应用, 并且矩阵元间有很大的相关性, 则如何利用这一特性更准确地实现矩阵重建对距离估算问题的研究具有积极意义。
参考文献
[1]A.Y.Alfakih, A.Khandani, H.Wolkowcz.Solving Euclidean distance matrix completion problemsviase midefiniteprogram ming[J], Computational Optimization and Applications, 1999, 12 (1) :13-30.
[2]Ivan Dokmanic, Reza Parhizkar, Juri Ranieri, Martin Vetterli.Euclidean distance matrices essential theory, algorithms and applications[J], Magazine, ISSN:1053-5888.
预处理矩阵 篇7
一、对称矩阵
定义:设A= (aij) n为n阶方阵, 如果满足TA=A, 即aij=aji (i, j=1, 2, ⋅⋅⋅, n) , 那么称A为对称矩阵。
由于对称矩阵形式的特殊性, 使其具有一般矩阵没有的性质, 下面列举出对称矩阵一系列的性质, 并运用对称矩阵的定义和转置运算的性质对每个性质进行了证明。
性质1:A为n阶对称矩阵, 则Am (m为正整数) 也是对称矩阵。
证明:因为A为n阶对称矩阵, 所以TA=A。则 (Am) T= (AT) m=Am, 所以由定义可知Am (m为正整数) 也是对称矩阵。
性质2:A为n阶对称矩阵, 则A+AT也是对称矩阵。
证明:因为 (A+AT) T=AT+ (AT) T=A+AT, 所以A+AT也是对称矩阵。
性质3:A为n阶对称矩阵且A可逆, 则A1-也是对称矩阵。
证明:因为 (A-1) T= (AT) -1=A-1, 所以A-1也是对称矩阵。
性质4:A为m×n阶的矩阵, 则AAT为m阶对称阵, AT A为n阶对称阵。
证明:显然AAT为m阶矩阵, AT A为n阶矩阵, 又由于 (AAT) T= (AT) TAT=AAT, (AT A) T=AT (AT) T=AT A, 所以AAT为m阶对称阵, TA A为n阶对称阵。
性质5:A, B都为n阶对称矩阵, 则A+B也是对称矩阵。
证明:因为 (A+B) T=AT+BT=A+B, 所以A+B也是对称矩阵。
性质6:A, B都为n阶对称矩阵, 则AB也是对称矩阵的充分必要条件是AB=BA。
证明:必要性:设AB为对称矩阵, 则 (AB) T=AB, 而 (AB) T=BT AT=BA, 所以AB=BA。
充分性:设AB=BA, 则 (AB) T= (BA) T=AT BT=AB, 所以AB为对称矩阵。
二、反对称矩阵
定义:设A= (aij) n为n阶方阵, 如果满足TA=-A, 即aij=-aji (i, j=1, 2, ⋅⋅⋅, n) , 那么称A为反对称矩阵。
由于反对称矩阵形式的特殊性, 使其具有了与对称矩阵不同的一些性质。
性质7:设A为n阶反对称矩阵, 则A的主对角线上的元素都为0。
证明:因为A为n阶反对称矩阵, 所以A的主对角线上的元素有aii=-aii (i=1, 2, ⋅⋅⋅, n) , 所以aii=0 (i=1, 2, ⋅⋅⋅, n) 。
性质8:设A为n阶反对称矩阵, n为奇数, 则A的行列式值为0。
证明:因为aij=-aji (i, j=1, 2, ⋅⋅⋅, n) , 所以将A的每一行提出一个公因子-1, 由于n为奇数, 则:。而根据行列式的性质有
性质9:设A为n阶对称矩阵, B为n阶反对称矩阵, 则 (1) AB-BA为对称矩阵。 (2) AB+BA为反对称矩阵。
证明: (1) 因为 (AB-BA) T= (A B) T- (B A) T=BT AT-AT BT=AB-BA, 所以AB-BA为对称矩阵。
(2) 同 (1) , 因为 (AB+BA) T= (AB) T+ (BA) T=BTAT+ATBT=- (AB+BA) , 所以AB+BA为反对称矩阵。
性质10:任一n阶方阵都可以表示为一个对称矩阵和一个反对称矩阵之和。
证明:假设n阶方阵A=B+C, 其中B为对称矩阵, C为反对称矩阵, 则AT= (B+C) T=BT+CT=B-C。由
而, 则B为对称矩阵, C为反对称矩阵, 且A=B+C。
性质11:设A为n阶反对称矩阵, B为n阶对称矩阵, 则AB为反对称矩阵的充分必要条件为AB=BA。
证明:必要性:设AB为反对称矩阵, 则 (AB) T=-AB, 而 (AB) T=BT AT=-BA, 所以AB=BA。
充分性:设AB=BA, 则 (AB) T= (BA) T=ATBT=-AB, 所以AB为反对称矩阵。
三、结束语
对称矩阵与反对称矩阵在高等代数和线形代数中的性质还有很多, 比如对称矩阵的特征值均为实数, 对应不同特征值得的特征向量必正交等等, 由于篇幅所限, 本文只介绍一些基本的性质, 方便读者参考。
参考文献
[1]同济大学应用数学系:《线性代数》.高等教育出版社, 2004
[2]肖马成、周概容:《线性代数、概率论与数理统计证明题500例解析》.高等教育出版社, 2008