Mentions légales du service

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

Add err_period argument to pyfaust eigtj.

parent 95e3c3ba
Branches
Tags
No related merge requests found
......@@ -85,7 +85,7 @@ void GivensFGFTGen<FPP,DEVICE,FPP2,FPP4>::compute_facts()
while(! stopping && (J == 0 || ite < facts.size())) // when J == 0 the stopping criterion is the error against Lap
{
next_step();
if(stopping = (ite > 1 && stoppingCritIsError && errs.size() > 2 && errs[ite-1]-errs[ite-2] > FLT_EPSILON))
if(stopping = (ite > 1 && stoppingCritIsError && errs.size() > 2 && errs[errs.size()-1]-errs[errs.size()-2] > FLT_EPSILON))
/*if(verbosity>0) */cerr << "WARNING: the eigtj algorithm stopped because the last error is greater than the previous one." << endl;
if(stopping || stoppingCritIsError && errs.size() > 0 && (*(errs.end()-1) - stoppingError ) < FLT_EPSILON)
{
......@@ -348,7 +348,7 @@ FPP2 GivensFGFTGen<FPP,DEVICE,FPP2,FPP4>::calc_err()
template<typename FPP, FDevice DEVICE, typename FPP2, typename FPP4>
void GivensFGFTGen<FPP, DEVICE, FPP2, FPP4>::update_err()
{
if(!(this->ite % this->err_period) && this->stoppingCritIsError || this->verbosity > 1)
if(!((this->ite+1) % this->err_period) && this->stoppingCritIsError || this->verbosity > 1)
{
auto err = this->calc_err();
if(this->verbosity)
......
......@@ -329,7 +329,7 @@ def svdtj(M, nGivens=None, tol=0, err_period=100, relerr=True,
V = Faust(core_obj=Vcore)
return U, S, V
def eigtj(M, nGivens=None, tol=0, order='ascend', relerr=True,
def eigtj(M, nGivens=None, tol=0, err_period=100, order='ascend', relerr=True,
nGivens_per_fac=None, verbosity=0, enable_large_Faust=False):
"""
Performs an approximate eigendecomposition of M and returns the eigenvalues in W along with the corresponding normalized right eigenvectors (as the columns of the Faust object V).
......@@ -351,6 +351,10 @@ def eigtj(M, nGivens=None, tol=0, order='ascend', relerr=True,
tol: (float) the tolerance error at which the algorithm stops. The
default value is zero so that stopping is based on reaching the
targeted nGivens (this argument is optional only if nGivens is set).
err_period: (int) it defines the period, in number of factors of V
the error is compared to tol (reducing the period spares
some factors but increases slightly the computational cost because the error
is computed more often).
order: order of eigenvalues, possible choices are ‘ascend,
'descend' or 'undef' (to avoid a sorting operation and save some time).
nGivens_per_fac: (int) targeted number of Givens rotations per factor
......@@ -431,19 +435,17 @@ def eigtj(M, nGivens=None, tol=0, order='ascend', relerr=True,
M_dtype = _sanitize_dtype(M.dtype)
eigtj_args = [M, nGivens, tol, relerr,
nGivens_per_fac, verbosity, order,
enable_large_Faust,
err_period]
if M_dtype == 'float32':
D, core_obj = _FaustCorePy.FaustAlgoGenGivensFlt.eigtj(M, nGivens, tol, relerr,
nGivens_per_fac, verbosity, order,
enable_large_Faust)
D, core_obj = _FaustCorePy.FaustAlgoGenGivensFlt.eigtj(*eigtj_args)
elif M_dtype == 'float64':
D, core_obj = _FaustCorePy.FaustAlgoGenGivensDbl.eigtj(M, nGivens, tol, relerr,
nGivens_per_fac, verbosity, order,
enable_large_Faust)
D, core_obj = _FaustCorePy.FaustAlgoGenGivensDbl.eigtj(*eigtj_args)
else: # M_dtype == 'complex'
D, core_obj = _FaustCorePy.FaustAlgoGenGivensCplxDbl.eigtj(M, nGivens, tol, relerr,
nGivens_per_fac, verbosity, order,
enable_large_Faust)
D, core_obj = _FaustCorePy.FaustAlgoGenGivensCplxDbl.eigtj(*eigtj_args)
return D, Faust(core_obj=core_obj)
def _check_fact_mat(funcname, M, is_real=None):
......
......@@ -164,7 +164,8 @@ cdef extern from "FaustFactGivensFGFT.h":
const FPP2 stoppingError,
const bool errIsRel,
const int order,
const bool enable_large_Faust)
const bool enable_large_Faust,
const int err_period)
cdef FaustCoreCppCPU[FPP]* fact_givens_fgft_sparse[FPP,FPP2](const FPP* data, int* row_ptr,
int* id_col, int nnz, unsigned int num_rows,
......@@ -174,17 +175,18 @@ cdef extern from "FaustFactGivensFGFT.h":
const FPP2 stoppingError,
const bool errIsRel,
const int order,
const bool enable_large_Faust)
const bool enable_large_Faust,
const int err_period)
cdef FaustCoreCppCPU[FPP]* fact_givens_fgft_cplx[FPP,FPP2](const
FPP* Lap, unsigned int num_rows,
cdef FaustCoreCppCPU[FPP]* fact_givens_fgft_cplx[FPP,FPP2](const FPP* Lap, unsigned int num_rows,
unsigned int num_cols, unsigned int J,
unsigned int t,
FPP2* D, unsigned int verbosity,
const FPP2 stoppingError,
const bool errIsRel,
const int order,
const bool enable_large_Faust)
const bool enable_large_Faust,
const int err_period)
cdef FaustCoreCppCPU[FPP]* fact_givens_fgft_sparse_cplx[FPP, FPP2](const FPP* data, int* row_ptr,
int* id_col, int nnz, unsigned int num_rows,
......@@ -195,7 +197,8 @@ cdef extern from "FaustFactGivensFGFT.h":
const FPP2 stoppingError,
const bool errIsRel,
const int order,
const bool enable_large_Faust)
const bool enable_large_Faust,
const int err_period)
cdef void svdtj[FPP, FPP2](FaustCoreCppCPU[FPP]** U, FaustCoreCppCPU[FPP] **V, FPP* S,
const FPP* M_data,
......
......@@ -3,26 +3,26 @@
#include "faust_GivensFGFT.h"
template<typename FPP, typename FPP2 = float>
FaustCoreCpp<FPP>* fact_givens_fgft(const FPP* Lap, unsigned int num_rows, unsigned int num_cols, unsigned int J, unsigned int t, FPP* D, unsigned int verbosity = 0, const FPP2 stoppingError = 0.0, const bool errIsRel = true, const int order = 1 /* ascendant order for eigenvalues */);
FaustCoreCpp<FPP>* fact_givens_fgft(const FPP* Lap, unsigned int num_rows, unsigned int num_cols, unsigned int J, unsigned int t, FPP* D, unsigned int verbosity = 0, const FPP2 stoppingError = 0.0, const bool errIsRel = true, const int order = 1 /* ascendant order for eigenvalues */, const bool enable_large_Faust=false, const int err_period=100);
template<typename FPP, typename FPP2 = float>
FaustCoreCpp<FPP>* fact_givens_fgft_sparse(const FPP* data, int* row_ptr, int* id_col, int nnz, int nrows, int ncols,
unsigned int J, unsigned int t /* end of input parameters*/, FPP* D, unsigned int verbosity = 0, const FPP2 stoppingError = 0.0, const bool errIsRel = true, const int order = 1 /* ascendant order for eigenvalues */, const bool enable_large_Faust=false);
unsigned int J, unsigned int t /* end of input parameters*/, FPP* D, unsigned int verbosity = 0, const FPP2 stoppingError = 0.0, const bool errIsRel = true, const int order = 1 /* ascendant order for eigenvalues */, const bool enable_large_Faust=false, const int err_period=100);
template<typename FPP, typename FPP2 = float>
FaustCoreCpp<FPP>* fact_givens_fgft_generic(GivensFGFT<FPP, Cpu, FPP2>* algo, FPP* D, const int order, const bool enable_large_Faust=false);
FaustCoreCpp<FPP>* fact_givens_fgft_generic(GivensFGFT<FPP, Cpu, FPP2>* algo, FPP* D, const int order);
template<typename FPP, typename FPP2 = float>
FaustCoreCpp<FPP>* fact_givens_fgft_cplx(const FPP* Lap, unsigned int num_rows, unsigned int num_cols, unsigned int J, unsigned int t, FPP2* D, unsigned int verbosity = 0, const FPP2 stoppingError = 0.0, const bool errIsRel = true, const int order = 1 /* ascendant order for eigenvalues */, const bool enable_large_Faust=false);
FaustCoreCpp<FPP>* fact_givens_fgft_cplx(const FPP* Lap, unsigned int num_rows, unsigned int num_cols, unsigned int J, unsigned int t, FPP2* D, unsigned int verbosity = 0, const FPP2 stoppingError = 0.0, const bool errIsRel = true, const int order = 1 /* ascendant order for eigenvalues */, const bool enable_large_Faust=false, const int err_period=100);
template<typename FPP, typename FPP2 = float>
FaustCoreCpp<FPP>* fact_givens_fgft_sparse_cplx(const FPP* data, int* row_ptr, int* id_col, int nnz, int nrows, int ncols,
unsigned int J, unsigned int t, FPP2* D, unsigned int verbosity = 0, const FPP2 stoppingError = 0.0, const bool errIsRel = true, const int order = 1 /* ascendant order for eigenvalues */, const bool enable_large_Faust=false);
unsigned int J, unsigned int t, FPP2* D, unsigned int verbosity = 0, const FPP2 stoppingError = 0.0, const bool errIsRel = true, const int order = 1 /* ascendant order for eigenvalues */, const bool enable_large_Faust=false, const int err_period=100);
template<typename FPP, typename FPP2 = float>
FaustCoreCpp<FPP>* fact_givens_fgft_generic_cplx(GivensFGFTComplex<FPP, Cpu, FPP2>* algo, FPP2* D, const int order, const bool enable_large_Faust=false);
FaustCoreCpp<FPP>* fact_givens_fgft_generic_cplx(GivensFGFTComplex<FPP, Cpu, FPP2>* algo, FPP2* D, const int order);
template<typename FPP, typename FPP2 = float>
void svdtj( FaustCoreCpp<FPP>** U, FaustCoreCpp<FPP> **V, FPP* S /** end of output parameters */, const FPP* M, unsigned int num_rows, unsigned int num_cols, unsigned int J1, unsigned int J2, unsigned int t1, unsigned int t2, unsigned int verbosity = 0, const FPP2 stoppingError = 0.0, const bool errIsRel = true, const bool enable_large_Faust = false /* end of input parameters*/, const int err_period=100);
......
......@@ -8,23 +8,23 @@ using namespace Faust;
template<typename FPP, typename FPP2>
FaustCoreCpp<FPP>* fact_givens_fgft_sparse(FPP* data, int* row_ptr, int* id_col, int nnz, int nrows, int ncols,
unsigned int J, unsigned int t /* end of input parameters*/, FPP* D, unsigned int verbosity, const FPP2 stoppingError, const bool errIsRel, const int order, const bool enable_large_Faust)
unsigned int J, unsigned int t /* end of input parameters*/, FPP* D, unsigned int verbosity, const FPP2 stoppingError, const bool errIsRel, const int order, const bool enable_large_Faust, const int err_period)
{
Faust::MatSparse<FPP, Cpu> mat_Lap(nnz, nrows, ncols, data, id_col, row_ptr);
GivensFGFT<FPP, Cpu, FPP2>* algo;
if(t <= 1)
{
algo = new GivensFGFT<FPP, Cpu, FPP2>(mat_Lap, (int)J, verbosity, stoppingError, errIsRel, enable_large_Faust);
algo = new GivensFGFT<FPP, Cpu, FPP2>(mat_Lap, (int)J, verbosity, stoppingError, errIsRel, enable_large_Faust, err_period);
}
else
{
algo = new GivensFGFTParallel<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity, stoppingError, errIsRel, enable_large_Faust);
algo = new GivensFGFTParallel<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity, stoppingError, errIsRel, enable_large_Faust, err_period);
}
return fact_givens_fgft_generic(algo, D, order);
}
template<typename FPP, typename FPP2>
FaustCoreCpp<FPP>* fact_givens_fgft(const FPP* Lap, unsigned int num_rows, unsigned int num_cols, unsigned int J, unsigned int t /* end of input parameters*/, FPP* D, unsigned int verbosity, const FPP2 stoppingError, const bool errIsRel, const int order, const bool enable_large_Faust)
FaustCoreCpp<FPP>* fact_givens_fgft(const FPP* Lap, unsigned int num_rows, unsigned int num_cols, unsigned int J, unsigned int t /* end of input parameters*/, FPP* D, unsigned int verbosity, const FPP2 stoppingError, const bool errIsRel, const int order, const bool enable_large_Faust, const int err_period)
{
//TODO: optimization possible here by avoiding Lap copy in MatDense (by
//just using the data in Lap as underlying pointer of MatDense)
......@@ -33,17 +33,17 @@ FaustCoreCpp<FPP>* fact_givens_fgft(const FPP* Lap, unsigned int num_rows, unsig
GivensFGFT<FPP, Cpu, FPP2>* algo;
if(t <= 1)
{
algo = new GivensFGFT<FPP, Cpu, FPP2>(mat_Lap, (int)J, verbosity, stoppingError, errIsRel, enable_large_Faust);
algo = new GivensFGFT<FPP, Cpu, FPP2>(mat_Lap, (int)J, verbosity, stoppingError, errIsRel, enable_large_Faust, err_period);
}
else
{
algo = new GivensFGFTParallel<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity, stoppingError, errIsRel, enable_large_Faust);
algo = new GivensFGFTParallel<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity, stoppingError, errIsRel, enable_large_Faust, err_period);
}
return fact_givens_fgft_generic(algo, D, order);
}
template<typename FPP, typename FPP2>
FaustCoreCpp<FPP>* fact_givens_fgft_generic(GivensFGFT<FPP, Cpu, FPP2>* algo, FPP* D, const int order, const bool enable_large_Faust)
FaustCoreCpp<FPP>* fact_givens_fgft_generic(GivensFGFT<FPP, Cpu, FPP2>* algo, FPP* D, const int order)
{
FaustCoreCpp<FPP>* fc = nullptr;
......@@ -70,23 +70,23 @@ FaustCoreCpp<FPP>* fact_givens_fgft_generic(GivensFGFT<FPP, Cpu, FPP2>* algo, FP
template<typename FPP, typename FPP2>
FaustCoreCpp<FPP>* fact_givens_fgft_sparse_cplx(FPP* data, int* row_ptr, int* id_col, int nnz, int nrows, int ncols,
unsigned int J, unsigned int t /* end of input parameters*/, FPP2* D, unsigned int verbosity, const FPP2 stoppingError, const bool errIsRel, const int order, const bool enable_large_Faust)
unsigned int J, unsigned int t /* end of input parameters*/, FPP2* D, unsigned int verbosity, const FPP2 stoppingError, const bool errIsRel, const int order, const bool enable_large_Faust, const int err_period)
{
Faust::MatSparse<FPP, Cpu> mat_Lap(nnz, nrows, ncols, data, id_col, row_ptr);
GivensFGFTComplex<FPP, Cpu, FPP2>* algo;
if(t <= 1)
{
algo = new GivensFGFTComplex<FPP, Cpu, FPP2>(mat_Lap, (int)J, verbosity, stoppingError, errIsRel, enable_large_Faust);
algo = new GivensFGFTComplex<FPP, Cpu, FPP2>(mat_Lap, (int)J, verbosity, stoppingError, errIsRel, enable_large_Faust, err_period);
}
else
{
algo = new GivensFGFTParallelComplex<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity, stoppingError, errIsRel, enable_large_Faust);
algo = new GivensFGFTParallelComplex<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity, stoppingError, errIsRel, enable_large_Faust, err_period);
}
return fact_givens_fgft_generic_cplx(algo, D, order);
}
template<typename FPP, typename FPP2>
FaustCoreCpp<FPP>* fact_givens_fgft_cplx(const FPP* Lap, unsigned int num_rows, unsigned int num_cols, unsigned int J, unsigned int t /* end of input parameters*/, FPP2* D, unsigned int verbosity, const FPP2 stoppingError, const bool errIsRel, const int order, const bool enable_large_Faust)
FaustCoreCpp<FPP>* fact_givens_fgft_cplx(const FPP* Lap, unsigned int num_rows, unsigned int num_cols, unsigned int J, unsigned int t /* end of input parameters*/, FPP2* D, unsigned int verbosity, const FPP2 stoppingError, const bool errIsRel, const int order, const bool enable_large_Faust, const int err_period)
{
//TODO: optimization possible here by avoiding Lap copy in MatDense (by
//just using the data in Lap as underlying pointer of MatDense)
......@@ -95,17 +95,17 @@ FaustCoreCpp<FPP>* fact_givens_fgft_cplx(const FPP* Lap, unsigned int num_rows,
GivensFGFTComplex<FPP, Cpu, FPP2>* algo;
if(t <= 1)
{
algo = new GivensFGFTComplex<FPP, Cpu, FPP2>(mat_Lap, (int)J, verbosity, stoppingError, errIsRel, enable_large_Faust);
algo = new GivensFGFTComplex<FPP, Cpu, FPP2>(mat_Lap, (int)J, verbosity, stoppingError, errIsRel, enable_large_Faust, err_period);
}
else
{
algo = new GivensFGFTParallelComplex<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity, stoppingError, errIsRel, enable_large_Faust);
algo = new GivensFGFTParallelComplex<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity, stoppingError, errIsRel, enable_large_Faust, err_period);
}
return fact_givens_fgft_generic_cplx(algo, D, order);
}
template<typename FPP, typename FPP2>
FaustCoreCpp<FPP>* fact_givens_fgft_generic_cplx(GivensFGFTComplex<FPP, Cpu, FPP2>* algo, FPP2* D, const int order, const bool enable_large_Faust)
FaustCoreCpp<FPP>* fact_givens_fgft_generic_cplx(GivensFGFTComplex<FPP, Cpu, FPP2>* algo, FPP2* D, const int order)
{
FaustCoreCpp<FPP>* fc = nullptr;
......
......@@ -11,7 +11,8 @@ cdef class FaustAlgoGenGivens@TYPE_NAME@:
@staticmethod
def fact_givens_fgft@FUNC_TYPE_SUFFIX@(Lap, J, t, verbosity=0, stoppingError = 0.0,
errIsRel=True, order=1, enable_large_Faust=False):
errIsRel=True, order=1, enable_large_Faust=False,
err_period=100):
if(order == 'ascend'): order = 1
elif(order == 'descend'): order = -1
elif(order == 'undef'): order = 0
......@@ -36,7 +37,8 @@ cdef class FaustAlgoGenGivens@TYPE_NAME@:
stoppingError,
errIsRel,
int(order),
enable_large_Faust)
enable_large_Faust,
err_period)
#from scipy.sparse import spdiags
#D_spdiag = spdiags(D, [0], Lap.shape[0], Lap.shape[0])
......@@ -50,7 +52,7 @@ cdef class FaustAlgoGenGivens@TYPE_NAME@:
@staticmethod
def fact_givens_fgft_sparse(Lap, J, t, verbosity=0, stoppingError = 0.0,
errIsRel=True, order='ascend',
enable_large_Faust=False):
enable_large_Faust=False, err_period=100):
from scipy.sparse import spdiags
cdef @TYPE@[:] data1d #only for csr mat factor
cdef int [:] indices # only for csr mat
......@@ -83,7 +85,8 @@ cdef class FaustAlgoGenGivens@TYPE_NAME@:
stoppingError,
errIsRel,
int(order),
enable_large_Faust)
enable_large_Faust,
err_period)
#D_spdiag = spdiags(D, [0], Lap.shape[0], Lap.shape[0])
#return core, D_spdiag
......@@ -94,8 +97,8 @@ cdef class FaustAlgoGenGivens@TYPE_NAME@:
return core, D
@staticmethod
def eigtj(M, nGivens=None, tol=0, relerr=True, nGivens_per_fac=None, verbosity=0,
order='ascend', enable_large_Faust=False):
def eigtj(M, nGivens=None, tol=0, relerr=True, nGivens_per_fac=None, verbosity=0,
order='ascend', enable_large_Faust=False, err_period=100):
if nGivens == None:
if tol == 0:
raise Exception("You must specify nGivens or tol argument"
......@@ -114,11 +117,11 @@ cdef class FaustAlgoGenGivens@TYPE_NAME@:
if isinstance(M, np.ndarray):
core_obj,D = FaustAlgoGenGivens@TYPE_NAME@.fact_givens_fgft@FUNC_TYPE_SUFFIX@(M, nGivens, nGivens_per_fac,
verbosity, tol,
relerr, order, enable_large_Faust)
relerr, order, enable_large_Faust, err_period)
elif isinstance(M, csr_matrix):
core_obj,D = FaustAlgoGenGivens@TYPE_NAME@.fact_givens_fgft_sparse@FUNC_TYPE_SUFFIX@(M, nGivens, nGivens_per_fac,
verbosity, tol,
relerr, order, enable_large_Faust)
relerr, order, enable_large_Faust, err_period)
else:
raise TypeError("The matrix to diagonalize must be a"
" scipy.sparse.csr_matrix or a numpy array.")
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment