Modify arrays pow : now computes A.^B instead of A.^b, A,B=matrices,b=scalar

10 files changed, 476 insertions, 154 deletions

diff --git a/src/matrixOperations/includes/matrixPow.h b/src/matrixOperations/includes/matrixPow.h
index d9d5684b..71a2d5d0 100644
--- a/src/matrixOperations/includes/matrixPow.h
+++ b/src/matrixOperations/includes/matrixPow.h
@@ -21,9 +21,9 @@
    powm is only working on square matrix
    so the size is limited to rows
 */
-void spowma(float* in, int rows, float expand, float* out);
+void spowma(float* in, int rows, float expand, floatComplex* out);
 
-void dpowma(double* in, int rows, double expand, double* out);
+void dpowma(double* in, int rows, double expand, doubleComplex* out);
 
 void cpowma(floatComplex* in, int rows, floatComplex expand, floatComplex* out);
 
diff --git a/src/matrixOperations/logm/zlogma.c b/src/matrixOperations/logm/zlogma.c
index 3a629bfd..ea0cd4e7 100644
--- a/src/matrixOperations/logm/zlogma.c
+++ b/src/matrixOperations/logm/zlogma.c
@@ -25,9 +25,8 @@
 void zlogma (doubleComplex* in, int rows, doubleComplex* out){
 
 	/* Algo : */
-	/* trouver les valeurs propres vp */
-	/* en déduire les vecteurs propres Vp */
-	/* utiliser la formule suivante
+	/* find eigenvalues and eigenvectors */
+	/* use formula
 	*	logm = Vp * diag(log(diag(vp)) * inv(Vp) */
 	
 	
@@ -45,7 +44,7 @@ void zlogma (doubleComplex* in, int rows, doubleComplex* out){
 	
 	
 	
-	/* regarde si in est hermitienne */
+	/* hermitian test */
 	for (i=0;i<rows;i++){
 		for(j=0;j<rows;j++)
 			if ( (zreals(in[i+j*rows])!=zreals(in[j+i*rows])) ||  
@@ -54,36 +53,26 @@ void zlogma (doubleComplex* in, int rows, doubleComplex* out){
 	}
 	if ((i==rows) && (j==rows)) hermitienne=1;
 		
-	/* trouver les valeurs propres vp ainsi que les vecteurs propres*/
+	/* find eigenvalues and eigenvectors */
 	zspec2a(in,rows,eigenvalues,eigenvectors);
 
 
-	/* utiliser la formule suivante
-	*	logm = Vp * diag(log(diag(vp)) * inv(Vp) */
-	
-
-	/* diag(log(diag(vp)) */
-	for (i=0;i<rows*rows;i++){
-		if ((i%(rows+1))==0) /* teste si i est sur la diagonale */
-			eigenvalues[i] = zlogs(eigenvalues[i]);
-		else    eigenvalues[i] = DoubleComplex(0,0);
+	/* make operation on eigenvalues and eigenvectors */
+	for (i=0;i<rows;i++){
+		eigenvalues[i+i*rows] = zlogs(eigenvalues[i+i*rows]);
 	}
 	
 	
-	
-	/* Vp * diag(log(diag(vp)) */
 	zmulma(eigenvectors, rows, rows, eigenvalues, rows, rows, tmp);
 	
-	/* Vp' ou inv(Vp) */
 	if (hermitienne) {
 		/* we use eigenvalues as a temporary matrix cause it's useless now*/
 		ztransposea(eigenvectors,rows,rows,eigenvalues);
-		zconja(eigenvalues,rows*rows,eigenvectors);
+		zconja(eigenvalues,rows*rows,eigenvalues);
 	}
-	else	zinverma(eigenvectors, eigenvectors, rows);
+	else	zinverma(eigenvectors, eigenvalues, rows);
 	
-	/* Vp * diag(log(diag(vp))*inv(Vp) */
-	zmulma(tmp, rows, rows, eigenvectors, rows, rows, out);
+	zmulma(tmp, rows, rows, eigenvalues, rows, rows, out);
 	
 
 	free(eigenvalues);
diff --git a/src/matrixOperations/multiplication/testDoubleMatrixMultiplication.c b/src/matrixOperations/multiplication/testDoubleMatrixMultiplication.c
index d1c210ea..ecee5428 100644
--- a/src/matrixOperations/multiplication/testDoubleMatrixMultiplication.c
+++ b/src/matrixOperations/multiplication/testDoubleMatrixMultiplication.c
@@ -214,10 +214,133 @@ static void dmulmaTest(void) {
 
 }
 
+
+
+
+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;
 }
diff --git a/src/matrixOperations/powm/Makefile.am b/src/matrixOperations/powm/Makefile.am
index 4215c21e..d3ce9bc7 100644
--- a/src/matrixOperations/powm/Makefile.am
+++ b/src/matrixOperations/powm/Makefile.am
@@ -30,7 +30,7 @@ libMatrixPow_la_SOURCES = 	$(HEAD) \
 				cpowma.c\
 				zpowma.c
 
-check_PROGRAMS = 	testDoubleMatrixPow
+check_PROGRAMS = 	testDoubleMatrixPow	testFloatMatrixPow
 
 check_LDADD =	$(top_builddir)/type/libDoubleComplex.la \
 			$(top_builddir)/type/libFloatComplex.la \
@@ -50,6 +50,8 @@ check_LDADD =	$(top_builddir)/type/libDoubleComplex.la \
 			$(top_builddir)/elementaryFunctions/cosh/libCosh.la \
 			$(top_builddir)/elementaryFunctions/sin/libSin.la \
 			$(top_builddir)/elementaryFunctions/sinh/libSinh.la \
+			$(top_builddir)/matrixOperations/spec2/libSpec2.la \
+			$(top_builddir)/matrixOperations/zeros/libMatrixZeros.la \
 			$(top_builddir)/matrixOperations/inversion/libMatrixInversion.la \
 			$(top_builddir)/matrixOperations/transpose/libMatrixTranspose.la \
 			$(top_builddir)/matrixOperations/cat/libMatrixConcatenation.la \
@@ -77,4 +79,8 @@ testDoubleMatrixPow_SOURCES = testDoublePowm.c
 testDoubleMatrixPow_LDADD = $(check_LDADD)
 testDoubleMatrixPow_CFLAGS = $(check_INCLUDES)
 
-TESTS =		testDoubleMatrixPow
+testFloatMatrixPow_SOURCES = testFloatPowm.c
+testFloatMatrixPow_LDADD = $(check_LDADD)
+testFloatMatrixPow_CFLAGS = $(check_INCLUDES)
+
+TESTS =	testDoubleMatrixPow	testFloatMatrixPow
diff --git a/src/matrixOperations/powm/Makefile.in b/src/matrixOperations/powm/Makefile.in
index 32b7f77b..0f801bb7 100644
--- a/src/matrixOperations/powm/Makefile.in
+++ b/src/matrixOperations/powm/Makefile.in
@@ -32,8 +32,9 @@ PRE_UNINSTALL = :
 POST_UNINSTALL = :
 build_triplet = @build@
 host_triplet = @host@
-check_PROGRAMS = testDoubleMatrixPow$(EXEEXT)
-TESTS = testDoubleMatrixPow$(EXEEXT)
+check_PROGRAMS = testDoubleMatrixPow$(EXEEXT) \
+	testFloatMatrixPow$(EXEEXT)
+TESTS = testDoubleMatrixPow$(EXEEXT) testFloatMatrixPow$(EXEEXT)
 subdir = matrixOperations/powm
 DIST_COMMON = $(srcdir)/Makefile.am $(srcdir)/Makefile.in
 ACLOCAL_M4 = $(top_srcdir)/aclocal.m4
@@ -69,6 +70,14 @@ testDoubleMatrixPow_LINK = $(LIBTOOL) --tag=CC $(AM_LIBTOOLFLAGS) \
 	$(LIBTOOLFLAGS) --mode=link $(CCLD) \
 	$(testDoubleMatrixPow_CFLAGS) $(CFLAGS) $(AM_LDFLAGS) \
 	$(LDFLAGS) -o $@
+am_testFloatMatrixPow_OBJECTS =  \
+	testFloatMatrixPow-testFloatPowm.$(OBJEXT)
+testFloatMatrixPow_OBJECTS = $(am_testFloatMatrixPow_OBJECTS)
+testFloatMatrixPow_DEPENDENCIES = $(check_LDADD)
+testFloatMatrixPow_LINK = $(LIBTOOL) --tag=CC $(AM_LIBTOOLFLAGS) \
+	$(LIBTOOLFLAGS) --mode=link $(CCLD) \
+	$(testFloatMatrixPow_CFLAGS) $(CFLAGS) $(AM_LDFLAGS) \
+	$(LDFLAGS) -o $@
 DEFAULT_INCLUDES = -I.@am__isrc@ -I$(top_builddir)/includes
 depcomp = $(SHELL) $(top_srcdir)/config/depcomp
 am__depfiles_maybe = depfiles
@@ -81,9 +90,10 @@ CCLD = $(CC)
 LINK = $(LIBTOOL) --tag=CC $(AM_LIBTOOLFLAGS) $(LIBTOOLFLAGS) \
 	--mode=link $(CCLD) $(AM_CFLAGS) $(CFLAGS) $(AM_LDFLAGS) \
 	$(LDFLAGS) -o $@
-SOURCES = $(libMatrixPow_la_SOURCES) $(testDoubleMatrixPow_SOURCES)
+SOURCES = $(libMatrixPow_la_SOURCES) $(testDoubleMatrixPow_SOURCES) \
+	$(testFloatMatrixPow_SOURCES)
 DIST_SOURCES = $(libMatrixPow_la_SOURCES) \
-	$(testDoubleMatrixPow_SOURCES)
+	$(testDoubleMatrixPow_SOURCES) $(testFloatMatrixPow_SOURCES)
 ETAGS = etags
 CTAGS = ctags
 DISTFILES = $(DIST_COMMON) $(DIST_SOURCES) $(TEXINFOS) $(EXTRA_DIST)
@@ -231,6 +241,8 @@ check_LDADD = $(top_builddir)/type/libDoubleComplex.la \
 			$(top_builddir)/elementaryFunctions/cosh/libCosh.la \
 			$(top_builddir)/elementaryFunctions/sin/libSin.la \
 			$(top_builddir)/elementaryFunctions/sinh/libSinh.la \
+			$(top_builddir)/matrixOperations/spec2/libSpec2.la \
+			$(top_builddir)/matrixOperations/zeros/libMatrixZeros.la \
 			$(top_builddir)/matrixOperations/inversion/libMatrixInversion.la \
 			$(top_builddir)/matrixOperations/transpose/libMatrixTranspose.la \
 			$(top_builddir)/matrixOperations/cat/libMatrixConcatenation.la \
@@ -257,6 +269,9 @@ check_INCLUDES = -I $(top_builddir)/type \
 testDoubleMatrixPow_SOURCES = testDoublePowm.c
 testDoubleMatrixPow_LDADD = $(check_LDADD)
 testDoubleMatrixPow_CFLAGS = $(check_INCLUDES)
+testFloatMatrixPow_SOURCES = testFloatPowm.c
+testFloatMatrixPow_LDADD = $(check_LDADD)
+testFloatMatrixPow_CFLAGS = $(check_INCLUDES)
 all: all-am
 
 .SUFFIXES:
@@ -329,6 +344,9 @@ clean-checkPROGRAMS:
 testDoubleMatrixPow$(EXEEXT): $(testDoubleMatrixPow_OBJECTS) $(testDoubleMatrixPow_DEPENDENCIES) 
 	@rm -f testDoubleMatrixPow$(EXEEXT)
 	$(testDoubleMatrixPow_LINK) $(testDoubleMatrixPow_OBJECTS) $(testDoubleMatrixPow_LDADD) $(LIBS)
+testFloatMatrixPow$(EXEEXT): $(testFloatMatrixPow_OBJECTS) $(testFloatMatrixPow_DEPENDENCIES) 
+	@rm -f testFloatMatrixPow$(EXEEXT)
+	$(testFloatMatrixPow_LINK) $(testFloatMatrixPow_OBJECTS) $(testFloatMatrixPow_LDADD) $(LIBS)
 
 mostlyclean-compile:
 	-rm -f *.$(OBJEXT)
@@ -341,6 +359,7 @@ distclean-compile:
 @AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/libMatrixPow_la-spowma.Plo@am__quote@
 @AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/libMatrixPow_la-zpowma.Plo@am__quote@
 @AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/testDoubleMatrixPow-testDoublePowm.Po@am__quote@
+@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/testFloatMatrixPow-testFloatPowm.Po@am__quote@
 
 .c.o:
 @am__fastdepCC_TRUE@	$(COMPILE) -MT $@ -MD -MP -MF $(DEPDIR)/$*.Tpo -c -o $@ $<
@@ -405,6 +424,20 @@ testDoubleMatrixPow-testDoublePowm.obj: testDoublePowm.c
 @AMDEP_TRUE@@am__fastdepCC_FALSE@	DEPDIR=$(DEPDIR) $(CCDEPMODE) $(depcomp) @AMDEPBACKSLASH@
 @am__fastdepCC_FALSE@	$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(testDoubleMatrixPow_CFLAGS) $(CFLAGS) -c -o testDoubleMatrixPow-testDoublePowm.obj `if test -f 'testDoublePowm.c'; then $(CYGPATH_W) 'testDoublePowm.c'; else $(CYGPATH_W) '$(srcdir)/testDoublePowm.c'; fi`
 
+testFloatMatrixPow-testFloatPowm.o: testFloatPowm.c
+@am__fastdepCC_TRUE@	$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(testFloatMatrixPow_CFLAGS) $(CFLAGS) -MT testFloatMatrixPow-testFloatPowm.o -MD -MP -MF $(DEPDIR)/testFloatMatrixPow-testFloatPowm.Tpo -c -o testFloatMatrixPow-testFloatPowm.o `test -f 'testFloatPowm.c' || echo '$(srcdir)/'`testFloatPowm.c
+@am__fastdepCC_TRUE@	mv -f $(DEPDIR)/testFloatMatrixPow-testFloatPowm.Tpo $(DEPDIR)/testFloatMatrixPow-testFloatPowm.Po
+@AMDEP_TRUE@@am__fastdepCC_FALSE@	source='testFloatPowm.c' object='testFloatMatrixPow-testFloatPowm.o' libtool=no @AMDEPBACKSLASH@
+@AMDEP_TRUE@@am__fastdepCC_FALSE@	DEPDIR=$(DEPDIR) $(CCDEPMODE) $(depcomp) @AMDEPBACKSLASH@
+@am__fastdepCC_FALSE@	$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(testFloatMatrixPow_CFLAGS) $(CFLAGS) -c -o testFloatMatrixPow-testFloatPowm.o `test -f 'testFloatPowm.c' || echo '$(srcdir)/'`testFloatPowm.c
+
+testFloatMatrixPow-testFloatPowm.obj: testFloatPowm.c
+@am__fastdepCC_TRUE@	$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(testFloatMatrixPow_CFLAGS) $(CFLAGS) -MT testFloatMatrixPow-testFloatPowm.obj -MD -MP -MF $(DEPDIR)/testFloatMatrixPow-testFloatPowm.Tpo -c -o testFloatMatrixPow-testFloatPowm.obj `if test -f 'testFloatPowm.c'; then $(CYGPATH_W) 'testFloatPowm.c'; else $(CYGPATH_W) '$(srcdir)/testFloatPowm.c'; fi`
+@am__fastdepCC_TRUE@	mv -f $(DEPDIR)/testFloatMatrixPow-testFloatPowm.Tpo $(DEPDIR)/testFloatMatrixPow-testFloatPowm.Po
+@AMDEP_TRUE@@am__fastdepCC_FALSE@	source='testFloatPowm.c' object='testFloatMatrixPow-testFloatPowm.obj' libtool=no @AMDEPBACKSLASH@
+@AMDEP_TRUE@@am__fastdepCC_FALSE@	DEPDIR=$(DEPDIR) $(CCDEPMODE) $(depcomp) @AMDEPBACKSLASH@
+@am__fastdepCC_FALSE@	$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(testFloatMatrixPow_CFLAGS) $(CFLAGS) -c -o testFloatMatrixPow-testFloatPowm.obj `if test -f 'testFloatPowm.c'; then $(CYGPATH_W) 'testFloatPowm.c'; else $(CYGPATH_W) '$(srcdir)/testFloatPowm.c'; fi`
+
 mostlyclean-libtool:
 	-rm -f *.lo
 
diff --git a/src/matrixOperations/powm/cpowma.c b/src/matrixOperations/powm/cpowma.c
index 51095302..5fb8e365 100644
--- a/src/matrixOperations/powm/cpowma.c
+++ b/src/matrixOperations/powm/cpowma.c
@@ -11,15 +11,53 @@
  */
 
 #include "matrixPow.h"
-#include "logm.h"
-#include "matrixExponential.h"
-#include "multiplication.h"
+#include "spec.h"
+#include "pow.h"
+#include "matrixTranspose.h"
+#include "conj.h"
+#include "matrixInversion.h"
+#include "matrixMultiplication.h"
 
-void cpowma(floatComplex* in, int rows, floatComplex expand, floatComplex* out){
-	int i=0;
-	/* use the formula a^b=exp(b*ln(a)) */
-	clogma(in,rows,out);
-	for(i=0;i<rows*rows;i++) out[i]= cmuls(expand,out[i]);
-	cexpma(out,out,(int)rows);
-					
+void cpowma(floatComplex* in, int rows, floatComplex power, floatComplex* out){
+	int i=0, j=0;
+	int hermitian=0;
+	floatComplex *eigenvalues,*eigenvectors,*tmp;
+	
+	/* Data initialization */
+	eigenvalues = malloc((uint)(rows*rows)*sizeof(floatComplex));
+	eigenvectors = malloc((uint)(rows*rows)*sizeof(floatComplex));
+	tmp = malloc((uint)(rows*rows)*sizeof(floatComplex));
+	
+	/* symmetric test*/
+	for(i=0;i<rows;i++) {
+		for (j=0;j<rows;j++)
+			if ((creals(in[i*rows+j])!=creals(in[j*rows+i])) || (cimags(in[i*rows+j])!=-cimags(in[j*rows+i]))) break;
+		
+		if (j!=rows) break;
+	}
+	
+	if ((i==rows)&&(j==rows)) hermitian=1;
+	
+	
+	/* find eigenvalues and eigenvectors */
+	cspec2a(in, rows, eigenvalues,eigenvectors);
+	
+	/* make operation on eigenvalues and eigenvectors */
+	for (i=0;i<rows;i++)
+			eigenvalues[i+i*rows]=cpows(eigenvalues[i+i*rows],power);
+			
+	cmulma(eigenvectors, rows, rows, eigenvalues, rows, rows, tmp);
+	
+	if (hermitian){		
+		ctransposea(eigenvectors, rows,rows, eigenvalues);
+		cconja(eigenvalues, rows*rows, eigenvalues);
+	}
+	else cinverma(eigenvectors, eigenvalues, rows);
+	
+	cmulma(tmp, rows, rows, eigenvalues, rows, rows, out);
+	
+	free(eigenvalues);
+	free(eigenvectors);
+	free(tmp);
+	
 }
diff --git a/src/matrixOperations/powm/dpowma.c b/src/matrixOperations/powm/dpowma.c
index 076aa405..5d67189b 100644
--- a/src/matrixOperations/powm/dpowma.c
+++ b/src/matrixOperations/powm/dpowma.c
@@ -11,65 +11,54 @@
  */
 
 #include "matrixPow.h"
-#include "lapack.h"
-#include "eye.h"
+#include "spec.h"
+#include "pow.h"
+#include "matrixTranspose.h"
+#include "conj.h"
+#include "matrixInversion.h"
 #include "matrixMultiplication.h"
+#include <stdio.h>
 
-
-void dpowma(double* in, int size, double expand, double* out){
-
-#ifndef WITHOUT_BLAS		
-	switch ((int)expand){
-		case 0 : 
-			deyea(out,size,size);
-			break;
-		case 1 : 
-			{
-				int i;
-				for (i=0;i<size*size;i++) out[i]=in[i];
-			}
-			break;
-		default : 
-			{
-				int i=0,j=0,k=0;
-				int One=1;
-				
-				for (i=1; i<expand; i++){
-					for (j=0;j<size;j++) /*column*/ {
-						for (k=0;k<size;k++)/*row*/{
-							out[k+j*size]=ddot_(&size,in+k,&size,in+j*size,&One);
-						}
-					}
-				}
-			}
-			break;
+void dpowma(double* in, int rows, double power, doubleComplex* out){
+	int i=0, j=0;
+	int symmetric=0;
+	doubleComplex *eigenvalues,*eigenvectors,*tmp;
+	
+	/* Data initialization */
+	eigenvalues = malloc((uint)(rows*rows)*sizeof(doubleComplex));
+	eigenvectors = malloc((uint)(rows*rows)*sizeof(doubleComplex));
+	tmp = malloc((uint)(rows*rows)*sizeof(doubleComplex));
+	
+	/* symmetric test*/
+	for(i=0;i<rows;i++) {
+		for (j=0;j<rows;j++)
+			if (in[i*rows+j]!=in[j*rows+i]) break;
+		
+		if (j!=rows) break;
 	}
-#else
-	switch ((int)expand){
-		case 0 : 
-			deyea(out,size,size);
-			break;
-		case 1 : 
-			{
-				int i;
-				for (i=0;i<size*size;i++) out[i]=in[i];
-			}
-			break;
-		default : 
-			{
-				int i=0,j=0;	
-				double* Pow;
-				Pow=malloc((uint)(size*size)*sizeof(double));
-				for (i=0;i<size*size;i++) out[i]=in[i];
-				for (i=1; i<expand; i++){
-					for (j=0;j<size*size;j++) Pow[j]=out[j];
-					dmulma(Pow,size,size,in,size,size,out);
-				}
-			}
-			break;
+	
+	if ((i==rows)&&(j==rows)) symmetric=1;
+	
+	
+	/* find eigenvalues and eigenvectors */
+	dspec2a(in, rows, eigenvalues,eigenvectors);
+/*	for (i=0;i<rows*rows;i++) printf("%f+%f*i\n",zreals(eigenvalues[i]),zimags(eigenvalues[i])); */
+	/* make operation on eigenvalues and eigenvectors */
+	for (i=0;i<rows;i++)
+			eigenvalues[i+i*rows]=zpows(eigenvalues[i+i*rows],DoubleComplex(power,0));
+			
+	zmulma(eigenvectors, rows, rows, eigenvalues, rows, rows, tmp);
+	
+	if (symmetric){		
+		ztransposea(eigenvectors, rows,rows, eigenvalues);
+		zconja(eigenvalues, rows*rows, eigenvalues);
 	}
-					
-
-#endif
-					
+	else zinverma(eigenvectors, eigenvalues, rows);
+	
+	zmulma(tmp, rows, rows, eigenvalues, rows, rows, out);
+	
+	free(eigenvalues);
+	free(eigenvectors);
+	free(tmp);
+	
 }
diff --git a/src/matrixOperations/powm/spowma.c b/src/matrixOperations/powm/spowma.c
index aa0f7043..634cd449 100644
--- a/src/matrixOperations/powm/spowma.c
+++ b/src/matrixOperations/powm/spowma.c
@@ -11,35 +11,53 @@
  */
 
 #include "matrixPow.h"
-#include "eye.h"
+#include "spec.h"
+#include "pow.h"
+#include "matrixTranspose.h"
+#include "conj.h"
+#include "matrixInversion.h"
 #include "matrixMultiplication.h"
 
-
-void spowma(float* in, int size, float expand, float* out){
-
-
-	switch ((int)expand){
-		case 0 : 
-			seyea(out,size,size);
-			break;
-		case 1 : 
-			{
-				int i;
-				for (i=0;i<size*size;i++) out[i]=in[i];
-			}
-			break;
-		default : 
-			{
-				int i=0,j=0;	
-				float* Pow;
-				Pow=malloc((uint)(size*size)*sizeof(float));
-				for (i=0;i<size*size;i++) out[i]=in[i];
-				for (i=1; i<expand; i++){
-					for (j=0;j<size*size;j++) Pow[j]=out[j];
-					smulma(Pow,size,size,in,size,size,out);
-				}
-			}
-			break;
+void spowma(float* in, int rows, float power, floatComplex* out){
+	int i=0, j=0;
+	int symmetric=0;
+	floatComplex *eigenvalues,*eigenvectors,*tmp;
+	
+	/* Data initialization */
+	eigenvalues = malloc((uint)(rows*rows)*sizeof(floatComplex));
+	eigenvectors = malloc((uint)(rows*rows)*sizeof(floatComplex));
+	tmp = malloc((uint)(rows*rows)*sizeof(floatComplex));
+	
+	/* symmetric test*/
+	for(i=0;i<rows;i++) {
+		for (j=0;j<rows;j++)
+			if (in[i*rows+j]!=in[j*rows+i]) break;
+		
+		if (j!=rows) break;
 	}
+	
+	if ((i==rows)&&(j==rows)) symmetric=1;
+	
+	
+	/* find eigenvalues and eigenvectors */
+	sspec2a(in, rows, eigenvalues,eigenvectors);
+	
+	/* make operation on eigenvalues and eigenvectors */
+	for (i=0;i<rows;i++)
+			eigenvalues[i+i*rows]=cpows(eigenvalues[i+i*rows],FloatComplex(power,0));
 			
+	cmulma(eigenvectors, rows, rows, eigenvalues, rows, rows, tmp);
+	
+	if (symmetric){		
+		ctransposea(eigenvectors, rows,rows, eigenvalues);
+		cconja(eigenvalues, rows*rows, eigenvalues);
+	}
+	else cinverma(eigenvectors, eigenvalues, rows);
+	
+	cmulma(tmp, rows, rows, eigenvalues, rows, rows, out);
+	
+	free(eigenvalues);
+	free(eigenvectors);
+	free(tmp);
+	
 }
diff --git a/src/matrixOperations/powm/testDoublePowm.c b/src/matrixOperations/powm/testDoublePowm.c
index 6fdbec10..05cad5e5 100644
--- a/src/matrixOperations/powm/testDoublePowm.c
+++ b/src/matrixOperations/powm/testDoublePowm.c
@@ -17,38 +17,126 @@
 
 
 static void dpowmaTest(void){
-	double in[4]={1,5,4,2};
-	double out[4];
+	double in1[4]={1,5,4,2};
+	double expand1=2.2;
+	double result1R[4]={  27.93459280052221771484 ,    23.580294119266994812278  ,
+	18.864235295413593007652 ,   32.650651624375619519469   };
+	double result1I[4]={  3.6611113731522362257920 , - 3.6611113731522362257920 , 
+  - 2.9288890985217883589087  ,  2.9288890985217883589087 };
+	doubleComplex out1[4];
 	int i;
 	
-	dpowma(in, 2, 2, out);
+	double in2[16]={    2.5358983855694532394409 ,   9.0725262500345706939697,    0.0026536155492067337036,    3.9639251008629798889160 , 
+	7.9845732506364583969116,    7.5407014600932598114014,   10.196942830458283424377 ,   8.2287722378969192504883 , 
+	10.538597775623202323914,    0.8204884417355060577393,    6.7301832754164934158325,    7.9482832476496696472168,  
+	8.7162081208080053329468 ,   2.3821726106107234954834 ,   6.5310877952724695205688,    2.784897476434707641602   };
+	double expand2 = 3.4683557949028909206390;
+	double result2R[16]={13801.893971410685480805 ,   9622.6108799100766191259 ,   10325.586569611912636901,   10694.791005280343597406  ,
+	24728.411825244897045195 ,   18392.823733925368287601  ,  18631.05868385956637212   ,  19357.84707477861229563   ,
+	16169.682243927050876664 ,   12258.542785024719705689  ,  12630.164466338968850323  ,  12827.915677254180991440  ,
+	13742.841851328515986097 ,   10198.0420642120679986    ,  10658.784670951883526868  ,  10839.51135004585739807  };
+	double result2I[16]={ - 7.1981835972120027378196  ,  1.9386514637886893552832,  - 17.692616672339234185074  ,  24.561537532538231687340 , 
+	- 2.2418859631076406557781 ,   0.6037961445855435371755,  - 5.5103941755046683681485,    7.649730724813480264857  , 
+	- 4.865855522250573272913  ,   1.310496989059492634056 ,  - 11.95992230200565309417 ,    16.603201547139228466676  ,
+	16.00935601900000193609 ,  - 4.3117212921047043394651  ,  39.34984366402868971591  , - 54.626892107189902958453 };
+	doubleComplex out2[16];
 	
-	for (i=0;i<4;i++) printf("out[%d] = %f\n",i,out[i]);
+	
+	dpowma(in1, 2, expand1, out1);
+	dpowma(in2, 4, expand2, out2);
+	
+	for (i=0;i<4;i++) {
+		assert( fabs(zreals(out1[i])-result1R[i]) / fabs(zreals(out1[i])) <1e-15);
+		assert( fabs(zimags(out1[i])-result1I[i]) / fabs(zimags(out1[i])) <1e-15);
+	}
+	
+/*	
+	FIXME : assert 1e-11 maybe due to spec2
+*/	
+	for (i=0;i<16;i++) {
+		printf("out[%d] = %1.16f+%1.16f*i --- result = %1.16f+%1.16f*i\n",i,zreals(out2[i]),zimags(out2[i]),result2R[i],result2I[i]);		
+		assert( fabs(zreals(out2[i])-result2R[i]) / fabs(zreals(out2[i])) <1e-14);
+		assert( fabs(zimags(out2[i])-result2I[i]) / fabs(zimags(out2[i])) <1e-11);
+	}	
 	
 }
 
-/* FIXME : assert à10^-12 */
+/* FIXME : assert 1e-14 */
 static void zpowmaTest(void){
-	double inR[9]={1,2,3,4,5,6,7,8,9};
-	double inI[9]={1,2,3,4,5,6,7,8,9};
-	double resultR[9]={- 4.7115011361608578610571,- 2.0782061409646632732517,0.5550888542315330909105,
-	- 2.3202132490900626571317,- 2.4412168031527574640904,- 2.5622203572154611528333,
-	0.0710746379807356554181,- 2.80422746534086453352,- 5.6795295686624518438634};
-	double resultI[9]={- 12.188702380084603049681,- 4.0827818504168584823333,4.0231386792508754268738,
-	- 3.0919079733956360556135,- 2.5964710348850239540752,- 2.1010340963744131848046,
-	6.0048864332933264975622,- 1.1101602193531934226201,- 8.2252068719997026846613};
-	doubleComplex *in,out[9];
-	int i;
+	/* Tests 1 */
+	{
+		double inR[9]={1,2,3,4,5,6,7,8,9};
+		double inI[9]={1,2,3,4,5,6,7,8,9};
+		double resultR[9]={- 4.7115011361608578610571,- 2.0782061409646632732517,0.5550888542315330909105,
+		- 2.3202132490900626571317,- 2.4412168031527574640904,- 2.5622203572154611528333,
+		0.0710746379807356554181,- 2.80422746534086453352,- 5.6795295686624518438634};
+		double resultI[9]={- 12.188702380084603049681,- 4.0827818504168584823333,4.0231386792508754268738,
+		- 3.0919079733956360556135,- 2.5964710348850239540752,- 2.1010340963744131848046,
+		6.0048864332933264975622,- 1.1101602193531934226201,- 8.2252068719997026846613};
+		doubleComplex *in,out[9];
+		int i;
+	
+		in=DoubleComplexMatrix(inR,inI,9);
+	
+		zpowma(in, 3, DoubleComplex(1,1), out);
+	
+		for (i=0;i<9;i++) printf("out[%d] = %f+%f*i\n",i,zreals(out[i]),zimags(out[i]));
+	
+		for (i=0;i<9;i++){
+			assert( (fabs(zreals(out[i])-resultR[i])/ fabs(zreals(out[i])) ) <1e-14);
+			assert( (fabs(zimags(out[i])-resultI[i])/ fabs(zimags(out[i])) ) <1e-14);
+		}
+	}
+	
+	/* Tests 2 and 3 */
+	{
+		double in1R[4]={1,5,4,2};
+		double in1I[4]={0};
+		double expand1=2.2;
+		double result1R[4]={  27.93459280052221771484 ,    23.580294119266994812278  ,
+		18.864235295413593007652 ,   32.650651624375619519469   };
+		double result1I[4]={  3.6611113731522362257920 , - 3.6611113731522362257920 , 
+	  - 2.9288890985217883589087  ,  2.9288890985217883589087 };
+		doubleComplex out1[4];
+		int i;
+	
+		double in2R[16]={    2.5358983855694532394409 ,   9.0725262500345706939697,    0.0026536155492067337036,    3.9639251008629798889160 , 
+		7.9845732506364583969116,    7.5407014600932598114014,   10.196942830458283424377 ,   8.2287722378969192504883 , 
+		10.538597775623202323914,    0.8204884417355060577393,    6.7301832754164934158325,    7.9482832476496696472168,  
+		8.7162081208080053329468 ,   2.3821726106107234954834 ,   6.5310877952724695205688,    2.784897476434707641602   };
+		double in2I[16]={0};
+		double expand2 = 3.4683557949028909206390;
+		double result2R[16]={13801.893971410685480805 ,   9622.6108799100766191259 ,   10325.586569611912636901,   10694.791005280343597406  ,
+		24728.411825244897045195 ,   18392.823733925368287601  ,  18631.05868385956637212   ,  19357.84707477861229563   ,
+		16169.682243927050876664 ,   12258.542785024719705689  ,  12630.164466338968850323  ,  12827.915677254180991440  ,
+		13742.841851328515986097 ,   10198.0420642120679986    ,  10658.784670951883526868  ,  10839.51135004585739807  };
+		double result2I[16]={ - 7.1981835972120027378196  ,  1.9386514637886893552832,  - 17.692616672339234185074  ,  24.561537532538231687340 , 
+		- 2.2418859631076406557781 ,   0.6037961445855435371755,  - 5.5103941755046683681485,    7.649730724813480264857  , 
+		- 4.865855522250573272913  ,   1.310496989059492634056 ,  - 11.95992230200565309417 ,    16.603201547139228466676  ,
+		16.00935601900000193609 ,  - 4.3117212921047043394651  ,  39.34984366402868971591  , - 54.626892107189902958453 };
+		doubleComplex out2[16];
+		doubleComplex *in1,*in2;
 	
-	in=DoubleComplexMatrix(inR,inI,9);
 	
-	zpowma(in, 3, DoubleComplex(1,1), out);
+		in1=DoubleComplexMatrix(in1R,in1I,4);
+		in2=DoubleComplexMatrix(in2R,in2I,16);
+		
+		zpowma(in1, 2, DoubleComplex(expand1,0), out1);
+		zpowma(in2, 4, DoubleComplex(expand2,0), out2);
 	
-	for (i=0;i<9;i++) printf("out[%d] = %f+%f*i\n",i,zreals(out[i]),zimags(out[i]));
+		for (i=0;i<4;i++) {
+			assert( fabs(zreals(out1[i])-result1R[i]) / fabs(zreals(out1[i])) <1e-15);
+			assert( fabs(zimags(out1[i])-result1I[i]) / fabs(zimags(out1[i])) <1e-15);
+		}
 	
-	for (i=0;i<9;i++){
-		assert( (fabs(zreals(out[i])-resultR[i])/ fabs(zreals(out[i])) ) <1e-12);
-		assert( (fabs(zimags(out[i])-resultI[i])/ fabs(zimags(out[i])) ) <1e-13);
+	/*	
+		FIXME : assert 1e-11 maybe due to spec2
+	*/	
+		for (i=0;i<16;i++) {
+			printf("out[%d] = %1.16f+%1.16f*i --- result = %1.16f+%1.16f*i\n",i,zreals(out2[i]),zimags(out2[i]),result2R[i],result2I[i]);		
+			assert( fabs(zreals(out2[i])-result2R[i]) / fabs(zreals(out2[i])) <1e-14);
+			assert( fabs(zimags(out2[i])-result2I[i]) / fabs(zimags(out2[i])) <1e-11);
+		}	
 	}
 	
 }
diff --git a/src/matrixOperations/powm/zpowma.c b/src/matrixOperations/powm/zpowma.c
index c8cac267..2def5f40 100644
--- a/src/matrixOperations/powm/zpowma.c
+++ b/src/matrixOperations/powm/zpowma.c
@@ -11,15 +11,53 @@
  */
 
 #include "matrixPow.h"
-#include "logm.h"
-#include "matrixExponential.h"
-#include "multiplication.h"
+#include "spec.h"
+#include "pow.h"
+#include "matrixTranspose.h"
+#include "conj.h"
+#include "matrixInversion.h"
+#include "matrixMultiplication.h"
 
-void zpowma(doubleComplex* in, int rows, doubleComplex expand, doubleComplex* out){
-	int i=0;
-	/* use the formula a^b=exp(b*ln(a)) */
-	zlogma(in,rows,out);
-	for(i=0;i<rows*rows;i++) out[i]= zmuls(expand,out[i]);
-	zexpma(out,out,(int)rows);
-					
+void zpowma(doubleComplex* in, int rows, doubleComplex power, doubleComplex* out){
+	int i=0, j=0;
+	int hermitian=0;
+	doubleComplex *eigenvalues,*eigenvectors,*tmp;
+	
+	/* Data initialization */
+	eigenvalues = malloc((uint)(rows*rows)*sizeof(doubleComplex));
+	eigenvectors = malloc((uint)(rows*rows)*sizeof(doubleComplex));
+	tmp = malloc((uint)(rows*rows)*sizeof(doubleComplex));
+	
+	/* hermitian test*/
+	for(i=0;i<rows;i++) {
+		for (j=0;j<rows;j++)
+			if ((zreals(in[i*rows+j])!=zreals(in[j*rows+i])) || (zimags(in[i*rows+j])!=-zimags(in[j*rows+i]))) break;
+		
+		if (j!=rows) break;
+	}
+	
+	if ((i==rows)&&(j==rows)) hermitian=1;
+	
+	
+	/* find eigenvalues and eigenvectors */
+	zspec2a(in, rows, eigenvalues,eigenvectors);
+	
+	/* make operation on eigenvalues and eigenvectors */
+	for (i=0;i<rows;i++)
+			eigenvalues[i+i*rows]=zpows(eigenvalues[i+i*rows],power);
+			
+	zmulma(eigenvectors, rows, rows, eigenvalues, rows, rows, tmp);
+	
+	if (hermitian){		
+		ztransposea(eigenvectors, rows,rows, eigenvalues);
+		zconja(eigenvalues, rows*rows, eigenvalues);
+	}
+	else zinverma(eigenvectors, eigenvalues, rows);
+	
+	zmulma(tmp, rows, rows, eigenvalues, rows, rows, out);
+	
+	free(eigenvalues);
+	free(eigenvectors);
+	free(tmp);
+	
 }