C語(yǔ)言科學(xué)計(jì)算入門之矩陣乘法的相關(guān)計(jì)算
更新時(shí)間：2015年12月02日 10:58:39 作者：Denlee
這篇文章主要介紹了C語(yǔ)言科學(xué)計(jì)算入門之矩陣乘法的相關(guān)計(jì)算,文章中還介紹了矩陣相關(guān)的斯特拉森算法的實(shí)現(xiàn),需要的朋友可以參考下
1.矩陣相乘
矩陣相乘應(yīng)滿足的條件：
(1) 矩陣A的列數(shù)必須等于矩陣B的行數(shù)，矩陣A與矩陣B才能相乘；
(2) 矩陣C的行數(shù)等于矩陣A的行數(shù)，矩陣C的列數(shù)等于矩陣B的列數(shù)；
(3) 矩陣C中第i行第j列的元素等于矩陣A的第i行元素與矩陣B的第j列元素對(duì)應(yīng)乘積之和，即
如：
則：
2. 常用矩陣相乘算法
用A的第i行分別和B的第j列的各個(gè)元素相乘求和，求得C的第i行j列的元素，這種算法中，B的訪問(wèn)是按列進(jìn)行訪問(wèn)的，代碼如下：
void arymul(int a[4][5], int b[5][3], int c[4][3])
{
int i, j, k;
int temp;
for(i = 0; i < 4; i++){
for(j = 0; j < 3; j++){
temp = 0;
for(k = 0; k < 5; k++){
temp += a[i][k] * b[k][j];
}
c[i][j] = temp;
printf("%d/t", c[i][j]);
}
printf("%d/n");
}
}
3. 改進(jìn)的算法
矩陣A、B、C都按行（數(shù)據(jù)的存儲(chǔ)順序）訪問(wèn)，以提高存儲(chǔ)器訪問(wèn)效率，對(duì)于A的第i行中，第j列的元素分別和B的第j行的元素相乘，對(duì)于B中相同的列k在上述計(jì)算過(guò)程中求和，從而得到C第i行k列的數(shù)據(jù)，代碼如下：
void arymul1(int a[4][5], int b[5][3], int c[4][3])
{
int i, j, k;
int temp[3] = {0};
for(i = 0; i < 4; i++){
for(k = 0; k < 3; k ++)
temp[k] = 0;
for(j = 0; j < 5; j++){//當(dāng)前行的每個(gè)元素
for(k = 0; k < 3; k++){
temp[k] += a[i][j] * b[j][k];
}
}
for(k = 0; k < 3; k++){
c[i][k] = temp[k];
printf("%d/t", c[i][k]);
}
printf("%d/n");
}
}
這種算法很容易轉(zhuǎn)到稀疏矩陣的相乘算法。
PS：斯特拉森算法的實(shí)現(xiàn)
斯特拉森方法，是由v.斯特拉森在1969年提出的一個(gè)方法。
我們先討論二階矩陣的計(jì)算方法。
對(duì)于二階矩陣
a11 a12 b11 b12 A = a21 a22 B = b21 b22
先計(jì)算下面7個(gè)量(1)
x1 = (a11 + a22) * (b11 + b22); x2 = (a21 + a22) * b11; x3 = a11 * (b12 - b22); x4 = a22 * (b21 - b11); x5 = (a11 + a12) * b22; x6 = (a21 - a11) * (b11 + b12); x7 = (a12 - a22) * (b21 + b22);
再設(shè)C = AB。根據(jù)矩陣相乘的規(guī)則，C的各元素為(2)
c11 = a11 * b11 + a12 * b21 c12 = a11 * b12 + a12 * b22 c21 = a21 * b11 + a22 * b21 c22 = a21 * b12 + a22 * b22
比較(1)(2)，C的各元素可以表示為(3)
c11 = x1 + x4 - x5 + x7 c12 = x3 + x5 c21 = x2 + x4 c22 = x1 + x3 - x2 + x6
根據(jù)以上的方法，我們就可以計(jì)算4階矩陣了，先將4階矩陣A和B劃分成四塊2階矩陣，分別利用公式計(jì)算它們的乘積，再使用(1)(3)來(lái)計(jì)算出最后結(jié)果。
ma11 ma12 mb11 mb12 A4 = ma21 ma22 B4 = mb21 mb22
其中
a11 a12 a13 a14 b11 b12 b13 b14 ma11 = a21 a22 ma12 = a23 a24 mb11 = b21 b22 mb12 = b23 b24 a31 a32 a33 a34 b31 b32 b33 b34 ma21 = a41 a42 ma22 = a43 a44 mb21 = b41 b42 mb22 = b43 b44
實(shí)現(xiàn)
// 計(jì)算2X2矩陣
void Multiply2X2(float& fOut_11, float& fOut_12, float& fOut_21, float& fOut_22,
float f1_11, float f1_12, float f1_21, float f1_22,
float f2_11, float f2_12, float f2_21, float f2_22)
{
const float x1((f1_11 + f1_22) * (f2_11 + f2_22));
const float x2((f1_21 + f1_22) * f2_11);
const float x3(f1_11 * (f2_12 - f2_22));
const float x4(f1_22 * (f2_21 - f2_11));
const float x5((f1_11 + f1_12) * f2_22);
const float x6((f1_21 - f1_11) * (f2_11 + f2_12));
const float x7((f1_12 - f1_22) * (f2_21 + f2_22));
fOut_11 = x1 + x4 - x5 + x7;
fOut_12 = x3 + x5;
fOut_21 = x2 + x4;
fOut_22 = x1 - x2 + x3 + x6;
}
// 計(jì)算4X4矩陣
void Multiply(CLAYMATRIX& mOut, const CLAYMATRIX& m1, const CLAYMATRIX& m2)
{
float fTmp[7][4];
// (ma11 + ma22) * (mb11 + mb22)
Multiply2X2(fTmp[0][0], fTmp[0][1], fTmp[0][2], fTmp[0][3],
m1._11 + m1._33, m1._12 + m1._34, m1._21 + m1._43, m1._22 + m1._44,
m2._11 + m2._33, m2._12 + m2._34, m2._21 + m2._43, m2._22 + m2._44);
// (ma21 + ma22) * mb11
Multiply2X2(fTmp[1][0], fTmp[1][1], fTmp[1][2], fTmp[1][3],
m1._31 + m1._33, m1._32 + m1._34, m1._41 + m1._43, m1._42 + m1._44,
m2._11, m2._12, m2._21, m2._22);
// ma11 * (mb12 - mb22)
Multiply2X2(fTmp[2][0], fTmp[2][1], fTmp[2][2], fTmp[2][3],
m1._11, m1._12, m1._21, m1._22,
m2._13 - m2._33, m2._14 - m2._34, m2._23 - m2._43, m2._24 - m2._44);
// ma22 * (mb21 - mb11)
Multiply2X2(fTmp[3][0], fTmp[3][1], fTmp[3][2], fTmp[3][3],
m1._33, m1._34, m1._43, m1._44,
m2._31 - m2._11, m2._32 - m2._12, m2._41 - m2._21, m2._42 - m2._22);
// (ma11 + ma12) * mb22
Multiply2X2(fTmp[4][0], fTmp[4][1], fTmp[4][2], fTmp[4][3],
m1._11 + m1._13, m1._12 + m1._14, m1._21 + m1._23, m1._22 + m1._24,
m2._33, m2._34, m2._43, m2._44);
// (ma21 - ma11) * (mb11 + mb12)
Multiply2X2(fTmp[5][0], fTmp[5][1], fTmp[5][2], fTmp[5][3],
m1._31 - m1._11, m1._32 - m1._12, m1._41 - m1._21, m1._42 - m1._22,
m2._11 + m2._13, m2._12 + m2._14, m2._21 + m2._23, m2._22 + m2._24);
// (ma12 - ma22) * (mb21 + mb22)
Multiply2X2(fTmp[6][0], fTmp[6][1], fTmp[6][2], fTmp[6][3],
m1._13 - m1._33, m1._14 - m1._34, m1._23 - m1._43, m1._24 - m1._44,
m2._31 + m2._33, m2._32 + m2._34, m2._41 + m2._43, m2._42 + m2._44);
// 第一塊
mOut._11 = fTmp[0][0] + fTmp[3][0] - fTmp[4][0] + fTmp[6][0];
mOut._12 = fTmp[0][1] + fTmp[3][1] - fTmp[4][1] + fTmp[6][1];
mOut._21 = fTmp[0][2] + fTmp[3][2] - fTmp[4][2] + fTmp[6][2];
mOut._22 = fTmp[0][3] + fTmp[3][3] - fTmp[4][3] + fTmp[6][3];
// 第二塊
mOut._13 = fTmp[2][0] + fTmp[4][0];
mOut._14 = fTmp[2][1] + fTmp[4][1];
mOut._23 = fTmp[2][2] + fTmp[4][2];
mOut._24 = fTmp[2][3] + fTmp[4][3];
// 第三塊
mOut._31 = fTmp[1][0] + fTmp[3][0];
mOut._32 = fTmp[1][1] + fTmp[3][1];
mOut._41 = fTmp[1][2] + fTmp[3][2];
mOut._42 = fTmp[1][3] + fTmp[3][3];
// 第四塊
mOut._33 = fTmp[0][0] - fTmp[1][0] + fTmp[2][0] + fTmp[5][0];
mOut._34 = fTmp[0][1] - fTmp[1][1] + fTmp[2][1] + fTmp[5][1];
mOut._43 = fTmp[0][2] - fTmp[1][2] + fTmp[2][2] + fTmp[5][2];
mOut._44 = fTmp[0][3] - fTmp[1][3] + fTmp[2][3] + fTmp[5][3];
}
比較
在標(biāo)準(zhǔn)的定義算法中我們需要進(jìn)行n * n * n次乘法運(yùn)算，新算法中我們需要進(jìn)行7log2n次乘法，對(duì)于最常用的4階矩陣： 原算法 新算法
加法次數(shù) 48 72(48次加法，24次減法)
乘法次數(shù) 64 49
需要額外空間 16 * sizeof(float) 28 * sizeof(float)
新算法要比原算法多了24次減法運(yùn)算，少了15次乘法。但因?yàn)楦↑c(diǎn)乘法的運(yùn)算速度要遠(yuǎn)遠(yuǎn)慢于加/減法運(yùn)算，所以新算法的整體速度有所提高。
相關(guān)文章
C++詳細(xì)講解模擬實(shí)現(xiàn)位圖和布隆過(guò)濾器的方法
位圖（bitset）是一種常用的數(shù)據(jù)結(jié)構(gòu)，常用在給一個(gè)很大范圍的數(shù)，判斷其中的一個(gè)數(shù)是不是在其中。在索引、數(shù)據(jù)壓縮方面有很大的應(yīng)用。布隆過(guò)濾器是由布隆提出的，它實(shí)際上是一個(gè)很長(zhǎng)的二進(jìn)制向量和一系列隨機(jī)映射函數(shù)。布隆過(guò)濾器可以用于檢索一個(gè)元素是否在一個(gè)集合中2022-06-06
這篇文章主要介紹了C++之拼接長(zhǎng)字符串問(wèn)題,具有很好的參考價(jià)值,希望對(duì)大家有所幫助,如有錯(cuò)誤或未考慮完全的地方,望不吝賜教2023-08-08
C/C++多參數(shù)函數(shù)參數(shù)的計(jì)算順序與壓棧順序的示例代碼
這篇文章主要介紹了C/C++多參數(shù)函數(shù)參數(shù)的計(jì)算順序與壓棧順序,本文通過(guò)實(shí)例代碼給大家介紹的非常詳細(xì)，對(duì)大家的學(xué)習(xí)或工作具有一定的參考借鑒價(jià)值，需要的朋友可以參考下2020-06-06
QWebEngineView是Qt框架中的一個(gè)組件,它是基于Chromium內(nèi)核的Web瀏覽器引擎,用于在Qt應(yīng)用程序中嵌入網(wǎng)頁(yè)內(nèi)容和實(shí)現(xiàn)各種Web應(yīng)用功能,這篇文章主要給大家介紹了關(guān)于Qt與QWebEngineView交互完整參考的相關(guān)資料,需要的朋友可以參考下2024-07-07
C語(yǔ)言實(shí)現(xiàn)隨機(jī)生成6位數(shù)密碼
這篇文章主要為大家詳細(xì)介紹了如何使用C語(yǔ)言實(shí)現(xiàn)一個(gè)簡(jiǎn)單而實(shí)用的隨機(jī)密碼生成器,該生成器將生成包含字母、數(shù)字和特殊字符的隨機(jī)密碼,有需要的小伙伴可以參考下2023-11-11
簡(jiǎn)單了解設(shè)計(jì)模式中的裝飾者模式及C++版代碼實(shí)現(xiàn)
這篇文章主要介紹了簡(jiǎn)單了解設(shè)計(jì)模式中的裝飾者模式及C++版代碼實(shí)現(xiàn),ConcreteComponent的引用(指針)也可以達(dá)到修飾的功能,需要的朋友可以參考下2016-03-03
這篇文章主要介紹了C語(yǔ)言中 printf 函數(shù)簡(jiǎn)介,通過(guò)實(shí)例代碼給大家介紹Printf輸出格式的相關(guān)知識(shí)，需要的朋友可以參考下2021-08-08
這篇文章介紹了C++引用使用與const修飾符,引用是給已經(jīng)定義的變量一個(gè)別名,可以簡(jiǎn)單理解成同一個(gè)變量的昵稱,既然是昵稱或者是別名,顯然它和原本的變量名有著同樣的效力,所以我們對(duì)別名進(jìn)行修改,原本的變量值也一樣會(huì)發(fā)生變化,下面來(lái)看看詳細(xì)內(nèi)容,需要的朋友可以參考下2021-11-11
C++使用遞歸和非遞歸算法實(shí)現(xiàn)的二叉樹葉子節(jié)點(diǎn)個(gè)數(shù)計(jì)算方法
這篇文章主要介紹了C++使用遞歸和非遞歸算法實(shí)現(xiàn)的二叉樹葉子節(jié)點(diǎn)個(gè)數(shù)計(jì)算方法,涉及C++二叉樹的定義、遍歷、統(tǒng)計(jì)相關(guān)操作技巧,需要的朋友可以參考下2017-05-05
C語(yǔ)言不用鏈表完成學(xué)生管理系統(tǒng)(完整代碼)
這篇文章主要介紹了C語(yǔ)言不用鏈表完成學(xué)生管理系統(tǒng)(完整代碼),本文通過(guò)實(shí)例代碼給大家介紹的非常詳細(xì)，對(duì)大家的學(xué)習(xí)或工作具有一定的參考借鑒價(jià)值，需要的朋友可以參考下2021-04-04