diff --git a/misc/test/src/Python/test_FaustPy.py b/misc/test/src/Python/test_FaustPy.py
index 1dad4bf12ceb9196a076a8be27a7d7cb80985f40..4e5a63137696ac551b4c5e1f8ad2cf67a28fd405 100644
--- a/misc/test/src/Python/test_FaustPy.py
+++ b/misc/test/src/Python/test_FaustPy.py
@@ -906,7 +906,7 @@ class TestFaustFactory(unittest.TestCase):
L = L.astype(np.float64)
J = \
int(loadmat(sys.path[0]+"/../../../misc/data/mat/test_GivensDiag_Lap_U_J.mat")['J'])
- F, D = pyfaust.fact.fgft_givens(L, J, 0, verbosity=0)
+ F, D = pyfaust.fact.fgft_givens(L, J, nGivens_per_fac=1, verbosity=0)
print("Lap norm:", norm(L, 'fro'))
err = norm((F*D.todense())*F.T.todense()-L,"fro")/norm(L,"fro")
print("err: ", err)
@@ -922,14 +922,14 @@ class TestFaustFactory(unittest.TestCase):
J = \
int(loadmat(sys.path[0]+"/../../../misc/data/mat/test_GivensDiag_Lap_U_J.mat")['J'])
t = int(L.shape[0]/2)
- F, D = fgft_givens(L, J, t, verbosity=0)
+ F, D = fgft_givens(L, J, nGivens_per_fac=t, verbosity=0)
print("Lap norm:", norm(L, 'fro'))
err = norm((F*D.todense())*F.T.todense()-L,"fro")/norm(L,"fro")
print("err: ", err)
# the error reference is from the C++ test,
# misc/test/src/C++/GivensFGFTParallel.cpp.in
self.assertAlmostEqual(err, 0.0410448, places=7)
- F2, D2 = eigtj(L, J, t, verbosity=0)
+ F2, D2 = eigtj(L, J, nGivens_per_fac=t, verbosity=0)
print("Lap norm:", norm(L, 'fro'))
err2 = norm((F2*D.todense())*F2.T.todense()-L,"fro")/norm(L,"fro")
print("err2: ", err2)
@@ -946,7 +946,7 @@ class TestFaustFactory(unittest.TestCase):
L = L.astype(np.float64)
J = \
int(loadmat(sys.path[0]+"/../../../misc/data/mat/test_GivensDiag_Lap_U_J.mat")['J'])
- F, D = pyfaust.fact.fgft_givens(csr_matrix(L), J, 0, verbosity=0)
+ F, D = pyfaust.fact.fgft_givens(csr_matrix(L), J, nGivens_per_fac=0, verbosity=0)
print("Lap norm:", norm(L, 'fro'))
err = norm((F*D.todense())*F.T.todense()-L,"fro")/norm(L,"fro")
print("err: ", err)
@@ -955,7 +955,7 @@ class TestFaustFactory(unittest.TestCase):
self.assertAlmostEqual(err, 0.0326529, places=7)
def testFGFTGivensParallelSparse(self):
- print("Test pyfaust.fact.fgft_givens() -- parallel")
+ print("Test pyfaust.fact.fgft_givens_sparse() -- parallel")
from pyfaust.fact import eigtj, fgft_givens
from scipy.sparse import csr_matrix
L = loadmat(sys.path[0]+"/../../../misc/data/mat/test_GivensDiag_Lap_U_J.mat")['Lap']
@@ -963,14 +963,14 @@ class TestFaustFactory(unittest.TestCase):
J = \
int(loadmat(sys.path[0]+"/../../../misc/data/mat/test_GivensDiag_Lap_U_J.mat")['J'])
t = int(L.shape[0]/2)
- F, D = fgft_givens(csr_matrix(L), J, t, verbosity=0)
+ F, D = fgft_givens(csr_matrix(L), J, nGivens_per_fac=t, verbosity=0)
print("Lap norm:", norm(L, 'fro'))
err = norm((F*D.todense())*F.T.todense()-L,"fro")/norm(L,"fro")
print("err: ", err)
# the error reference is from the C++ test,
# misc/test/src/C++/GivensFGFTParallel.cpp.in
self.assertAlmostEqual(err, 0.0410448, places=7)
- F2, D2 = eigtj(L, J, t, verbosity=0)
+ F2, D2 = eigtj(L, J, nGivens_per_fac=t, verbosity=0)
print("Lap norm:", norm(L, 'fro'))
err2 = norm((F2*D.todense())*F2.T.todense()-L,"fro")/norm(L,"fro")
print("err2: ", err2)
diff --git a/src/algorithm/factorization/faust_GivensFGFT.h b/src/algorithm/factorization/faust_GivensFGFT.h
index 85fffec2919c47a87983b600df9e7a206044e595..073e800fcad25c549add5931b484004ada8fb872 100644
--- a/src/algorithm/factorization/faust_GivensFGFT.h
+++ b/src/algorithm/factorization/faust_GivensFGFT.h
@@ -35,6 +35,7 @@ namespace Faust {
/** \brief Pivot candidates q coordinates. */
int* q_candidates; /* default IndexType for underlying eigen matrix is int. */
protected:
+ const static unsigned int ERROR_CALC_PERIOD = 100;
/** \brief Fourier matrix factorization matrices (Givens matrix). */
vector<Faust::MatSparse<FPP,DEVICE>> facts;
/** \brief Diagonalization approximate of Laplacian. */
@@ -88,6 +89,10 @@ namespace Faust {
/** \brief In calc_theta() two values are calculated for theta, this boolean is set to true to always choose theta2 (useful for GivensFGFTParallel). */
bool always_theta2;
+ bool stoppingCritIsError;
+ double stoppingError;
+ bool errIsRel;
+
/**
* Function pointer to any step of the algorithm (internal purpose only).
*/
@@ -100,8 +105,8 @@ namespace Faust {
* \param Lap The Laplacian matrix to approximate/diagonalize.
* \param J The number of iterations, Givens rotations factors.
* */
- GivensFGFT(Faust::MatSparse<FPP,DEVICE>& Lap, int J, unsigned int verbosity = 0);
- GivensFGFT(Faust::MatDense<FPP,DEVICE>& Lap, int J, unsigned int verbosity = 0);
+ GivensFGFT(Faust::MatSparse<FPP,DEVICE>& Lap, int J, unsigned int verbosity = 0, const double stoppingError = 0.0, const bool errIsRel = true);
+ GivensFGFT(Faust::MatDense<FPP,DEVICE>& Lap, int J, unsigned int verbosity = 0, const double stoppingError = 0.0, const bool errIsRel = true);
/** Destructor */
virtual ~GivensFGFT() {delete[] q_candidates; delete L;};
diff --git a/src/algorithm/factorization/faust_GivensFGFT.hpp b/src/algorithm/factorization/faust_GivensFGFT.hpp
index 55aada53957dfbe878c9023331ec4f112a447065..3bf3bc78d50f5c03282eb1de54115c9c6766af57 100644
--- a/src/algorithm/factorization/faust_GivensFGFT.hpp
+++ b/src/algorithm/factorization/faust_GivensFGFT.hpp
@@ -371,7 +371,7 @@ void GivensFGFT<FPP,DEVICE,FPP2>::update_err()
// % disp(['Number of edges: ' num2str(N_edges)])
// end
//
- if(!((ite+1)%100) || verbosity > 1)
+ if(!((ite+1)%GivensFGFT<FPP,DEVICE,FPP2>::ERROR_CALC_PERIOD) || stoppingCritIsError || verbosity > 1)
{
MatDense<FPP,DEVICE> tmp = this->get_Dspm(false);
MatDense<FPP,DEVICE>* dL;
@@ -388,11 +388,14 @@ void GivensFGFT<FPP,DEVICE,FPP2>::update_err()
}
FPP2 err = tmp.norm(), err_d;
err *= err;
- err_d = Lap.norm();
- err_d *= err_d;
- err /= err_d;
+ if(errIsRel)
+ {
+ err_d = Lap.norm();
+ err_d *= err_d;
+ err /= err_d;
+ }
if(verbosity)
- cout << "GivensFGFT ite. i: "<< ite << " err.: " << err << endl;
+ cout << "GivensFGFT factor : "<< ite << ", transform " << ((errIsRel)?"relative ":"absolute ") << "err.: " << err << endl;
errs.push_back(err);
}
}
@@ -430,11 +433,16 @@ void GivensFGFT<FPP,DEVICE,FPP2>::compute_facts()
{
next_step();
ite++;
+ if(stoppingCritIsError && *(errs.end()-1) <= stoppingError)
+ {
+ facts.resize(ite);
+ break;
+ }
}
}
template<typename FPP, Device DEVICE, typename FPP2>
-GivensFGFT<FPP,DEVICE,FPP2>::GivensFGFT(Faust::MatSparse<FPP,DEVICE>& Lap, int J, unsigned int verbosity /* deft val == 0 */) : Lap(Lap), facts(J), D(Lap.getNbRow()), C(Lap.getNbRow(), Lap.getNbCol()), errs(0), coord_choices(J), q_candidates(new int[Lap.getNbCol()]), is_D_ordered(false), always_theta2(false), verbosity(verbosity)
+GivensFGFT<FPP,DEVICE,FPP2>::GivensFGFT(Faust::MatSparse<FPP,DEVICE>& Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stopingError /* default to 0.0 */, const bool errIsRel) : Lap(Lap), facts(J), D(Lap.getNbRow()), C(Lap.getNbRow(), Lap.getNbCol()), errs(0), coord_choices(J), q_candidates(new int[Lap.getNbCol()]), is_D_ordered(false), always_theta2(false), verbosity(verbosity), stoppingCritIsError(stoppingError != 0.0), stoppingError(stoppingError), errIsRel(errIsRel)
{
/* Matlab ref. code:
* facts = cell(1,J);
@@ -466,7 +474,7 @@ GivensFGFT<FPP,DEVICE,FPP2>::GivensFGFT(Faust::MatSparse<FPP,DEVICE>& Lap, int J
}
template<typename FPP, Device DEVICE, typename FPP2>
-GivensFGFT<FPP,DEVICE,FPP2>::GivensFGFT(Faust::MatDense<FPP,DEVICE>& Lap, int J, unsigned int verbosity /* deft val == 0 */) : Lap(Lap), facts(J), D(Lap.getNbRow()), C(Lap.getNbRow(), Lap.getNbCol()), errs(0), coord_choices(J), q_candidates(new int[Lap.getNbCol()]), is_D_ordered(false), always_theta2(false), verbosity(verbosity)
+GivensFGFT<FPP,DEVICE,FPP2>::GivensFGFT(Faust::MatDense<FPP,DEVICE>& Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError, const bool errIsRel) : Lap(Lap), facts(J), D(Lap.getNbRow()), C(Lap.getNbRow(), Lap.getNbCol()), errs(0), coord_choices(J), q_candidates(new int[Lap.getNbCol()]), is_D_ordered(false), always_theta2(false), verbosity(verbosity), stoppingCritIsError(stoppingError != 0.0), stoppingError(stoppingError), errIsRel(errIsRel)
{
/* Matlab ref. code:
* facts = cell(1,J);
diff --git a/src/algorithm/factorization/faust_GivensFGFTParallel.h b/src/algorithm/factorization/faust_GivensFGFTParallel.h
index 34d2eb761b40b974c08f24746c1d6939fe612a5f..8e51c7f426973ebef5dbbb52cca5ed309b88a1ac 100644
--- a/src/algorithm/factorization/faust_GivensFGFTParallel.h
+++ b/src/algorithm/factorization/faust_GivensFGFTParallel.h
@@ -81,10 +81,11 @@ namespace Faust {
* \param Lap The Laplacian matrix to approximate/diagonalize.
* \param J round(J/t) is the number of factors to approximate the Fourier matrix (as a product of factors).
* \param t the maximum number of 2D rotation matrices (Givens matrices) to insert in each factor. The effective number can be less than t if there is not enough pivot candidates left in the matrix L of the current iteration.
+ * \param stoppingError defines a stopping criterion based on error (absolute relative error).
*
*/
- GivensFGFTParallel(Faust::MatDense<FPP,DEVICE>& Lap, int J, int t, unsigned int verbosity = 0);
- GivensFGFTParallel(Faust::MatSparse<FPP,DEVICE>& Lap, int J, int t, unsigned int verbosity = 0);
+ GivensFGFTParallel(Faust::MatDense<FPP,DEVICE>& Lap, int J, int t, unsigned int verbosity = 0, const double stoppingCritIsError = 0.0, const bool errIsRel = true);
+ GivensFGFTParallel(Faust::MatSparse<FPP,DEVICE>& Lap, int J, int t, unsigned int verbosity = 0, const double stoppingCritIsError = 0.0, const bool errIsRel = true);
};
diff --git a/src/algorithm/factorization/faust_GivensFGFTParallel.hpp b/src/algorithm/factorization/faust_GivensFGFTParallel.hpp
index 56eacdc99db052e96db7dca5f214f485c72cf922..46a7388ea7b1715f8d5b93d9dff14fc24d7eb2ca 100644
--- a/src/algorithm/factorization/faust_GivensFGFTParallel.hpp
+++ b/src/algorithm/factorization/faust_GivensFGFTParallel.hpp
@@ -6,7 +6,7 @@ using namespace Faust;
#endif
template<typename FPP, Device DEVICE, typename FPP2>
-GivensFGFTParallel<FPP,DEVICE,FPP2>::GivensFGFTParallel(Faust::MatDense<FPP,DEVICE>& Lap, int J, int t, unsigned int verbosity) : GivensFGFT<FPP,DEVICE,FPP2>(Lap, J, verbosity), t(t), fact_nrots(0)
+GivensFGFTParallel<FPP,DEVICE,FPP2>::GivensFGFTParallel(Faust::MatDense<FPP,DEVICE>& Lap, int J, int t, unsigned int verbosity, const double stoppingError, const bool errIsRel) : GivensFGFT<FPP,DEVICE,FPP2>(Lap, J, verbosity, stoppingError, errIsRel), t(t), fact_nrots(0)
{
this->facts.resize(round(J/(float)t));
this->always_theta2 = true;
@@ -14,7 +14,7 @@ GivensFGFTParallel<FPP,DEVICE,FPP2>::GivensFGFTParallel(Faust::MatDense<FPP,DEVI
}
template<typename FPP, Device DEVICE, typename FPP2>
-GivensFGFTParallel<FPP,DEVICE,FPP2>::GivensFGFTParallel(Faust::MatSparse<FPP,DEVICE>& Lap, int J, int t, unsigned int verbosity) : GivensFGFT<FPP,DEVICE,FPP2>(Lap, J, verbosity), t(t), fact_nrots(0)
+GivensFGFTParallel<FPP,DEVICE,FPP2>::GivensFGFTParallel(Faust::MatSparse<FPP,DEVICE>& Lap, int J, int t, unsigned int verbosity, const double stoppingError, const bool errIsRel) : GivensFGFT<FPP,DEVICE,FPP2>(Lap, J, verbosity, stoppingError, errIsRel), t(t), fact_nrots(0)
{
this->facts.resize(round(J/(float)t));
this->always_theta2 = true;
diff --git a/wrapper/python/pyfaust/fact.py b/wrapper/python/pyfaust/fact.py
index c81bc16ec3b683da20267db6c13e62beffe8265d..1d81184454fb2329849086343113d3f720cd7cf9 100644
--- a/wrapper/python/pyfaust/fact.py
+++ b/wrapper/python/pyfaust/fact.py
@@ -78,27 +78,27 @@ def svdtj(M, J, t=1):
V = W2[:,0:S.shape[0]]*Faust(Id[:,I])
return U,S,V
-def eigtj(M, J, t=1, verbosity=0):
+def eigtj(M, maxiter, tol=0, relerr=True, nGivens_per_fac=None, verbosity=0):
"""
- Computes the eigenvalues and the eigenvectors transform (as a Faust object) using the truncated Jacobi algorithm.
+ Runs the truncated Jacobi algorithm to compute the eigenvalues of M (returned in W) and the corresponding eigenvectors (in Faust V) using the truncated Jacobi algorithm.
+
+ The output is such that V*W.todense()*V.T approximates M.
- The eigenvalues and the eigenvectors are approximate.
The trade-off between accuracy and sparsity can be set through the
- parameters J and t.
+ parameters maxiter and nGivens_per_fac.
Args:
M: (numpy.ndarray) the matrix to diagonalize. Must be real and symmetric.
- J: (int) defines the number of factors in the eigenvector transform V.
- The number of factors is round(J/t). Note that the last permutation factor
- is not in count here (in fact, the total number of factors in V is rather
- round(J/t)+1).
- t: (int) the number of Givens rotations per factor. Note that t is
- forced to the value min(J,t). Besides, a value of t such that t >
- M.shape[0]/2 won't lead to the desired effect because the maximum
- number of rotation matrices per factor is anyway M.shape[0]/2.
- The parameter t is meaningful in the parallel version of the
- truncated Jacobi algorithm (cf. references below). If t <= 1 (by default)
- then the function runs the non-parallel algorithm.
+ maxiter: (int) defines the number of Givens rotations that are computed in
+ eigenvector transform V. The number of rotations per factor of V is
+ defined by nGivens_per_fac.
+ tol: (float) the tolerance error under what the algorithm stops. By default,
+ it's zero for not stopping on error criterion.
+ relErr: (bool) the type of error stopping criterion. Default to True to use
+ relative error, otherwise (False) the absolute error is used.
+ nGivens_per_fac: (int) the number of Givens rotations per factor of V, must be
+ an integer between 1 to M.shape[0]/2 which is the default value (when
+ nGivens_per_fac == None).
verbosity: (int) the level of verbosity, the greater the value the more info.
is displayed.
@@ -106,11 +106,11 @@ def eigtj(M, J, t=1, verbosity=0):
Returns:
The tuple (V,W):
- V the Faust object representing the approximate eigenvector
- transform. V has its last factor being a permutation matrix,
- the goal of this factor is to apply to the columns of V the same order as
+ transform. The last factor of V is a permutation matrix.
+ The goal of this factor is to apply to the columns of V the same order as
eigenvalues set in W.
- - W (scipy.sparse.dia.dia_matrix) the approximate diagonal matrix of the eigenvalues
- (in ascendant order along the rows/columns).
+ - W (scipy.sparse.dia.dia_matrix) the diagonal matrix of the
+ approximate eigenvalues (in ascendant order along the rows/columns).
References:
[1] Le Magoarou L., Gribonval R. and Tremblay N., "Approximate fast
@@ -129,8 +129,8 @@ def eigtj(M, J, t=1, verbosity=0):
demo_path = sep.join((get_data_dirpath(),'Laplacian_256_community.mat'))
data_dict = loadmat(demo_path)
Lap = data_dict['Lap'].astype(np.float)
- Uhat, Dhat = eigtj(Lap,J=Lap.shape[0]*100,t=int(Lap.shape[0]/2))
- # Uhat is the Fourier matrix/eigenvectors approximattion as a Faust
+ Uhat, Dhat = eigtj(Lap, maxiter=Lap.shape[0]*100)
+ # Uhat is the Fourier matrix/eigenvectors approximation as a Faust
# (200 factors + permutation mat.)
# Dhat the eigenvalues diagonal matrix approx.
print(Uhat)
@@ -139,16 +139,18 @@ def eigtj(M, J, t=1, verbosity=0):
See also:
fgft_givens, fgft_palm
"""
- return fgft_givens(M, J, t, verbosity)
+ return fgft_givens(M, maxiter, tol, relerr, nGivens_per_fac, verbosity)
-def fgft_givens(Lap, J, t=1, verbosity=0):
+def fgft_givens(Lap, maxiter, tol=0.0, relerr=True, nGivens_per_fac=None, verbosity=0):
"""
Diagonalizes the graph Laplacian matrix Lap using the Givens FGFT algorithm.
Args:
Lap: the Laplacian matrix as a numpy array. Must be real and symmetric.
- J: see eigtj
- t: see eigtj
+ maxiter: see eigtj
+ tol: see eigtj
+ relerr: see eigtj
+ nGivens_per_fac: see eigtj
verbosity: see eigtj
Returns:
@@ -170,17 +172,21 @@ def fgft_givens(Lap, J, t=1, verbosity=0):
"""
if((isinstance(Lap, np.ndarray) and (Lap.T != Lap).any()) or Lap.shape[0] != Lap.shape[1]):
raise ValueError(" the matrix/array must be square and should be symmetric.")
- if(not isinstance(J, int)): raise TypeError("J must be a int")
- if(J < 1): raise ValueError("J must be >= 1")
- if(not isinstance(t, int)): raise TypeError("t must be a int")
- if(t < 0): raise ValueError("t must be >= 0")
- t = min(t,J)
+ if(not isinstance(maxiter, int)): raise TypeError("maxiter must be a int")
+ if(maxiter < 1): raise ValueError("maxiter must be >= 1")
+ if(not isinstance(nGivens_per_fac, int)): raise TypeError("nGivens_per_fac must be a int")
+ nGivens_per_fac = max(nGivens_per_fac, 1)
+ nGivens_per_fac = min(nGivens_per_fac, maxiter)
if(isinstance(Lap, np.ndarray)):
- core_obj,D = _FaustCorePy.FaustFact.fact_givens_fgft(Lap, J, t,
- verbosity)
+ core_obj,D = _FaustCorePy.FaustFact.fact_givens_fgft(Lap, maxiter, nGivens_per_fac,
+ verbosity, tol,
+ relerr)
elif(isinstance(Lap, csr_matrix)):
- core_obj,D = _FaustCorePy.FaustFact.fact_givens_fgft_sparse(Lap, J, t,
- verbosity)
+ core_obj,D = _FaustCorePy.FaustFact.fact_givens_fgft_sparse(Lap, maxiter,
+ nGivens_per_fac,
+ verbosity,
+ tol,
+ relerr)
else:
raise TypeError("The matrix to diagonalize must be a"
" scipy.sparse.csr_matrix or a numpy array.")
@@ -429,6 +435,7 @@ def _prepare_hierarchical_fact(M, p, callee_name, ret_lambda, ret_params,
"with the first factor constraint defined in p.")
return p
+# experimental block start
def hierarchical_constends(M, p, A, B, ret_lambda=False, ret_params=False):
"""
Tries to approximate M by \f$ A \prod_j S_j B \f$ using the hierarchical.
@@ -510,7 +517,9 @@ def hierarchical_constends(M, p, A, B, ret_lambda=False, ret_params=False):
return ret_list
else:
return F
+# experimental block end
+# experimental block start
def palm4msa_constends(M, p, A, B=None, ret_lambda=False):
"""
Tries to approximate M by \f$ A \prod_j S_j B\f$ using palm4msa (B being optional).
@@ -564,6 +573,7 @@ def palm4msa_constends(M, p, A, B=None, ret_lambda=False):
return F, _lambda
else:
return F
+# experimental block end
def fgft_palm(U, Lap, p, init_D=None, ret_lambda=False, ret_params=False):
"""
diff --git a/wrapper/python/src/FaustCoreCy.pxd b/wrapper/python/src/FaustCoreCy.pxd
index 64a342eb582465e4c81d8b651039e4f52378bca8..410bc797a932af92b35efb9d71ec279bd2e76449 100644
--- a/wrapper/python/src/FaustCoreCy.pxd
+++ b/wrapper/python/src/FaustCoreCy.pxd
@@ -184,9 +184,14 @@ cdef extern from "FaustFactGivensFGFT.h":
cdef FaustCoreCpp[FPP]* fact_givens_fgft[FPP,FPP2](const FPP* Lap, unsigned int num_rows,
unsigned int num_cols, unsigned int J,
- unsigned int t, FPP* D, unsigned int verbosity)
+ unsigned int t, FPP* D, unsigned int verbosity,
+ const double stoppingError,
+ const bool errIsRel)
cdef FaustCoreCpp[FPP]* fact_givens_fgft_sparse[FPP,FPP2](const FPP* data, int* row_ptr,
int* id_col, int nnz, unsigned int num_rows,
unsigned int num_cols, unsigned int J,
- unsigned int t, FPP* D, unsigned int verbosity)
+ unsigned int t, FPP* D,
+ unsigned int verbosity,
+ const double stoppingError,
+ const bool errIsRel)
diff --git a/wrapper/python/src/FaustFactGivensFGFT.h b/wrapper/python/src/FaustFactGivensFGFT.h
index 013c590b832a55716e8e1cdb3ca5d255d1db6ef8..c5a86f8861b7e2e43a6a5653d9cb677994868bdb 100644
--- a/wrapper/python/src/FaustFactGivensFGFT.h
+++ b/wrapper/python/src/FaustFactGivensFGFT.h
@@ -3,11 +3,11 @@
#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 /* end of input parameters*/, FPP* D, unsigned int verbosity = 0);
+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 = 0, const double stoppingError = 0.0, const bool errIsRel = true);
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);
+ unsigned int J, unsigned int t /* end of input parameters*/, FPP* D, unsigned int verbosity = 0, const double stoppingError = 0.0, const bool errIsRel = true);
template<typename FPP, typename FPP2 = float>
diff --git a/wrapper/python/src/FaustFactGivensFGFT.hpp b/wrapper/python/src/FaustFactGivensFGFT.hpp
index 13a2568a7f4b326bdca66ecde1a146c44882ec46..4e04497890846582bf90c16e363de165b8485e21 100644
--- a/wrapper/python/src/FaustFactGivensFGFT.hpp
+++ b/wrapper/python/src/FaustFactGivensFGFT.hpp
@@ -6,23 +6,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)
+ unsigned int J, unsigned int t /* end of input parameters*/, FPP* D, unsigned int verbosity, const double stoppingError, const bool errIsRel)
{
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);
+ algo = new GivensFGFT<FPP, Cpu, FPP2>(mat_Lap, (int)J, verbosity, stoppingError, errIsRel);
}
else
{
- algo = new GivensFGFTParallel<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity);
+ algo = new GivensFGFTParallel<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity, stoppingError, errIsRel);
}
return fact_givens_fgft_generic(algo, D);
}
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)
+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 double stoppingError, const bool errIsRel)
{
//TODO: optimization possible here by avoiding Lap copy in MatDense (by
//just using the data in Lap as underlying pointer of MatDense)
@@ -31,11 +31,11 @@ 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);
+ algo = new GivensFGFT<FPP, Cpu, FPP2>(mat_Lap, (int)J, verbosity, stoppingError, errIsRel);
}
else
{
- algo = new GivensFGFTParallel<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity);
+ algo = new GivensFGFTParallel<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity, stoppingError, errIsRel);
}
return fact_givens_fgft_generic(algo, D);
}
diff --git a/wrapper/python/src/_FaustCorePy.pyx b/wrapper/python/src/_FaustCorePy.pyx
index f0a6e431d86cc3ef3e743946954d505cacdb3126..768f14d5db658c7621d29aeca014d00f90609dab 100644
--- a/wrapper/python/src/_FaustCorePy.pyx
+++ b/wrapper/python/src/_FaustCorePy.pyx
@@ -1350,7 +1350,8 @@ cdef class FaustFact:
return core, np.real(_out_buf[0])
@staticmethod
- def fact_givens_fgft_sparse(Lap, J, t, verbosity=0):
+ def fact_givens_fgft_sparse(Lap, J, t, verbosity=0, stoppingError = 0.0,
+ errIsRel=True):
from scipy.sparse import spdiags
cdef double [:] data1d #only for csr mat factor
cdef int [:] indices # only for csr mat
@@ -1374,14 +1375,18 @@ cdef class FaustFact:
Lap.nnz,
Lap_num_rows,
Lap_num_cols, J, t,
- &D_view[0], verbosity)
+ &D_view[0],
+ verbosity,
+ stoppingError,
+ errIsRel)
core._isReal = True
D_spdiag = spdiags(D, [0], Lap.shape[0], Lap.shape[0])
return core, D_spdiag
@staticmethod
- def fact_givens_fgft(Lap, J, t, verbosity=0):
+ def fact_givens_fgft(Lap, J, t, verbosity=0, stoppingError = 0.0,
+ errIsRel=True):
isReal = Lap.dtype in [ 'float', 'float128',
'float16', 'float32',
'float64', 'double']
@@ -1403,9 +1408,11 @@ cdef class FaustFact:
core = FaustCore(core=True)
core.core_faust_dbl = FaustCoreCy.fact_givens_fgft[double, double](&Lap_view[0,0],
- Lap_num_rows,
- Lap_num_cols, J, t,
- &D_view[0], verbosity)
+ Lap_num_rows,
+ Lap_num_cols, J, t,
+ &D_view[0], verbosity,
+ stoppingError,
+ errIsRel)
core._isReal = True
from scipy.sparse import spdiags