diff options
Diffstat (limited to '2.3-1/src/c/matrixOperations/powm/cpowma.c')
-rw-r--r-- | 2.3-1/src/c/matrixOperations/powm/cpowma.c | 64 |
1 files changed, 64 insertions, 0 deletions
diff --git a/2.3-1/src/c/matrixOperations/powm/cpowma.c b/2.3-1/src/c/matrixOperations/powm/cpowma.c new file mode 100644 index 00000000..d023a7aa --- /dev/null +++ b/2.3-1/src/c/matrixOperations/powm/cpowma.c @@ -0,0 +1,64 @@ +/* + * Scilab ( http://www.scilab.org/ ) - This file is part of Scilab + * Copyright (C) 2008 - INRIA - Arnaud TORSET + * + * 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 "matrixPow.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 power, floatComplex* out){ + int i=0, j=0; + int hermitian=0; + floatComplex *eigenvalues,*eigenvectors,*tmp; + + /* Data initialization */ + eigenvalues = (floatComplex*)malloc((unsigned int)(rows*rows)*sizeof(floatComplex)); + eigenvectors = (floatComplex*)malloc((unsigned int)(rows*rows)*sizeof(floatComplex)); + tmp = (floatComplex*)malloc((unsigned int)(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); + +} |