Mentions légales du service

Skip to content
Snippets Groups Projects
Commit ce364670 authored by hhakim's avatar hhakim
Browse files

Implement Faust over complex field into Python wrapper and add related unit tests.

parent 21b920f7
No related branches found
No related tags found
No related merge requests found
......@@ -168,7 +168,7 @@ class TestFaustPy(unittest.TestCase):
print("testMul()")
rmat = np.random.rand(self.F.get_nb_cols(),
self.r.randint(1,TestFaustPy.MAX_DIM_SIZE))
prod = self.mulFactors()*rmat
prod = self.mulFactors().dot(rmat)
test_prod = self.F*rmat
self.assertProdEq(prod, test_prod)
......@@ -233,10 +233,33 @@ class TestFaustPy(unittest.TestCase):
print("Test Faust.getH()")
test_Fct = self.F.getH().todense()
ref_Fct = self.F.todense().conj().T
print("test_Fct=", test_Fct)
print("ref_Fct=", ref_Fct)
ref_Fct[ref_Fct==0] = 1
test_Fct[test_Fct==0] = 1
self.assertTrue(((((test_Fct-ref_Fct)/ref_Fct) < 0.01)).all())
class TestFaustPyCplx(TestFaustPy):
def setUp(self):
""" Initializes the tests objects """
r = random.Random() # initialized from time or system
num_factors = r.randint(1, TestFaustPy.MAX_NUM_FACTORS)
factors = []
d2 = r.randint(1, TestFaustPy.MAX_DIM_SIZE)
for i in range(0, num_factors):
d1, d2 = d2, r.randint(1, TestFaustPy.MAX_DIM_SIZE)
factors += [np.random.rand(d1, d2).astype(np.complex)]
factors[i].imag = [np.random.rand(d1, d2)]
self.F = Faust(factors)
self.factors = factors
print("Tests on random Faust with dims=", self.F.get_nb_rows(),
self.F.get_nb_cols())
print("Num. factors:", num_factors)
self.r = r
self.num_factors = num_factors
if __name__ == "__main__":
if(len(sys.argv)> 1):
# argv[1] is for adding a directory in PYTHONPATH
......
......@@ -332,3 +332,8 @@ class Faust:
if(format == "Matlab"):
F.m_faust.save_mat_file(filepath)
def isReal(F):
"""
Returns True if F is a real Faust and False if it's a complex Faust.
"""
return F.m_faust.isReal()
......@@ -47,39 +47,57 @@ import copy
from libc.stdlib cimport malloc, free;
from libc.string cimport memcpy;
from libcpp cimport bool
from libcpp cimport complex
from cpython.mem cimport PyMem_Malloc, PyMem_Realloc, PyMem_Free
cdef class FaustCore:
#### ATTRIBUTE ########
# classe Cython
cdef FaustCoreCy.FaustCoreCpp[double]* m_faust
cdef FaustCoreCy.FaustCoreCpp[double]* core_faust_dbl
cdef FaustCoreCy.FaustCoreCpp[complex]* core_faust_cplx
cdef bool _isReal
#### CONSTRUCTOR ####
#def __cinit__(self,np.ndarray[double, mode="fortran", ndim=2] mat):
def __cinit__(self,list_factors=None, core=False):
#print 'inside cinit'
#self.m_faust = FaustCoreCy.FaustCoreCpp[double]();
cdef double [:,:] data
cdef complex [:,:] data_cplx
cdef unsigned int nbrow
cdef unsigned int nbcol
if(list_factors is not None):
#print 'longueur list'
#print len(list_factors)
#print 'avant boucle'
self.m_faust = new FaustCoreCy.FaustCoreCpp[double]()
self._isReal = True
for factor in list_factors:
if(not isinstance(factor, np.ndarray)):
raise ValueError('Faust factors must be numpy.ndarray')
if(isinstance(factor[0,0], np.complex)):
# str(factor.dtype) == complex128/64/complex_
self._isReal = False
break
if(self._isReal):
self.core_faust_dbl = new FaustCoreCy.FaustCoreCpp[double]()
else:
self.core_faust_cplx = new FaustCoreCy.FaustCoreCpp[complex]()
for factor in list_factors:
data=factor.astype(float,'F');
if(self._isReal):
data=factor.astype(float,'F')
else:
data_cplx=factor.astype(np.complex128,'F')
nbrow=factor.shape[0];
nbcol=factor.shape[1];
# print nbrow
# print nbcol
self.m_faust.push_back(&data[0,0],nbrow,nbcol);
if(self._isReal):
self.core_faust_dbl.push_back(&data[0,0], nbrow, nbcol)
else:
self.core_faust_cplx.push_back(&data_cplx[0,0], nbrow, nbcol)
#print(self.__dict__)
#print 'apres boucle'
elif(core): # trick to initialize a new FaustCoreCpp from C++ (see
# transpose)
# transpose, conj and adjoint)
pass
#else:
#TODO: raise error for undefined object here
......@@ -87,7 +105,7 @@ cdef class FaustCore:
#### METHOD ####
#~ def getNbRow(self):
#~ #return self.m_faust.getNbRow();
#~ #return self.core_faust_dbl.getNbRow();
#~ (dim1,dim2)=self.shape();
#~ return dim1
......@@ -99,8 +117,12 @@ cdef class FaustCore:
def shape(self):
cdef unsigned int nbrow = 0
cdef unsigned int nbcol = 0
nbrow = self.m_faust.getNbRow();
nbcol = self.m_faust.getNbCol();
if(self._isReal):
nbrow = self.core_faust_dbl.getNbRow();
nbcol = self.core_faust_dbl.getNbCol();
else:
nbrow = self.core_faust_cplx.getNbRow();
nbcol = self.core_faust_cplx.getNbCol();
return (nbrow,nbcol)
......@@ -122,9 +144,14 @@ cdef class FaustCore:
if not isinstance(M, (np.ndarray) ):
raise ValueError('input M must a numpy ndarray')
#transform into float F continous matrix
M=M.astype(float,'F')
if not M.dtype=='float':
raise ValueError('input M must be double array')
if(self._isReal):
M=M.astype(float,'F')
if not M.dtype=='float':
raise ValueError('input M must be double array')
else:
M=M.astype(complex,'F')
if(M.dtype not in ['complex', 'complex128', 'complex64'] ): #could fail if complex128 etc.
raise ValueError('input M must be complex array')
if not M.flags['F_CONTIGUOUS']:
raise ValueError('input M must be Fortran contiguous (Colmajor)')
......@@ -146,13 +173,21 @@ cdef class FaustCore:
cdef double[:] xview_1D
cdef double[:,:] xview_2D
cdef complex[:] xview_1D_cplx
cdef complex[:,:] xview_2D_cplx
if ndim_M == 1:
nbcol_x=1
xview_1D=M;
if(self._isReal):
xview_1D=M
else:
xview_1D_cplx=M
else:
nbcol_x=M.shape[1]
xview_2D=M;
if(self._isReal):
xview_2D=M
else:
xview_2D_cplx=M
if (nbrow_x != nbColThis):
raise ValueError('y=F*M multiplication with Faust: invalid dimension of the input matrix M');
......@@ -160,12 +195,27 @@ cdef class FaustCore:
#void multiply(FPP* value_y,int nbrow_y,int nbcol_y,FPP* value_x,int nbrow_x,int nbcol_x,bool isTranspose);
nbcol_y = nbcol_x;
cdef y = np.zeros([nbrow_y,nbcol_y], dtype='d',order='F')
cdef double[:,:] yview=y
cdef y
cdef double[:,:] yview
cdef complex[:,:] yview_cplx
if(self._isReal):
y = np.zeros([nbrow_y,nbcol_y], dtype='d',order='F')
yview = y
else:
y = np.zeros([nbrow_y, nbcol_y], dtype='complex', order='F')
yview_cplx = y
if ndim_M == 1:
self.m_faust.multiply(&yview[0,0],nbrow_y,nbcol_y,&xview_1D[0],nbrow_x,nbcol_x)
if(self._isReal):
self.core_faust_dbl.multiply(&yview[0,0],nbrow_y,nbcol_y,&xview_1D[0],nbrow_x,nbcol_x)
else:
self.core_faust_cplx.multiply(&yview_cplx[0,0], nbrow_y,
nbcol_y, &xview_1D_cplx[0], nbrow_x,nbcol_x)
else:
self.m_faust.multiply(&yview[0,0],nbrow_y,nbcol_y,&xview_2D[0,0],nbrow_x,nbcol_x)
if(self._isReal):
self.core_faust_dbl.multiply(&yview[0,0],nbrow_y,nbcol_y,&xview_2D[0,0],nbrow_x,nbcol_x)
else:
self.core_faust_cplx.multiply(&yview_cplx[0,0],nbrow_y,nbcol_y,&xview_2D_cplx[0,0],nbrow_x,nbcol_x)
return y
......@@ -173,32 +223,51 @@ cdef class FaustCore:
# print information about the faust (size, number of factor, type of factor (dense/sparse) ...)
def display(self):
self.m_faust.Display();
self.core_faust_dbl.Display();
def nnz(self):
cdef unsigned long long nnz = 0
nnz = self.m_faust.nnz()
if(self._isReal):
nnz = self.core_faust_dbl.nnz()
else:
nnz = self.core_faust_cplx.nnz()
return nnz
def norm(self):
cdef double norm
norm = self.m_faust.norm()
if(self._isReal):
norm = self.core_faust_dbl.norm()
else:
norm = self.core_faust_cplx.norm()
return norm
def get_nb_factors(self):
cdef int nb_factors
nb_factors = int(self.m_faust.get_nb_factors())
if(self._isReal):
nb_factors = int(self.core_faust_dbl.get_nb_factors())
else:
nb_factors = int(self.core_faust_cplx.get_nb_factors())
return nb_factors
def get_fact(self,i):
if(i >= self.get_nb_factors() or i < 0):
raise ValueError("factor index must be greater or equal 0 and "
"lower than "+str(self.get_nb_factors())+".")
cdef fact = np.zeros([self.m_faust.get_fact_nb_rows(i),
self.m_faust.get_fact_nb_cols(i)], dtype='d',
cdef fact
cdef double[:,:] fact_dbl_view
cdef complex[:,:] fact_cplx_view
if(self._isReal):
fact = np.zeros([self.core_faust_dbl.get_fact_nb_rows(i),
self.core_faust_dbl.get_fact_nb_cols(i)], dtype='d',
order='F')
cdef double[:,:] fact_view = fact
self.m_faust.get_fact(i, &fact_view[0, 0])
fact_dbl_view = fact
self.core_faust_dbl.get_fact(i, &fact_dbl_view[0, 0])
else:
fact = np.zeros([self.core_faust_cplx.get_fact_nb_rows(i),
self.core_faust_cplx.get_fact_nb_cols(i)],
dtype='complex', order='F')
fact_cplx_view = fact
self.core_faust_cplx.get_fact(i, &fact_cplx_view[0, 0])
return fact
def save_mat_file(self,filepath):
......@@ -208,25 +277,45 @@ cdef class FaustCore:
for i in range(0,len(filepath)):
cfilepath[i] = fparr[i]
cfilepath[i+1] = 0
self.m_faust.save_mat_file(cfilepath)
if(self._isReal):
self.core_faust_dbl.save_mat_file(cfilepath)
else:
self.core_faust_cplx.save_mat_file(cfilepath)
PyMem_Free(cfilepath)
def transpose(self):
core = FaustCore(core=True)
core.m_faust = self.m_faust.transpose()
if(self._isReal):
core.core_faust_dbl = self.core_faust_dbl.transpose()
else:
core.core_faust_cplx = self.core_faust_cplx.transpose()
core._isReal = self._isReal
return core
def conj(self):
core = FaustCore(core=True)
core.m_faust = self.m_faust.conjugate()
if(self._isReal):
core.core_faust_dbl = self.core_faust_dbl.conjugate()
else:
core.core_faust_cplx = self.core_faust_cplx.conjugate()
core._isReal = self._isReal
return core
def getH(self):
core = FaustCore(core=True)
core.m_faust = self.m_faust.adjoint()
if(self._isReal):
core.core_faust_dbl = self.core_faust_dbl.adjoint()
else:
core.core_faust_cplx = self.core_faust_cplx.adjoint()
core._isReal = self._isReal
return core
def isReal(self):
if(self._isReal): return True
return False
def __dealloc__(self):
del self.m_faust
if(self._isReal):
del self.core_faust_dbl
else:
del self.core_faust_cplx
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment