/*
* Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
* Copyright (C) 2008-2008 - INRIA - Bruno JOFRET
*
* This file must be used under the terms of the CeCILL.
* This source file is licensed as described in the file COPYING, which
* you should have received as part of this distribution. The terms
* are also available at
* http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
*
*/
#include <stdlib.h>
#include <stdio.h>
#include <assert.h>
#include <math.h>
#include "matrixMultiplication.h"
int testFloatMultiplication(void);
static void zmulmaTest(void) {
double realM1[4] = {1.0, 2.0, 3.0, 4.0};
double imagM1[4] = {1.0, 2.0, 3.0, 4.0};
double realM3[6] = {1.0, 2.0, 3.0, 4.0, 5.0, 6.0};
double imagM3[6] = {1.0, 2.0, 3.0, 4.0, 5.0, 6.0};
doubleComplex *M1;
doubleComplex *M2;
doubleComplex M1_mul_M2[4];
doubleComplex *M3;
doubleComplex *M4;
doubleComplex M3_mul_M4[4];
doubleComplex miscM3_mul_M4[9];
int i = 0;
printf("\n>>>> Matrix Complex Double Multiplication Tests\n");
M1 = DoubleComplexMatrix(realM1, imagM1, 4);
M2 = DoubleComplexMatrix(realM1, imagM1, 4);
/*
[ 1+1.%i 3+3.%i] * [ 1+1.%i 3+3.%i] = [ 14.%i 30.%i ]
[ 2+2.%i 4+4.%i] [ 2+2.%i 4+4.%i] [ 20.%i 44.%i ]
*/
zmulma(M1, 2, 2, M2, 2, 2, M1_mul_M2);
for (i = 0; i < 4; ++i)
{
printf("M1_mul_M2[%d] = %e + %e i\n", i, zreals(M1_mul_M2[i]), zimags(M1_mul_M2[i]));
}
for (i = 0; i < 4; ++i)
{
assert(zreals(M1_mul_M2[i]) == 0.0);
}
assert(zimags(M1_mul_M2[0]) == 14.0);
assert(zimags(M1_mul_M2[1]) == 20.0);
assert(zimags(M1_mul_M2[2]) == 30.0);
assert(zimags(M1_mul_M2[3]) == 44.0);
M3 = DoubleComplexMatrix(realM3, imagM3, 6);
M4 = DoubleComplexMatrix(realM3, imagM3, 6);
/*
[ 1+1.%i 3+3.%i 5+5.%i ] * [ 1+1.%i 4+4.%i ] = [ 44.%i 98.%i ]
[ 2+2.%i 4+4.%i 6+6.%i ] [ 2+2.%i 5+5.%i ] [ 56.%i 128.%i ]
[ 3+3.%i 6+6.%i ]
*/
zmulma(M3, 2, 3, M4, 3, 2, M3_mul_M4);
for (i = 0; i < 4; ++i)
{
printf("M3_mul_M4[%d] = %e + %e i\n", i, zreals(M3_mul_M4[i]), zimags(M3_mul_M4[i]));
}
for (i = 0; i < 4; ++i)
{
assert(zreals(M3_mul_M4[i]) == 0.0);
}
assert(zimags(M3_mul_M4[0]) == 44.0);
assert(zimags(M3_mul_M4[1]) == 56.0);
assert(zimags(M3_mul_M4[2]) == 98.0);
assert(zimags(M3_mul_M4[3]) == 128.0);
/*
[ 1+1.%i 4+4.%i ] * [ 1+1.%i 3+3.%i 5+5.%i ] = [ 18.%i 38.%i 58.%i ]
[ 2+2.%i 5+5.%i ] [ 2+2.%i 4+4.%i 6+6.%i ] [ 24.%i 52.%i 80.%i ]
[ 3+3.%i 6+6.%i ] [ 30.%i 66.%i 102.%i ]
*/
zmulma(M3, 3, 2, M4, 2, 3, miscM3_mul_M4);
for (i = 0; i < 9; ++i)
{
printf("miscM3_mul_M4[%d] = %e + %e i\n", i, zreals(miscM3_mul_M4[i]), zimags(miscM3_mul_M4[i]));
}
for (i = 0; i < 9; ++i)
{
assert(zreals(miscM3_mul_M4[i]) == 0.0);
}
assert(zimags(miscM3_mul_M4[0]) == 18.0);
assert(zimags(miscM3_mul_M4[1]) == 24.0);
assert(zimags(miscM3_mul_M4[2]) == 30.0);
assert(zimags(miscM3_mul_M4[3]) == 38.0);
assert(zimags(miscM3_mul_M4[4]) == 52.0);
assert(zimags(miscM3_mul_M4[5]) == 66.0);
assert(zimags(miscM3_mul_M4[6]) == 58.0);
assert(zimags(miscM3_mul_M4[7]) == 80.0);
assert(zimags(miscM3_mul_M4[8]) == 102.0);
}
static void dmulmaTest(void) {
double M1[4] = {1.0, 2.0, 3.0, 4.0};
double M2[4] = {1.0, 2.0, 3.0, 4.0};
double M1_by_M2[4];
double M3[4] = {1.0, 0.0, 1.0, 0.0};
double M4[4] = {0.0, 1.0, 0.0, 1.0};
double M3_by_M4[4];
double M5[4] = {1.0, 0.0, 0.0, 1.0};
double M6[4] = {42.0, 51.0, 69.0, 1664.0};
double M5_by_M6[4];
double M7[6] = {1.0, 4.0, 2.0, 5.0, 3.0, 6.0};
double M8[6] = {1.0, 3.0, 5.0, 2.0, 4.0, 6.0};
double M7_by_M8[4];
double miscM7_by_M8[9];
double M9[6] = {1, 4, 2, 5, 3, 6};
double M10[9] = {4, 8, 3, 2, 8, 4, 3, 4, 5};
double M9_by_M10[6];
int i = 0;
printf("\n>>>> Matrix Real Double Multiplication Tests\n");
/*
[ 1 3 ] * [ 1 3 ] = [ 7 15 ]
[ 2 4 ] [ 2 4 ] [10 22 ]
*/
dmulma(M1, 2, 2, M2, 2, 2, M1_by_M2);
for (i = 0; i < 4; ++i) {
printf("M1_by_M2[%d] = %e\n", i, M1_by_M2[i]);
}
assert(M1_by_M2[0] == 7.0);
assert(M1_by_M2[1] == 10.0);
assert(M1_by_M2[2] == 15.0);
assert(M1_by_M2[3] == 22.0);
/*
[ 1 1 ] * [ 0 0 ] = [ 1 1 ]
[ 0 0 ] [ 1 1 ] [ 0 0 ]
*/
dmulma(M3, 2, 2, M4, 2, 2, M3_by_M4);
for (i = 0; i < 4; ++i) {
printf("M3_by_M4[%d] = %e\n", i, M3_by_M4[i]);
}
assert(M3_by_M4[0] == 1.0);
assert(M3_by_M4[1] == 0.0);
assert(M3_by_M4[2] == 1.0);
assert(M3_by_M4[3] == 0.0);
/*
[ 1 0 ] * [ 42 69 ] = [ 42 69 ]
[ 0 1 ] [ 51 1664 ] [ 51 1664 ]
*/
dmulma(M5, 2, 2, M6, 2, 2, M5_by_M6);
for (i = 0; i < 4; ++i) {
printf("M5_by_M6[%d] = %e\n", i, M5_by_M6[i]);
}
assert(M5_by_M6[0] == 42.0);
assert(M5_by_M6[1] == 51.0);
assert(M5_by_M6[2] == 69.0);
assert(M5_by_M6[3] == 1664.0);
/*
[ 1 2 3 ] * [ 1 2 ] = [ 22 28 ]
[ 4 5 6 ] [ 3 4 ] [ 49 64 ]
[ 5 6 ]
*/
dmulma(M7, 2, 3, M8, 3, 2, M7_by_M8);
for (i = 0; i < 4; ++i) {
printf("M7_by_M8[%d] = %e\n", i, M7_by_M8[i]);
}
assert(M7_by_M8[0] == 22.0);
assert(M7_by_M8[1] == 49.0);
assert(M7_by_M8[2] == 28.0);
assert(M7_by_M8[3] == 64.0);
/*
[ 1 5 ] * [ 1 5 4 ] = [ 16 15 34 ]
[ 4 3 ] [ 3 2 6 ] [ 13 26 34 ]
[ 2 6 ] [ 20 22 44 ]
*/
dmulma(M7, 3, 2, M8, 2, 3, miscM7_by_M8);
for (i = 0; i < 9; ++i) {
printf("miscM7_by_M8[%d] = %e\n", i, miscM7_by_M8[i]);
}
assert(miscM7_by_M8[0] == 16.0);
assert(miscM7_by_M8[1] == 13.0);
assert(miscM7_by_M8[2] == 20.0);
assert(miscM7_by_M8[3] == 15.0);
assert(miscM7_by_M8[4] == 26.0);
assert(miscM7_by_M8[5] == 22.0);
assert(miscM7_by_M8[6] == 34.0);
assert(miscM7_by_M8[7] == 34.0);
assert(miscM7_by_M8[8] == 44.0);
/*
[ 1 2 3 ] * [ 4 2 3 ] = [ 29 30 26 ]
[ 4 5 6 ] [ 8 8 4 ] [ 74 72 62 ]
[ 3 4 5 ]
*/
dmulma(M9, 2, 3, M10, 3, 3, M9_by_M10);
for (i = 0; i < 6; ++i) {
printf("M9_by_M10[%d] = %e\n", i, M9_by_M10[i]);
}
assert(M9_by_M10[0] == 29.0);
assert(M9_by_M10[1] == 74.0);
assert(M9_by_M10[2] == 30.0);
assert(M9_by_M10[3] == 72.0);
assert(M9_by_M10[4] == 26.0);
assert(M9_by_M10[5] == 62.0);
}
static void dmulma2Test(void){
int i=0;
double in1[16]={0.2164632631465792655945 , 0.8833887814544141292572 , 0.6525134947150945663452 , 0.3076090742833912372589, 0.9329616213217377662659,
0.2146007861010730266571 , 0.3126419968903064727783 , 0.3616361008025705814362 , 0.2922266637906432151794 , 0.5664248815737664699554 ,
0.4826471973210573196411 , 0.3321718913502991199493 , 0.5935094701126217842102 , 0.5015341597609221935272 , 0.4368587583303451538086 ,
0.2693124809302389621735 };
double in2[16]={ 0.6325744865462183952332 , 0.4051954015158116817474 , 0.9184707831591367721558 , 0.0437334333546459674835 , 0.4818508932366967201233 ,
0.2639556000940501689911 , 0.4148103706538677215576 , 0.2806498021818697452545 , 0.1280058464035391807556 , 0.7783128595910966396332 ,
0.2119030449539422988892 , 0.1121354666538536548615 , 0.6856895955279469490051 , 0.1531216683797538280487 , 0.6970850601792335510254 ,
0.8415518426336348056793 };
double result1[1]={3.4777275993941634268936};
double result2[4]={1.9089008209228131018875 , 1.5406061588610213686223 ,
1.8239702765316110344429, 1.4540285665075025622883};
double result4[16]={0.8093187558996659536348 , 1.1879429718699099360890 , 1.0018471710504197602631 , 0.6579870739173818705581,
0.6383504274229201413959 , 0.8580211304925904336471 , 0.7197492031038768001139 , 0.4570484210841743166753 ,
0.8823217718625468997118 , 0.4563724043650834172325 , 0.4781206002867167681458 , 0.4214295036121353255076 ,
0.9944590770529683210199 , 1.4555038456021862636192 , 1.199379422070573131265 , 0.7244911422701945102887};
double result8[64]={0.2553380379093421193026,0.7883219621670439769545,0.6083298137928665472174,0.3291801751097247485944,0.83065502662064605310 ,
0.3389703173185232287779,0.3747825106430331398855,0.3378859496255101069195,
0.211595258152759591042, 0.8361384907451508974319,0.6204223995507252009674,0.2970569646365784355346,0.8828541978493160691244 ,
0.2190383628278015915036,0.3062578731251097141630,0.3439301521590905075243,
0.1814378811044508321704,0.5751726929430913681784,0.4418116408800580319216,0.2359003381269093035932,0.6062085389345149843976 ,
0.2357883305794587769366,0.2659581411763266567405,0.2453412157151133865529 ,
0.1718045618038900324009,0.5254058588311405486593,0.406124205051207498585, 0.2208234097178256027938,0.5535704713264770759906 ,
0.2297740943495505672178,0.2522914667188620452265,0.2255928995139138970583 ,
0.255152333558011423786, 0.5539346979946898619218,0.4591760625013173724440,0.2979094145476890442836,0.5813605948770583786711 ,
0.4178206413298739541062,0.3800327928667431298671,0.2559009023296956453208 ,
0.0786381978942892051476,0.2507090910879711254111,0.1923914651044380252909,0.1024315495419729910021,0.2642508697778033210923 ,
0.1017143270719664260859,0.1152371518487332324732,0.1068312716501549353154 ,
0.1931728416583101681781,0.6924624191277827245372,0.5213254583603328384811,0.2617870559185208612085,0.7306012370300282166014 ,
0.2239452735136546190908,0.2812679063146140134855,0.2892076880831473406630 ,
0.3968171941968219318397,1.09247302468375151463, 0.8610300469787290911228,0.4939695572954598823401,1.149822596137767938274 ,
0.5716619981842931963456,0.5855773582603787108525,0.4787315376620214779635 };
double result16[256]={0.1369291375410663369472,0.5588092048492155905493,0.4127633888838795339638,0.1945856522217737638592,0.5901677185749256704384 ,
0.1357509820803011191259,0.1977693506556700286936,0.2287617707817625745115,0.1848551318024805323326,0.3583059286285278921547 ,
0.3053103030283391694510,0.2101234636160017044126,0.3754389483168098506916,0.3172577135961544003706,0.2763457047440365799140 ,
0.1703602043449341518766 ,
0.0877099188241009936062,0.3579450719959849647189,0.2643954674855681918899,0.1246417823641658567668,0.3780317587503042031649 ,
0.0869552516898331046002,0.1266810994606728801859,0.1465332850673101428018,0.1184089003682757823555,0.2295127573178283864053 ,
0.1955664249089870132536,0.1345945228879510124909,0.2404873080457204104210,0.203219335238222109652, 0.1770131599873631389475 ,
0.1091241788437475462859 ,
0.1988151828274209720338,0.8113667859364313006410,0.5993145805128781145044,0.2825299473639234504319,0.8568979909927945115555 ,
0.1971045520768189285299,0.287152539732276168571, 0.3321521927227535364580,0.2684016527517738559538,0.5202447045798785918436 ,
0.4432973493130339392998,0.30509017719196090956, 0.5451211078267039766843,0.4606444724966738402117,0.4012420058936001510119 ,
0.2473556452745266642790 ,
0.0094666816925501173080,0.0386336243999785347580,0.0285366554341297228026,0.0134528009494570312266,0.0408016148885166698990 ,
0.0093852291768059120658,0.0136729079348656499399,0.0158155883130832522476,0.0127800753253386294162,0.0247717048087195516881 ,
0.0211078190378471211575,0.0145270172726550076542,0.0259562068565215887939,0.0219338107509826525654,0.0191053333928335386527 ,
0.0117779594363369693316 ,
0.1043030167001093577728,0.4256616734190865058984,0.3144142102774669100818,0.1482217072111654543321,0.4495483905894361309485 ,
0.1034055804720993287926,0.1506468254648987425970,0.174254678198354717722, 0.1408096789751012989367,0.272932335137809511938 ,
0.2325639831473396978989,0.1600573225552646139391,0.2859830683182053578939,0.2416646828695165061252,0.210500782919751050493 ,
0.1297684594960265003216 ,
0.0571366905221716223084,0.2331754159251516189677,0.1722345910669886293043,0.0811951377968477827629,0.2462604446206972941269 ,
0.0566450792759636348461,0.0825236059037830066432,0.0954558740030149316391,0.0771348644043414743976,0.1495110195240048345244 ,
0.1273974306025911218399,0.0876786309157438353834,0.1566601483450788256935,0.1323827501073594481440,0.1153113157114278863880 ,
0.0710865375167586688310 ,
0.0897912064187782466007,0.3664388278665740572571,0.2706693645994189267867,0.1275994341199866433101,0.3870021559463034388493 ,
0.08901863162519749184, 0.1296871426120334014343,0.1500104050157337609583,0.121218650721939893988, 0.2349589150731871967093 ,
0.2002070628157982257189,0.1377883453718138639932,0.2461938832839972723932,0.2080415707060042362908,0.1812135434663988964310 ,
0.1117136100363851081152 ,
0.0607503719817294868255,0.2479228867648642919974,0.1831277832127918003824,0.0863304258469818203148,0.2618354944672221518154 ,
0.0602276681673398900374,0.0877429145810092459001,0.1014931001520641418168,0.0820133553851117719224,0.1589670309645665469134 ,
0.1354548404517885995535,0.0932239755978389733615,0.1665683153801736371413,0.1407554627243530753766,0.1226043241068286160367 ,
0.0755824944981801360999 ,
0.0277085632143499079050,0.1130789286734633775078,0.0835255421807369752285,0.0393757599150546533373,0.1194245419993072432341 ,
0.0274701552637327192641,0.0400200034332363480116,0.0462915351733086607999,0.0374067214402037592769,0.0725056963898744188946 ,
0.0617816630073779321508,0.0425199441037594910719,0.0759726820702821964426,0.0641993046204846901093,0.0559204751188750043278 ,
0.0344735720685122454254 ,
0.1684761413360341408829,0.6875528486244794068583,0.5078596439934851547449,0.2394160982416763117087,0.7261360273796675368629 ,
0.1670265515008234058492,0.2433332866279651618413,0.2814660277470227844887,0.2274437703436615054731,0.440855769321226442781 ,
0.3756505203205803833377,0.2585336546326343909463,0.4619360528777511909482,0.3903504860661412556588,0.3400127894335067679776 ,
0.209609367156386977582 ,
0.0458692245814266291726,0.1871927726683430270871,0.1382695964036666780306,0.0651832994961140249623,0.1976974083832430950647 ,
0.0454745600243270306495,0.0662497911215368939786,0.0766317909253355217247,0.0619237198739692532024,0.1200271571431572731337 ,
0.1022744107508182648259,0.0703882352252384740909,0.1257664639258653793696,0.1062766156017563717873,0.0925717011049985438742 ,
0.0570681347532181540427 ,
0.0242732090263575867040,0.0990592132451698637041,0.0731699052278139971950,0.0344938870917280079875,0.1046180867770489658630 ,
0.0240643592937276537547,0.0350582562168871778030,0.0405522329223763147721,0.0327689733128625759750,0.0635163184196280983240 ,
0.0541218687007713472603,0.0372482500458589660397,0.066553461394560253783, 0.0562397670476394023753,0.0489873607271963315557 ,
0.0301994807248194101212 ,
0.1484266073536374686004,0.6057304962494031164155,0.4474217142677203162826,0.2109243417261047193190,0.6397220767671999430704 ,
0.1471495262216242216891,0.2143753643927638896649,0.247970111687618482987, 0.2003767828970874909089,0.3883916479432811819095 ,
0.3309461615137729850744,0.2277668098257397411022,0.4069632685035297625653,0.3438967551499154673422,0.2995495053023755072097 ,
0.1846647661196834933062 ,
0.0331452159959298997549,0.1352659640442576416408,0.0999139549510789648501,0.0471016146630244963989,0.1428566399910645989824 ,
0.0328600304034029899736,0.0478722641694215370056,0.0553743231012384373724,0.0447462343046726856999,0.0867319228783796081217 ,
0.0739037440926125499541,0.0508627141924161210729,0.0908791602628282874088,0.0767957472920304112796,0.0668925419218501238205 ,
0.0412375763955288299201 ,
0.1508933068171264824109,0.6157971218818100167525,0.4548574087312337033140,0.2144296900585161147479,0.6503536079439788952783 ,
0.1495950018937773073890,0.2179380652168350362974,0.2520911230909432809710,0.2037068415144771826864,0.3948463226588642127624 ,
0.3364461505898876381160,0.2315520628717731010937,0.4137265846904019661601,0.3496119699388837709364,0.3045277138405338979155 ,
0.1877337069762743115842 ,
0.1821650579634931654827,0.7434174567948435008447,0.5491239338208003806230,0.2588689832740145302381,0.7851355715297718695922 ,
0.1805976869739845391827,0.2631044485676965472010,0.3043355269932461748539,0.2459238873796954771844,0.4766759028019415533883 ,
0.4061726382474953322976,0.2795398672369437398366,0.4994689881937890429953,0.4220669962905158612010,0.3676392930435437023107 ,
0.2266404145710782247480 };
double out1[1],out2[4],out4[16],out8[64],out16[256];
dmulma(in1, 1, 16, in2, 16, 1, out1);
dmulma(in1, 2, 8, in2, 8, 2, out2);
dmulma(in1, 4, 4, in2, 4, 4, out4);
dmulma(in1, 8, 2, in2, 2, 8, out8);
dmulma(in1, 16, 1, in2, 1, 16, out16);
assert( (fabs(out1[0]-result1[0]) / fabs(out1[0])) <1e-16);
for (i=0;i<4;i++) assert( (fabs(out2[i]-result2[i]) / fabs(out2[i])) <1e-15);
for (i=0;i<16;i++) assert( (fabs(out4[i]-result4[i]) / fabs(out4[i])) <1e-15);
for (i=0;i<64;i++) assert( (fabs(out8[i]-result8[i]) / fabs(out8[i])) <1e-15);
for (i=0;i<256;i++) assert( (fabs(out16[i]-result16[i]) / fabs(out16[i])) <1e-16);
}
static int testDoubleMultiplication(void) {
printf("\n>>>> Matrix Double Multiplication Tests\n");
dmulmaTest();
dmulma2Test();
zmulmaTest();
return 0;
}
int main(void) {
assert(testFloatMultiplication() == 0);
assert(testDoubleMultiplication() == 0);
return 0;
}