這個(gè)程序?qū)崿F(xiàn)了將Matlab當(dāng)中的矩陣輸入到AutoCAD當(dāng)中并制成表格(完整版)實(shí)用資料

這個(gè)程序?qū)崿F(xiàn)了將Matlab當(dāng)中的矩陣輸入到AutoCAD當(dāng)中并制成表格（完整版）實(shí)用資料(可以直接使用，可編輯完整版實(shí)用資料，歡迎下載）

這個(gè)程序?qū)崿F(xiàn)了將Matlab當(dāng)中的矩陣輸入到AutoCAD當(dāng)中并制成表格。

軟件說(shuō)明這個(gè)程序?qū)崿F(xiàn)了將Matlab當(dāng)中的矩陣輸入到AutoCAD當(dāng)中并制成表格（完整版）實(shí)用資料(可以直接使用，可編輯完整版實(shí)用資料，歡迎下載）

這個(gè)程序是在AutoCAD當(dāng)中畫(huà)表格的Matlab程序，在AutoCAD2002和AutoCAD2005當(dāng)中測(cè)試過(guò)。表格是由單行文本和直線組成，文本和表格線分別在不同的層上，方便格式的設(shè)置。AutoCAD2005當(dāng)中已經(jīng)包含了表格對(duì)象，這里沒(méi)有使用這個(gè)對(duì)象。整個(gè)程序包括四個(gè)文件：

acadtalbe.m

畫(huà)表格的Matlab程序的主程序

acadstart.m

啟動(dòng)AutoCAD的Matlab程序

acadquit.m

退出AutoCAD

ACADTable.dvb

在AutoCAD當(dāng)中畫(huà)圖的vba程序使用說(shuō)明

主程序acadtable的調(diào)用方式：

acadtable(m,propname,propvalue)

m:

輸入矩陣，表格的內(nèi)容

propname:

輸入?yún)?shù)名稱

propvalue:

輸入?yún)?shù)的值

propname可以是：

GridHeight

單元格高度,默認(rèn)為20

GridWidth

單元格寬度,默認(rèn)為60

CharHeight

字符高度,默認(rèn)為10

Prec

保留小數(shù)的位數(shù),默認(rèn)為3

isNumber

是否添加行列編號(hào)

調(diào)用示例：

m=rand(10,10);

acadtable(m,'gridheight',10,'gridwidth',30,'charheight',6,'isnumber',true,'prec',3)

水箱水位模糊控制系統(tǒng)設(shè)計(jì)一．在MATLAB命令窗口中輸入sltank，便可打開(kāi)如圖所示的模型窗口。圖1sltank仿真圖（1）打開(kāi)MATLAB，輸入指令fuzzy，打開(kāi)模糊邏輯工具箱的圖形用戶界面窗口，新建一個(gè)Mamdani模糊推理系統(tǒng)。（2）增加一個(gè)輸入變量，將輸入變量命名為水位誤差、誤差變化，將輸出變量命名為閥門開(kāi)關(guān)速度。這樣就建立了一個(gè)兩輸入單輸出的模糊推理系統(tǒng)，保存為shuiwei1。圖2增加一個(gè)輸入變量（3）設(shè)計(jì)模糊化模塊；設(shè)水位誤差level的論域?yàn)閇2.953.05],誤差變化率rate的論域?yàn)閇-0.20.2]；兩個(gè)輸入量的模糊集為level設(shè)為為7個(gè)，rate設(shè)為5個(gè)：其中水位誤差level定為NB、NM、NS、ZE、PS、PM、PB；參數(shù)分別為[0.012.95]、[0.012.97],[0.012.99]、[0.013]、[0.013.01]、[0.013.03]、[0.013.05]，隸屬度均為高斯函數(shù)；圖3輸入量level的參數(shù)設(shè)定誤差變化率rate分別為負(fù)大，負(fù)小，不變，正小，正大。參數(shù)分別為，[0.03-0.2]、[0.03-0.1]、[0.030]、[0.030.1]、[0.03-0.2]，隸屬度函數(shù)均為高斯函數(shù)。圖4誤差變化率rate的參數(shù)設(shè)定閥門的開(kāi)關(guān)速度設(shè)為七個(gè)等級(jí)：快關(guān)，中關(guān)，慢關(guān)，不動(dòng)，慢開(kāi)，中開(kāi)，快開(kāi),其論域?yàn)閇2.953.05]。參數(shù)分別為；[2.942.952.96]、[2.9652.972.975]、[2.992.992.995]、[2.99933.001]、[3.0053.013.015]、[3.023.033.035]、[3.043.053.06]，隸屬函數(shù)為三角形函數(shù)。圖5輸出量valve的參數(shù)設(shè)定（4）設(shè)計(jì)模糊規(guī)則打開(kāi)RuelEditor窗口，通過(guò)選擇添加模糊規(guī)則；1）If(levelisNBand(rateis負(fù)大then(valveis快關(guān)(12）If(levelisNBand(rateis負(fù)小then(valveis快關(guān)(13）If(levelisNBand(rateis不變then(valveis快關(guān)(14）If(levelisNBand(rateis正小then(valveis中關(guān)(15）If(levelisNBand(rateis正大then(valveis不動(dòng)(16）If(levelisNMand(rateis負(fù)大then(valveis快關(guān)(17）If(levelisNMand(rateis負(fù)小then(valveis快關(guān)(18）If(levelisNMand(rateis不變then(valveis快關(guān)(19）If(levelisNMand(rateis正小then(valveis中關(guān)(110）If(levelisNMand(rateis正大then(valveis不動(dòng)(111）If(levelisNSand(rateis負(fù)大then(valveis中關(guān)(112）If(levelisNSand(rateis負(fù)小then(valveis中關(guān)(113）If(levelisNSand(rateis不變then(valveis中關(guān)(114）If(levelisNSand(rateis正小then(valveis不動(dòng)(115）If(levelisNSand(rateis正大then(valveis慢開(kāi)(116）If(levelisZEand(rateis負(fù)大then(valveis中關(guān)(117）If(levelisZEand(rateis負(fù)小then(valveis慢關(guān)(118）If(levelisZEand(rateis不變then(valveis不動(dòng)(119）If(levelisZEand(rateis正小then(valveis慢開(kāi)(120）If(levelisZEand(rateis正大then(valveis中開(kāi)(121）If(levelisPSand(rateis負(fù)大then(valveis慢關(guān)(122）If(levelisPSand(rateis負(fù)小then(valveis不動(dòng)(123）If(levelisPSand(rateis不變then(valveis中開(kāi)(124）If(levelisPSand(rateis正小then(valveis中開(kāi)(125）If(levelisPSand(rateis正大then(valveis中開(kāi)(126）If(levelisPMand(rateis負(fù)大then(valveis不動(dòng)(127）If(levelisPMand(rateis負(fù)小then(valveis中開(kāi)(128）If(levelisPMand(rateis不變then(valveis快開(kāi)(129）If(levelisPMand(rateis正小then(valveis快開(kāi)(130）If(levelisPMand(rateis正大then(valveis快開(kāi)(131）If(levelisPBand(rateis負(fù)大then(valveis不動(dòng)(132）If(levelisPBand(rateis負(fù)小then(valveis中開(kāi)(133）If(levelisPBand(rateis不變then(valveis快開(kāi)(134）If(levelisPBand(rateis正小then(valveis快開(kāi)(135）If(levelisPBand(rateis正大then(valveis快開(kāi)(1這35條模糊控制規(guī)則的權(quán)重都為1.圖6模糊控制規(guī)則的設(shè)定(5利用編輯器的File/SavetoWorkspace,將當(dāng)前的模糊推理系統(tǒng)，以shuiwei1保存到工作空間中。(6在如圖1所示的Simulink仿真系統(tǒng)中，打開(kāi)FuzzyLogicController模糊邏輯控制器模塊對(duì)話框，在其FISFileorStructure參數(shù)對(duì)話框中輸入：shuiwei1。(7在如圖1所示的Simulink系統(tǒng)中，打開(kāi)仿真參數(shù)設(shè)置窗口，正確設(shè)置仿真參數(shù)后，啟動(dòng)仿真便可看到水位變化曲線。圖7水位變化曲線通過(guò)曲面觀察器也可以清晰的看見(jiàn)水箱液位模糊推理的輸入輸出關(guān)系。圖8SurfaceViewer總結(jié)：隨著科學(xué)技術(shù)的發(fā)展，智能控制技術(shù)必會(huì)日趨完善，并且能夠在多領(lǐng)域應(yīng)用。#ifndef

_MATRIX_H

#define

_MATRIX_H

class

Matrix

{

private:

int

row;

//

矩陣的行數(shù)

int

col;

//

矩陣的列數(shù)

int

n;

//

矩陣元素個(gè)數(shù)

double*

mtx;

//

動(dòng)態(tài)分配用來(lái)存放數(shù)組的空間

public:

Matrix(int

row=1,

int

col=1);

//

帶默認(rèn)參數(shù)的構(gòu)造函數(shù)

Matrix(int

row,

int

col,

double

mtx[]);

//

用數(shù)組創(chuàng)建一個(gè)矩陣

Matrix(const

Matrix

&obj);

//

copy構(gòu)造函數(shù)

~Matrix()

{

delete[]

this->mtx;

}

void

print()const;

//

格式化輸出矩陣

int

getRow()const

{

return

this->row;

}//

訪問(wèn)矩陣行數(shù)

int

getCol()const

{

return

this->col;

}//

訪問(wèn)矩陣列數(shù)

int

getN()const

{

return

this->n;

}//

訪問(wèn)矩陣元素個(gè)數(shù)

double*

getMtx()const

{

return

this->mtx;

}//

獲取該矩陣的數(shù)組

//

用下標(biāo)訪問(wèn)矩陣元素

double

get(const

int

i,

const

int

j)const;

//

用下標(biāo)修改矩陣元素值

void

set(const

int

i,

const

int

j,

const

double

e);

//

重載了一些常用操作符，包括

+，-，x，=，負(fù)號(hào)，正號(hào)，

//

A

=

B

Matrix

&operator=

(const

Matrix

&obj);

//

+A

Matrix

operator+

()const

{

return

*this;

}

//

-A

Matrix

operator-

()const;

//

A

+

B

friend

Matrix

operator+

(const

Matrix

&A,

const

Matrix

&B);

//

A

-

B

friend

Matrix

operator-

(const

Matrix

&A,

const

Matrix

&B);

//

A

*

B

兩矩陣相乘

friend

Matrix

operator*

(const

Matrix

&A,

const

Matrix

&B);

//

a

*

B

實(shí)數(shù)與矩陣相乘

friend

Matrix

operator*

(const

double

&a,

const

Matrix

&B);

//

A

的轉(zhuǎn)置

friend

Matrix

trv(const

Matrix

&A);

//

A

的行列式值，采用列主元消去法

//

求行列式須將矩陣化為三角陣,此處為了防止修改原矩陣，采用傳值調(diào)用

friend

double

det(Matrix

A);

//

A

的逆矩陣，采用高斯-若當(dāng)列主元消去法

friend

Matrix

inv(Matrix

A);

};

#endif

實(shí)現(xiàn)

#include<iostream.h>

#include<math.h>

#include<stdlib.h>

#include<iomanip.h>

#include"matrix.h"

//

帶默認(rèn)參數(shù)值的構(gòu)造函數(shù)

//

構(gòu)造一個(gè)row行col列的零矩陣

Matrix::Matrix(int

row,

int

col){

this->row

=

row;

this->col

=

col;

this->n

=

row

*

col;

this->mtx

=

new

double[n];

for(int

i=0;

i<n;

i++)

this->mtx[i]

=

0.0;

}

//

用一個(gè)數(shù)組初始化矩陣

Matrix::Matrix(int

row,

int

col,

double

mtx[])

{

this->row

=

row;

this->col

=

col;

this->n

=

row

*

col;

this->mtx

=

new

double[n];

for(int

i=0;

i<n;

i++)

this->mtx[i]

=

mtx[i];

}

//

拷貝構(gòu)造函數(shù)，因?yàn)槌蓡T變量含有動(dòng)態(tài)空間，防止傳遞參數(shù)

//

等操作發(fā)生錯(cuò)誤

Matrix::Matrix(const

Matrix

&obj){

this->row

=

obj.getRow();

this->col

=

obj.getCol();

this->n

=

obj.getN();

this->mtx

=

new

double[n];

for(int

i=0;

i<n;

i++)

this->mtx[i]

=

obj.getMtx()[i];

}

//

格式化輸出矩陣所有元素

void

Matrix::print()const{

for(int

i=0;

i<this->row;

i++){

for(int

j=0;

j<this->col;

j++)

if(fabs(this->get(i,j))

<=

1.0e-10)

cout

<<

setiosflags(ios::left)

<<

setw(12)

<<

0.0

<<

'

';

else

cout

<<

setiosflags(ios::left)

<<

setw(12)

<<

this->get(i,j)

<<

'

';

cout

<<endl;

}

}

//

獲取矩陣元素

//

注意這里矩陣下標(biāo)從(0,0)開(kāi)始

double

Matrix::get(const

int

i,

const

int

j)const{

return

this->mtx[i*this->col

+

j];

}

//

修改矩陣元素

void

Matrix::set(const

int

i,

const

int

j,

const

double

e){

this->mtx[i*this->col

+

j]

=

e;

}

//

重載賦值操作符，由于成員變量中含有動(dòng)態(tài)分配

Matrix

&Matrix::operator=

(const

Matrix

&obj){

if(this

==

&obj)

//

將一個(gè)矩陣賦給它自己時(shí)簡(jiǎn)單做返回即可

return

*this;

delete[]

this->mtx;

//

首先刪除目的操作數(shù)的動(dòng)態(tài)空間

this->row

=

obj.getRow();

this->col

=

obj.getCol();

this->n

=

obj.getN();

this->mtx

=

new

double[n];

//

重新分配空間保存obj的數(shù)組

for(int

i=0;

i<n;

i++)

this->mtx[i]

=

obj.getMtx()[i];

return

*this;

}

//

負(fù)號(hào)操作符，返回值為該矩陣的負(fù)矩陣，原矩陣不變

Matrix

Matrix::operator-

()const{

//

為了不改變?cè)瓉?lái)的矩陣，此處從新構(gòu)造一個(gè)矩陣

Matrix

_A(this->row,

this->col);

for(int

i=0;

i<_A.n;

i++)

_A.mtx[i]

=

-(this->mtx[i]);

return

_A;

}

//

矩陣求和，對(duì)應(yīng)位置元素相加

Matrix

operator+

(const

Matrix

&A,

const

Matrix

&B){

Matrix

AB(A.row,

A.col);

if(A.row!=B.row

||

A.col!=B.col){

cout

<<

"Can't

do

A+B\n";

//

如果矩陣A和B行列數(shù)不一致則不可相加

exit(0);

}

for(int

i=0;

i<AB.n;

i++)

AB.mtx[i]

=

A.mtx[i]

+

B.mtx[i];

return

AB;

}

//

矩陣減法，用加上一個(gè)負(fù)矩陣來(lái)實(shí)現(xiàn)

Matrix

operator-

(const

Matrix

&A,

const

Matrix

&B){

return

(A

+

(-B));

}

//

矩陣乘法

Matrix

operator*

(const

Matrix

&A,

const

Matrix

&B){

if(A.col

!=

B.row){

//

A的列數(shù)必須和B的行數(shù)一致

cout

<<

"Can't

multiply\n";

exit(0);

}

Matrix

AB(A.row,

B.col);

//

AB用于保存乘積

for(int

i=0;

i<AB.row;

i++)

for(int

j=0;

j<AB.col;

j++)

for(int

k=0;

k<A.col;

k++)

AB.set(i,

j,

AB.get(i,j)

+

A.get(i,k)*B.get(k,j));

return

AB;

}

//

矩陣與實(shí)數(shù)相乘

Matrix

operator*

(const

double

&a,

const

Matrix

&B){

Matrix

aB(B);

for(int

i=0;

i<aB.row;

i++)

for(int

j=0;

j<aB.col;

j++)

aB.set(i,j,

a*B.get(i,j));

return

aB;

}

//

矩陣的轉(zhuǎn)置

將(i,j)與(j,i)互換

//

此函數(shù)返回一個(gè)矩陣的轉(zhuǎn)置矩陣，并不改變?cè)瓉?lái)的矩陣

Matrix

trv(const

Matrix

&A){

Matrix

AT(A.col,

A.row);

for(int

i=0;

i<AT.row;

i++)

for(int

j=0;

j<AT.col;

j++)

AT.set(i,

j,

A.get(j,i));

return

AT;

}

//

矩陣行列式值，采用列主元消去法

double

det(Matrix

A){

if(A.row

!=

A.col)

{

//

矩陣必須為n*n的才可進(jìn)行行列式求值

cout

<<

"error"

<<

endl;

return

0.0;

//如果不滿足行列數(shù)相等返回0.0

}

double

detValue

=

1.0;

//

用于保存行列式值

for(int

i=0;

i<A.getRow()-1;

i++){

//

需要n-1步列化零操作

//------------------

選主元

---------------------------------

double

max

=

fabs(A.get(i,i));

//

主元初始默認(rèn)為右下方矩陣首個(gè)元素

int

ind

=

i;

//

主元行號(hào)默認(rèn)為右下方矩陣首行

for(int

j=i+1;

j<A.getRow();

j++){

//

選擇列主元

if(fabs(A.get(j,i))

>

max){

//

遇到絕對(duì)值更大的元素

max

=

fabs(A.get(j,i));

//

更新主元值

ind

=

j;

//

更新主元行號(hào)

}

}//loop

j

//-------------------

移動(dòng)主元行

-----------------------------

if(max

<=

1.0e-10)

return

0.0;

//

右下方矩陣首行為零，顯然行列式值為零

if(ind

!=

i){//

主元行非右下方矩陣首行

for(int

k=i;

k<A.getRow();

k++){

//

將主元行與右下方矩陣首行互換

double

temp

=

A.get(i,k);

A.set(i,k,A.get(ind,k));

A.set(ind,k,temp);

}

detValue

=

-detValue;

//

互換行列式兩行，行列式值反號(hào)

}

//-------------------

消元

----------------------------------

for(j=i+1;

j<A.getRow();

j++){

//

遍歷行

double

temp

=

A.get(j,i)/A.get(i,i);

for(int

k=i;

k<A.getRow();

k++)

//

遍歷行中每個(gè)元素，行首置0

A.set(j,

k,

A.get(j,k)-A.get(i,k)*temp);

}

detValue

*=

A.get(i,i);

//

每步消元都會(huì)產(chǎn)生一個(gè)對(duì)角線上元素，將其累乘

}//

loop

i

//

注意矩陣最后一個(gè)元素在消元的過(guò)程中沒(méi)有被累乘到

return

detValue

*

A.get(A.getRow()-1,

A.getRow()-1);}//det()

//

A的逆矩陣

高斯-若當(dāng)消去法，按列選主元

Matrix

inv(Matrix

A){

if(A.row

!=

A.col){

//

只可求狹義逆矩陣，即行列數(shù)相同

cout

<<

"Matrix

should

be

N

x

N\n";

exit(0);

}

//

構(gòu)造一個(gè)與A行列相同的單位陣B

Matrix

B(A.row,A.col);

for(int

r=0;

r<A.row;

r++)

for(int

c=0;

c<A.col;

c++)

if(r

==

c)

B.set(r,c,1.0);

//

對(duì)矩陣A進(jìn)行A.row次消元運(yùn)算，每次保證第K列只有對(duì)角線上非零

//

同時(shí)以同樣的操作施與矩陣B，結(jié)果A變?yōu)閱挝魂嘊為所求逆陣

for(int

k=0;

k<A.row;

k++){

//------------------

選主元

--------------------------------------

double

max

=

fabs(A.get(k,k));

//

主元初始默認(rèn)為右下方矩陣首個(gè)元素

int

ind

=

k;

//

主元行號(hào)默認(rèn)為右下方矩陣首行

//

結(jié)果第ind行為列主元行

for(int

n=k+1;

n<A.getRow();

n++){

if(fabs(A.get(n,k))

>

max){//

遇到絕對(duì)值更大的元素

max

=

fabs(A.get(n,k));

//

更新主元值

ind

=

n;//

更新主元行號(hào)

}

}

//-------------------

移動(dòng)主元行

--------------------------------

if(ind

!=

k){//

主元行不是右下方矩陣首行

for(int

m=k;

m<A.row;

m++){//

將主元行與右下方矩陣首行互換

double

tempa

=

A.get(k,m);

A.set(k,

m,

A.get(ind,m));

A.set(ind,

m,

tempa);

}

for(m=0;

m<B.row;

m++){

double

tempb

=

B.get(k,m);

//

對(duì)矩陣B施以相同操作

B.set(k,

m,

B.get(ind,m));

//

B與A階數(shù)相同，可在一個(gè)循環(huán)中

B.set(ind,

m,

tempb);

}

}

//---------------------

消元

-----------------------------------

//

第k次消元操作，以第k行作為主元行，將其上下各行的第k列元素化為零

//

同時(shí)以同樣的參數(shù)對(duì)B施以同樣的操作，此時(shí)可以將B看作A矩陣的一部分

for(int

i=0;

i<A.col;

i++){

if(i

!=

k){

double

Mik

=

-A.get(i,k)/A.get(k,k);

for(int

j=k+1;

j<A.row;

j++)

A.set(i,

j,

A.get(i,j)

+

Mik*A.get(k,j));

for(j=0;

j<B.row;

j++)

B.set(i,

j,

B.get(i,j)

+

Mik*B.get(k,j));

}//end

if

}//loop

i

double

Mkk

=

1.0/A.get(k,k);

for(int

j=0;

j<A.row;

j++)

A.set(k,

j,

A.get(k,j)

*

Mkk);

for(j=0;

j<B.row;

j++)

B.set(k,

j,

B.get(k,j)

*

Mkk);

}//loop

k

return

B;

}//inv()9矩陣特征值計(jì)算在實(shí)際的工程計(jì)算中，經(jīng)常會(huì)遇到求n階方陣A的特征值(Eigenvalue)與特征向量(Eigenvector)的問(wèn)題。對(duì)于一個(gè)方陣A，如果數(shù)值λ使方程組Ax=λx即(A-λIn)x=0有非零解向量(SolutionVector)x，則稱λ為方陣A的特征值，而非零向量x為特征值λ所對(duì)應(yīng)的特征向量，其中In為n階單位矩陣。由于根據(jù)定義直接求矩陣特征值的過(guò)程比較復(fù)雜，因此在實(shí)際計(jì)算中，往往采取一些數(shù)值方法。本章主要介紹求一般方陣絕對(duì)值最大的特征值的乘冪(Power)法、求對(duì)稱方陣特征值的雅可比法和單側(cè)旋轉(zhuǎn)(One-sideRotation)法以及求一般矩陣全部特征值的QR方法及一些相關(guān)的并行算法。求解矩陣最大特征值的乘冪法乘冪法及其串行算法在許多實(shí)際問(wèn)題中，只需要計(jì)算絕對(duì)值最大的特征值，而并不需要求矩陣的全部特征值。乘冪法是一種求矩陣絕對(duì)值最大的特征值的方法。記實(shí)方陣A的n個(gè)特征值為λii=(1,2,…,n)，且滿足：│λ1│≥│λ2│≥│λ3│≥…≥│λn│特征值λi對(duì)應(yīng)的特征向量為xi。乘冪法的做法是：①取n維非零向量v0作為初始向量；②對(duì)于k=1,2,…,做如下迭代：uk=Avk-1vk=uk/║uk║∞直至ε為止，這時(shí)vk+1就是A的絕對(duì)值最大的特征值λ1所對(duì)應(yīng)的特征向量x1。若vk-1與vk的各個(gè)分量同號(hào)且成比例，則λ1=║uk║∞；若vk-1與vk的各個(gè)分量異號(hào)且成比例，則λ1=-║uk║∞。若各取一次乘法和加法運(yùn)算時(shí)間、一次除法運(yùn)算時(shí)間、一次比較運(yùn)算時(shí)間為一個(gè)單位時(shí)間，則因?yàn)橐惠営?jì)算要做一次矩陣向量相乘、一次求最大元操作和一次規(guī)格化操作，所以下述乘冪法串行算法21.1的一輪計(jì)算時(shí)間為n2+2n=O(n2)。算法21.1單處理器上乘冪法求解矩陣最大特征值的算法輸入：系數(shù)矩陣An×n，初始向量vn×1，ε輸出：最大的特征值maxBegin while(│diff│>ε)do(1)fori=1tondo(1.1)sum=0(1.2)forj=1tondosum=sum+a[i,j]*x[j]endfor(1.3)b[i]=sumendfor(2)max=│b[1]│(3)fori=2tondoif(│b[i]│>max)thenmax=│b[i]│endifendfor(4)fori=1tondox[i]=b[i]/maxendfor(5)diff=max-oldmax,oldmax=maxendwhileEnd乘冪法并行算法乘冪法求矩陣特征值由反復(fù)進(jìn)行矩陣向量相乘來(lái)實(shí)現(xiàn)，因而可以采用矩陣向量相乘的數(shù)據(jù)劃分方法。設(shè)處理器個(gè)數(shù)為p，對(duì)矩陣A按行劃分為p塊，每塊含有連續(xù)的m行向量，這里，編號(hào)為i的處理器含有A的第im至第(i+1)m-1行數(shù)據(jù)，(i=0,1,…,p-1)，初始向量v被廣播給所有處理器。各處理器并行地對(duì)存于局部存儲(chǔ)器中a的行塊和向量v做乘積操作，并求其m維結(jié)果向量中的最大值localmax，然后通過(guò)歸約操作的求最大值運(yùn)算得到整個(gè)n維結(jié)果向量中的最大值max并廣播給所有處理器，各處理器利用max進(jìn)行規(guī)格化操作。最后通過(guò)擴(kuò)展收集操作將各處理器中的維結(jié)果向量按處理器編號(hào)連接起來(lái)并廣播給所有處理器，以進(jìn)行下一次迭代計(jì)算，直至收斂。具體算法框架描述如下：算法21.2乘冪法求解矩陣最大特征值的并行算法輸入：系數(shù)矩陣An×n，初始向量vn×1，ε輸出：最大的特征值maxBegin 對(duì)所有處理器my_rank(my_rank=0,…,p-1)同時(shí)執(zhí)行如下的算法: while(│diff│>ε)do/*diff為特征向量的各個(gè)分量新舊值之差的最大值*/(1)fori=0tom-1do/*對(duì)所存的行計(jì)算特征向量的相應(yīng)分量*/(1.1)sum=0(1.2)forj=0ton-1dosum=sum+a[i,j]*x[j]endfor(1.3)b[i]=sumendfor(2)localmax=│b[0]│/*對(duì)所計(jì)算的特征向量的相應(yīng)分量求新舊值之差的最大值localmax*/(3)fori=1tom-1doif(│b[i]│>localmax)thenlocalmax=│b[i]│endifendfor(4)用Allreduce操作求出所有處理器中l(wèi)ocalmax值的最大值max并廣播到所有處理器中 (5)fori=0tom-1do/*對(duì)所計(jì)算的特征向量歸一化*/b[i]=b[i]/maxendfor(6)用Allgather操作將各個(gè)處理器中計(jì)算出的特征向量的分量的新值組合并廣播到所有處理器中(7)diff=max-oldmax,oldmax=maxendwhileEnd若各取一次乘法和加法運(yùn)算時(shí)間、一次除法運(yùn)算時(shí)間、一次比較運(yùn)算時(shí)間為一個(gè)單位時(shí)間，一輪迭代的計(jì)算時(shí)間為mn+2m；一輪迭代中，各處理器做一次歸約操作，通信量為1，一次擴(kuò)展收集操作，通信量為m，則通信時(shí)間為。因此乘冪法的一輪并行計(jì)算時(shí)間為。MPI源程序請(qǐng)參見(jiàn)所附光盤。求對(duì)稱矩陣特征值的雅可比法雅可比法求對(duì)稱矩陣特征值的串行算法若矩陣A=[aij]是n階實(shí)對(duì)稱矩陣，則存在一個(gè)正交矩陣U，使得UTAU=D其中D是一個(gè)對(duì)角矩陣，即這里λi(i=1,2,…,n)是A的特征值，U的第i列是與λi對(duì)應(yīng)的特征向量。雅可比算法主要是通過(guò)正交相似變換將一個(gè)實(shí)對(duì)稱矩陣對(duì)角化，從而求出該矩陣的全部特征值和對(duì)應(yīng)的特征向量。因此可以用一系列的初等正交變換逐步消去A的非對(duì)角線元素，從而使矩陣A對(duì)角化。設(shè)初等正交矩陣為R(p,q,θ)，其(p,p)(q,q)位置的數(shù)據(jù)是cosθ，(p,q)(q,p)位置的數(shù)據(jù)分別是-sinθ和sinθ(p≠q)，其它位置的數(shù)據(jù)和一個(gè)同階數(shù)的單位矩陣相同。顯然可以得到：R(p,q,θ)TR(p,q,θ)=In不妨記B=R(p,q,θ)TAR(p,q,θ)，如果將右端展開(kāi)，則可知矩陣B的元素與矩陣A的元素之間有如下關(guān)系：bpp=appcos2θ+aqqsin2θ+apqsin2θbqq=appsin2θ+aqqcos2θ-apqsin2θbpq=((aqq-app)sin2θ)/2+apqcos2θbqp=bpqbpj=apjcosθ+aqjsinθbqj=-apjsinθ+aqjcosθbip=aipcosθ+aiqsinθbiq=-aipsinθ+aiqcosθbij=aij其中1≤i,j≤n且i,j≠p,q。因?yàn)锳為對(duì)稱矩陣，R為正交矩陣，所以B亦為對(duì)稱矩陣。若要求矩陣B的元素bpq=0，則只需令((aqq-app)sin2θ)/2+apqcos2θ=0，即：在實(shí)際應(yīng)用時(shí)，考慮到并不需要解出θ，而只需要求出sin2θ，sinθ和cosθ就可以了。如果限制θ值在-π/2<2θ≤π/2，則可令容易推出：,,利用sin2θ，sinθ和cosθ的值，即得矩陣B的各元素。矩陣A經(jīng)過(guò)旋轉(zhuǎn)變換，選定的非主對(duì)角元素apq及aqp（一般是絕對(duì)值最大的）就被消去，且其主對(duì)角元素的平方和增加了，而非主對(duì)角元素的平方和減少了2a2pq，矩陣元素總的平方和不變。通過(guò)反復(fù)選取主元素apq，并作旋轉(zhuǎn)變換，就逐步將矩陣A變?yōu)閷?duì)角矩陣。在對(duì)稱矩陣中共有(n2-n)/2個(gè)非主對(duì)角元素要被消去,而每消去一個(gè)非主對(duì)角元素需要對(duì)2n個(gè)元素進(jìn)行旋轉(zhuǎn)變換,對(duì)一個(gè)元素進(jìn)行旋轉(zhuǎn)變換需要2次乘法和1次加法，若各取一次乘法運(yùn)算時(shí)間或一次加法運(yùn)算時(shí)間為一個(gè)單位時(shí)間，則消去一個(gè)非主對(duì)角元素需要3個(gè)單位運(yùn)算時(shí)間，所以下述算法21.3的一輪計(jì)算時(shí)間復(fù)雜度為(n2-n)/2*2n*3=3n2(n-1)=O(n3)。算法21.3單處理器上雅可比法求對(duì)稱矩陣特征值的算法輸入：矩陣An×n，ε輸出：矩陣特征值EigenvalueBegin (1)while(max>ε)do (1.1)max=a[1,2] (1.2)fori=1tondo forj=i+1tondo if(│a[i,j])│>max)thenmax=│a[i,j]│，p=i，q=jendif endfor endfor(1.3)Compute:m=-a[p,q]，n=(a[q,q]-a[p,p])/2，w=sgn(n)*m/sqrt(m*m+n*n)，si2=w，si1=w/sqrt(2(1+sqrt(1-w*w)))，co1=sqrt(1-si1*si1)，b[p,p]=a[p,p]*co1*co1+a[q,q]*si1*si1+a[p,q]*si2，b[q,q]=a[p,p]*si1*si1+a[q,q]*co1*co1-a[p,q]*si2，b[q,p]=0，b[p,q]=0(1.4)forj=1tondoif((j≠p)and(j≠q))then (i)b[p,j]=a[p,j]*co1+a[q,j]*si1 (ii)b[q,j]=-a[p,j]*si1+a[q,j]*co1endifendfor(1.5)fori=1tondoif((i≠p)and(i≠q))then(i)b[i,p]=a[i,p]*co1+a[i,q]*si1 (ii)b[i,q]=-a[i,p]*si1+a[i,q]*co1endifendfor(1.6)fori=1tondoforj=1tondoa[i,j]=b[i,j]endforendforendwhile(2)fori=1tondo Eigenvalue[i]=a[i,i]endforEnd雅可比法求對(duì)稱矩陣特征值的并行算法串行雅可比算法逐次尋找非主對(duì)角元絕對(duì)值最大的元素的方法并不適合于并行計(jì)算。因此，在并行處理中，我們每次并不尋找絕對(duì)值最大的非主對(duì)角元消去，而是按一定的順序?qū)中的所有上三角元素全部消去一遍，這樣的過(guò)程稱為一輪。由于對(duì)稱性，在一輪中，A的下三角元素也同時(shí)被消去一遍。經(jīng)過(guò)若干輪，可使A的非主對(duì)角元的絕對(duì)值減少，收斂于一個(gè)對(duì)角方陣。具體算法如下：設(shè)處理器個(gè)數(shù)為p，對(duì)矩陣A按行劃分為2p塊，每塊含有連續(xù)的m/2行向量，記，這些行塊依次記為A0,A1,…,A2p-1，并將A2i與A2i+1存放與標(biāo)號(hào)為i的處理器中。每輪計(jì)算開(kāi)始，各處理器首先對(duì)其局部存儲(chǔ)器中所存的兩個(gè)行塊的所有行兩兩配對(duì)進(jìn)行旋轉(zhuǎn)變換，消去相應(yīng)的非對(duì)角線元素。然后按圖21.1所示規(guī)律將數(shù)據(jù)塊在不同處理器之間傳送，以消去其它非主對(duì)角元素。圖STYLEREF1\s21.SEQ圖\*ARABIC\s11p=4時(shí)的雅可比算法求對(duì)稱矩陣特征值的數(shù)據(jù)交換圖這里長(zhǎng)方形框中兩個(gè)方格內(nèi)的整數(shù)被看成是所移動(dòng)的行塊的編號(hào)。在要構(gòu)成新的行塊配對(duì)時(shí)，只要將每個(gè)處理器所對(duì)應(yīng)的行塊按箭頭方向移至相鄰的處理器即可，這樣的傳送可以在行塊之間實(shí)現(xiàn)完全配對(duì)。當(dāng)編號(hào)為i和j的兩個(gè)行塊被送至同一處理器時(shí)，令編號(hào)為i的行塊中的每行元素和編號(hào)為j的行塊中的每行元素配對(duì)，以消去相應(yīng)的非主對(duì)角元素，這樣在所有的行塊都兩兩配對(duì)之后，可以將所有的非主對(duì)角元素消去一遍，從而完成一輪計(jì)算。由REF_Ref37669518\h圖21.1可以看出，在一輪計(jì)算中，處理器之間要2p-1次交換行塊。為保證不同行塊配對(duì)時(shí)可以將原矩陣A的非主對(duì)角元素消去，引入變量b記錄每個(gè)處理器中的行塊數(shù)據(jù)在原矩陣A中的實(shí)際行號(hào)。在行塊交換時(shí)，變量b也跟著相應(yīng)的行塊在各處理器之間傳送。處理器之間的數(shù)據(jù)塊交換存在如下規(guī)律：0號(hào)處理器前一個(gè)行塊(簡(jiǎn)稱前數(shù)據(jù)塊，后一個(gè)行塊簡(jiǎn)稱后數(shù)據(jù)塊)始終保持不變，將后數(shù)據(jù)塊發(fā)送給1號(hào)處理器，作為1號(hào)處理器的前數(shù)據(jù)塊。同時(shí)接收1號(hào)處理器發(fā)送的后數(shù)據(jù)塊作為自己的后數(shù)據(jù)塊。p-1號(hào)處理器首先將后數(shù)據(jù)塊發(fā)送給編號(hào)為p-2的處理器，作為p-2號(hào)處理器的后數(shù)據(jù)塊。然后將前數(shù)據(jù)塊移至后數(shù)據(jù)塊的位置上，最后接收p-2號(hào)處理器發(fā)送的前數(shù)據(jù)塊作為自己的前數(shù)據(jù)塊。編號(hào)為my_rank的其余處理器將前數(shù)據(jù)塊發(fā)送給編號(hào)為my_rank+1的處理器，作為my_rank+1號(hào)處理器的前數(shù)據(jù)塊。將后數(shù)據(jù)塊發(fā)送給編號(hào)為my_rank-1的處理器，作為my_rank-1號(hào)處理器的后數(shù)據(jù)塊。為了避免在通信過(guò)程中發(fā)生死鎖，奇數(shù)號(hào)處理器和偶數(shù)號(hào)處理器的收發(fā)順序被錯(cuò)開(kāi)。使偶數(shù)號(hào)處理器先發(fā)送后接收；而奇數(shù)號(hào)處理器先將數(shù)據(jù)保存在緩沖區(qū)中，然后接收偶數(shù)號(hào)處理器發(fā)送的數(shù)據(jù)，最后再將緩沖區(qū)中的數(shù)據(jù)發(fā)送給偶數(shù)號(hào)處理器。數(shù)據(jù)塊傳送的具體過(guò)程描述如下：算法21.4雅可比法求對(duì)稱矩陣特征值的并行過(guò)程中處理器之間的數(shù)據(jù)塊交換算法輸入：矩陣A的行數(shù)據(jù)塊和向量b的數(shù)據(jù)塊分布于各個(gè)處理器中輸出：在處理器陣列中傳送后的矩陣A的行數(shù)據(jù)塊和向量b的數(shù)據(jù)塊Begin 對(duì)所有處理器my_rank(my_rank=0,…,p-1)同時(shí)執(zhí)行如下的算法: /*矩陣A和向量b為要傳送的數(shù)據(jù)塊*/ (1)if(my-rank=0)then/*0號(hào)處理器*/ (1.1)將后數(shù)據(jù)塊發(fā)送給1號(hào)處理器 (1.2)接收1號(hào)處理器發(fā)送來(lái)的后數(shù)據(jù)塊作為自己的后數(shù)據(jù)塊endif(2)if((my-rank=p-1)and(my-rankmod2≠0))then/*處理器p-1且其為奇數(shù)*/ (2.1)fori=m/2tom-1do/*將后數(shù)據(jù)塊保存在緩沖區(qū)buffer中*/ forj=0ton-1do buffer[i-m/2,j]=a[i,j] endfor endfor (2.2)fori=m/2tom-1do buf[i-m/2]=b[i] endfor (2.3)fori=0tom/2-1do/*將前數(shù)據(jù)塊移至后數(shù)據(jù)塊的位置上*/ forj=0ton-1do a[i+m/2,j]=a[i,j] endfor endfor (2.4)fori=0tom/2-1do b[i+m/2]=b[i]endfor (2.5)接收p-2號(hào)處理器發(fā)送的數(shù)據(jù)塊作為自己的前數(shù)據(jù)塊 (2.6)將buffer中的后數(shù)據(jù)塊發(fā)送給編號(hào)為p-2的處理器endif(3)if((my-rank=p-1)and(my-rankmod2=0))then/*處理器p-1且其為偶數(shù)*/ (3.1)將后數(shù)據(jù)塊發(fā)送給編號(hào)為p-2的處理器 (3.2)fori=0tom/2-1do/*將前數(shù)據(jù)塊移至后數(shù)據(jù)塊的位置上*/ forj=0ton-1do a[i+m/2,j]=a[i,j] endfor endfor (3.3)fori=0tom/2-1do b[i+m/2]=b[i] endfor (3.4)接收p-2號(hào)處理器發(fā)送的數(shù)據(jù)塊作為自己的前數(shù)據(jù)塊endif(4)if((my-rank≠p-1)and(my-rank≠0))then/*其它的處理器*/ (4.1)if(my-rankmod2=0)then/*偶數(shù)號(hào)處理器*/(i)將前數(shù)據(jù)塊發(fā)送給編號(hào)為my_rank+1的處理器(ii)將后數(shù)據(jù)塊發(fā)送給編號(hào)為my_rank-1的處理器(ii)接收編號(hào)為my_rank-1的處理器發(fā)送的數(shù)據(jù)塊作為自己的前數(shù)據(jù)塊(iv)接收編號(hào)為my_rank+1的處理器發(fā)送的數(shù)據(jù)塊作為自己的后數(shù)據(jù)塊else/*奇數(shù)號(hào)處理器*/(v)fori=0tom-1do/*將前后數(shù)據(jù)塊保存在緩沖區(qū)buffer中*/ forj=0ton-1do buffer[i,j]=a[i,j] endforendfor(vi)fori=0tom-1do buf[i]=b[i]endfor(vii)接收編號(hào)為my_rank-1的處理器發(fā)送的數(shù)據(jù)塊作為自己的前數(shù)據(jù)塊(viii)接收編號(hào)為my_rank+1的處理器發(fā)送的數(shù)據(jù)塊作為自己的后數(shù)據(jù)塊(ix)將存于buffer中的前數(shù)據(jù)塊發(fā)送給編號(hào)為my_rank+1的處理器(x)將存于buffer中的后數(shù)據(jù)塊發(fā)送給編號(hào)為my_rank-1的處理器endifendifEnd各處理器并行地對(duì)其局部存儲(chǔ)器中的非主對(duì)角元素aij進(jìn)行消去，首先計(jì)算旋轉(zhuǎn)參數(shù)并對(duì)第i行和第j行兩行元素進(jìn)行旋轉(zhuǎn)行變換。然后通過(guò)擴(kuò)展收集操作將相應(yīng)的旋轉(zhuǎn)參數(shù)及第i列和第j列的列號(hào)按處理器編號(hào)連接起來(lái)并廣播給所有處理器。各處理器在收到這些旋轉(zhuǎn)參數(shù)和列號(hào)之后，按0,1,…,p-1的順序依次讀取旋轉(zhuǎn)參數(shù)及列號(hào)并對(duì)其m行中的第i列和第j列元素進(jìn)行旋轉(zhuǎn)列變換。經(jīng)過(guò)一輪計(jì)算的2p-1次數(shù)據(jù)交換之后，原矩陣A的所有非主對(duì)角元素都被消去一次。此時(shí)，各處理器求其局部存儲(chǔ)器中的非主對(duì)角元素的最大元localmax，然后通過(guò)歸約操作的求最大值運(yùn)算求得將整個(gè)n階矩陣非主對(duì)角元素的最大元max，并廣播給所有處理器以決定是否進(jìn)行下一輪迭代。具體算法框架描述如下：算法21.5雅可比法求對(duì)稱矩陣特征值的并行算法輸入：矩陣An×n，ε輸出：矩陣A的主對(duì)角元素即為A的特征值Begin 對(duì)所有處理器my_rank(my_rank=0,…,p-1)同時(shí)執(zhí)行如下的算法: (a)fori=0tom-1do b[i]=myrank*m+i/*b記錄處理器中的行塊數(shù)據(jù)在原矩陣A中的實(shí)際行號(hào)*/endfor (b)while(│max│>ε)do/*max為A中所有非對(duì)角元最大的絕對(duì)值*/ (1)fori=my_rank*mto(my_rank+1)*m-2do/*對(duì)本處理器內(nèi)部所有行兩兩配對(duì)進(jìn)行旋轉(zhuǎn)變換*/ forj=i+1to(my_rank+1)*m-1do(1.1)r=imodm,t=jmodm/*i,j為進(jìn)行旋轉(zhuǎn)變換行(稱為主行)的實(shí)際行號(hào),r,t為它們?cè)趬K內(nèi)的相對(duì)行號(hào)*/ (1.2)if(a[r,j]≠0)then/*對(duì)四個(gè)主元素的旋轉(zhuǎn)變換*/ (i)Compute:f=-a[r,j],g=(a[t,j]-a[r,i])/2,h=sgn(g)*f/sqrt(f*f+g*g),sin2=h,sin1=h/sqrt(2*(1+sqrt(1-h*h))),cos1=sqrt(1-sin1*sin1)，bpp=a[r,i]*cos1*cos1+a[t,j]*sin1*sin1+a[r,j]*sin2，bqq=a[r,i]*sin1*sin1+a[t,j]*cos1*cos1-a[r,j]*sin2，bpq=0,bqp=0 (ii)forv=0ton-1do/*對(duì)兩個(gè)主行其余元素的旋轉(zhuǎn)變換*/ if((v≠i)and(v≠j))then br[v]=a[r,v]*cos1+a[t,v]*sin1 a[t,v]=-a[r,v]*sin1+a[t,v]*cos1 endifendfor(iii)forv=0ton-1do if((v≠i)and(v≠j))then a[r,v]=br[v] endifendfor (iv)forv=0tom-1do/*對(duì)兩個(gè)主列在本處理器內(nèi)的其余元素的旋轉(zhuǎn)變換*/if((v≠r)and(v≠t))then bi[v]=a[v,i]*cos1+a[v,j]*sin1 a[v,j]=-a[v,i]*sin1+a[v,j]*cos1endifendfor (v)forv=0tom-1doif((v≠r)and(v≠t))then a[v,i]=bi[v]endifendfor (vi)Compute: a[r,i]=bpp,a[t,j]=bqq,a[r,j]=bpq,a[t,i]=bqp，/*用temp1保存本處理器主行的行號(hào)和旋轉(zhuǎn)參數(shù)*/temp1[0]=sin1,temp1[1]=cos1，temp1[2]=(float)i,temp1[3]=(float)jelse(vii)Compute: temp1[0]=0,temp1[1]=0，temp1[2]=0,temp1[3]=0endif (1.3)將所有處理器temp1中的旋轉(zhuǎn)參數(shù)及主行的行號(hào) 按處理器編號(hào)連接起來(lái)并廣播給所有處理器,存于temp2中 (1.4)current=0 (1.5) forv=1topdo/*根據(jù)temp2中的其它處理器的旋轉(zhuǎn)參數(shù)及主行的行號(hào)對(duì)相關(guān)的列在本處理器的部分進(jìn)行旋轉(zhuǎn)變換*/ (i)Compute: s1=temp2[(v-1)*4+0],c1=temp2[(v-1)*4+1]， i1=(int)temp2[(v-1)*4+2],j1=(int)temp2[(v-1)*4+3] (ii)if(s1、c1、i1、j1中有一不為0)then if(my-rank≠current)then forz=0tom-1do zi[z]=a[z,i1]*c1+a[z,j1]*s1 a[z,j1]=-a[z,i1]*s1+a[z,j1]*c1 endfor forz=0tom-1do a[z,i1]=zi[z] endfor endifendif (iii)current=current+1endfor endforendfor (2)forcounter=1to2p-2do/*進(jìn)行2p-2次處理器間的數(shù)據(jù)交換,并對(duì)交換后處理器中所有行兩兩配對(duì)進(jìn)行旋轉(zhuǎn)變換*/ (2.1)Data_exchange()/*處理器間的數(shù)據(jù)交換*/ (2.2)fori=0tom/2-1doforj=m/2tom-1do (i)if(a[i,b[j]]≠0)then/*對(duì)四個(gè)主元素的旋轉(zhuǎn)變換*/ ①Compute:f=-a[i,b[j]]，g=(a[j,b[j]]-a[i,b[i]])/2，h=sgn(g)*f/sqrt(f*f+g*g)，sin2=h,sin1=h/sqrt(2*(1+sqrt(1-h*h)))，cos1=sqrt(1-sin1*sin1)，bpp=a[i,b[i]]*cos1*cos1+a[j,b[j]]*sin1*sin1+a[i,b[j]]*sin2，bqq=a[i,b[i]]*sin1*sin1+a[j,b[j]]*cos1*cos1-a[i,b[j]]*sin2，bpq=0,bqp=0 ②forv=0ton-1do/*對(duì)兩個(gè)主行其余元素的旋轉(zhuǎn)變換*/ if((v≠b[i])and(v≠b[j]))then br[v]=a[i,v]*cos1+a[j,v]*sin1 a[j,v]=-a[i,v]*sin1+a[j,v]*cos1 endifendfor ③forv=0ton-1do if((v≠b[i])and(v≠b[j]))then a[i,v]=br[v] endifendfor ④forv=0tom-1do/*對(duì)本處理器內(nèi)兩個(gè)主列的其余元素旋轉(zhuǎn)變換*/ if((v≠i)and(v≠j))thenbi[v]=a[v,b[i]]*cos1+a[v,b[j]]*sin1a[v,b[j]]=-a[v,b[i]]*sin1+a[v,b[j]]*cos1 endifendfor ⑤forv=0tom-1do if((v≠i)and(v≠j))then a[v,b[i]]=bi[v] endifendfor ⑥Compute:a[i,b[i]]=bpp,a[j,b[j]]=bqq,a[i,b[j]]=bpq,a[j,b[i]]=bqp/*用temp1保存本處理器主行的行號(hào)和旋轉(zhuǎn)參數(shù)*/temp1[0]=sin1,temp1[1]=cos1,temp1[2]=(float)b[i],temp1[3]=(float)b[j] else ⑦Compute:temp1[0]=0,temp1[1]=0,temp1[2]=0,temp1[3]=0 endif(ii)將所有處理器temp1中的旋轉(zhuǎn)參數(shù)及主行的行號(hào)按處理器編號(hào)連接起來(lái)并廣播給所有處理器,存于temp2中 (iii)current=0 (iv)forv=1topdo/*根據(jù)temp2中的其它處理器的旋轉(zhuǎn)參數(shù)及主行的行號(hào)對(duì)相關(guān)的列在本處理器的部分進(jìn)行旋轉(zhuǎn)變換*/ ①Compute: s1=temp2[(v-1)*4+0],c1=temp2[(v-1)*4+1], i1=(int)temp2[(v-1)*4+2],j1=(int)temp2[(v-1)*4+3] ②if(s1、c1、i1、j1中有一不為0)then if(my-rank≠current)then forz=0tom-1do zi[z]=a[z,i1]*c1+a[z,j1]*s1 a[z,j1]=-a[z,i1]s1+a[z,j1]*c1 endfor forz=0tom-1do a[z,i1]=zi[z] endfor endif endif③current=current+1endforendforendforendfor (3)Data-exchange()/*進(jìn)行一輪中的最后一次處理器間的數(shù)據(jù)交換,使數(shù)據(jù)回到原來(lái)的位置*/ (4)localmax=max/*localmax為本處理器中非對(duì)角元最大的絕對(duì)值*/ (5)fori=0tom-1do forj=0ton-1do if((m*my-rank+i)≠j)then if(│a[i,j]│>localmax)thenlocalmax=│a[i,j]│endif endif endforendfor (6)通過(guò)Allreduce操作求出所有處理器中l(wèi)ocalmax的最大值為新的max值endwhileEnd在上述算法中,每個(gè)處理器在一輪迭代中要處理2p個(gè)行塊對(duì),由于每一個(gè)行塊對(duì)含有m/2行,因而對(duì)每一個(gè)行塊對(duì)的處理要有(m/2)2=m2/4個(gè)行的配對(duì),即消去m2/4個(gè)非主對(duì)角元素.每消去一個(gè)非主對(duì)角元素要對(duì)同行n個(gè)元素和同列n個(gè)元素進(jìn)行旋轉(zhuǎn)變換.由于按行劃分?jǐn)?shù)據(jù)塊,對(duì)同行n個(gè)元素進(jìn)行旋轉(zhuǎn)變換的過(guò)程是在本處理器中進(jìn)行的.而對(duì)同列n個(gè)元素進(jìn)行旋轉(zhuǎn)變換的過(guò)程則分布在所有處理器中進(jìn)行.但由于所有處理器同時(shí)在為其它處理器的元素消去過(guò)程進(jìn)行列的旋轉(zhuǎn)變換，對(duì)每個(gè)處理器而言，對(duì)列元素進(jìn)行旋轉(zhuǎn)變換的處理總量仍然是n個(gè)元素。因此，一個(gè)處理器每消去一個(gè)非主對(duì)角元素共要對(duì)2n個(gè)元素進(jìn)行旋轉(zhuǎn)變換。而對(duì)一個(gè)元素進(jìn)行旋轉(zhuǎn)變換需要2次乘法和1次加法，若取一次乘法運(yùn)算時(shí)間或一次加法運(yùn)算時(shí)間為一個(gè)單位時(shí)間，則其需要3個(gè)單位運(yùn)算時(shí)間，即一輪迭代的計(jì)算時(shí)間為T1＝3×2p×2n×m2/4＝3n3/p；在每輪迭代中，各個(gè)處理器之間以點(diǎn)對(duì)點(diǎn)的通信方式相互錯(cuò)開(kāi)交換數(shù)據(jù)2p-1次，通信量為mn+m，擴(kuò)展收集操作n(n-1)/(2p)次，通信量為4，另外有歸約操作1次，通信量為1，從而不難得出上述算法求解過(guò)程中的總通信時(shí)間為：因此雅可比算法求對(duì)矩陣特征值的一輪并行計(jì)算時(shí)間為