diff --git a/misc/test/src/Matlab/FaustFactoryTest.m b/misc/test/src/Matlab/FaustFactoryTest.m
index 80134f293b0034a19d7cbadb3b46c81b456eee1b..eb479f24a2e66598f5650d863785f7fe204fafbc 100644
--- a/misc/test/src/Matlab/FaustFactoryTest.m
+++ b/misc/test/src/Matlab/FaustFactoryTest.m
@@ -163,7 +163,7 @@ classdef FaustFactoryTest < matlab.unittest.TestCase
 			import matfaust.*
 			load([this.faust_paths{1} '../../../misc/data/mat/test_GivensDiag_Lap_U_J.mat'])
 			% Lap and J available
-			[F,D] = matfaust.fact.fgft_givens(Lap, 'nGivens', J);%, 0, 'verbosity', 1);
+			[F,D] = matfaust.fact.fgft_givens(Lap, 'nGivens', J, 'enable_large_Faust', true, 'nGivens_per_fac', 1);%, 0, 'verbosity', 1);
 			this.verifyEqual(size(F), size(Lap))
 			%disp('norm F: ')
 			%norm(F, 'fro')
@@ -173,7 +173,7 @@ classdef FaustFactoryTest < matlab.unittest.TestCase
 			% misc/test/src/C++/GivensFGFT.cpp.in
 			this.verifyEqual(err, 0.0326529, 'AbsTol', 0.00001)
 			% verify it works the same using the eigtj() alias function
-			[F2,D2] = matfaust.fact.eigtj(Lap, 'nGivens', J);%, 0, 'verbosity', 2);
+			[F2,D2] = matfaust.fact.eigtj(Lap, 'nGivens', J, 'enable_large_Faust', true, 'nGivens_per_fac', 1);%, 0, 'verbosity', 2);
 			this.verifyEqual(full(F2),full(F))
 			this.verifyEqual(D,D2)
 		end
@@ -183,7 +183,7 @@ classdef FaustFactoryTest < matlab.unittest.TestCase
 			import matfaust.*
 			load([this.faust_paths{1} '../../../misc/data/mat/test_GivensDiag_Lap_U_J.mat'])
 			% Lap and J available
-			[F,D] = matfaust.fact.fgft_givens(sparse(Lap), 'nGivens', J);%, 0, 'verbosity', 1);
+			[F,D] = matfaust.fact.fgft_givens(sparse(Lap), 'nGivens', J, 'enable_large_Faust', true, 'nGivens_per_fac', 1);%, 0, 'verbosity', 1);
 			this.verifyEqual(size(F), size(Lap))
 			%disp('norm F: ')
 			%norm(F, 'fro')
@@ -193,18 +193,18 @@ classdef FaustFactoryTest < matlab.unittest.TestCase
 			% misc/test/src/C++/GivensFGFT.cpp.in
 			this.verifyEqual(err, 0.0326529, 'AbsTol', 0.00001)
 			% verify it works the same using the eigtj() alias function
-			[F2,D2] = matfaust.fact.eigtj(Lap, 'nGivens', J);%, 0, 'verbosity', 2);
+			[F2,D2] = matfaust.fact.eigtj(Lap, 'nGivens', J, 'nGivens_per_fac', 1, 'enable_large_Faust', true);%, 0, 'verbosity', 2);
 			this.verifyEqual(full(F2),full(F))
 			this.verifyEqual(D,D2)
 		end
 
 		function test_fgft_givens_parallel(this)
-			disp('Test matfaust.fact.fgft_givens_sparse() -- parallel')
+			disp('Test matfaust.fact.fgft_givens() -- parallel')
 			import matfaust.*
 			load([this.faust_paths{1} '../../../misc/data/mat/test_GivensDiag_Lap_U_J.mat'])
 			% Lap and J available
 			t = size(Lap,1)/2;
-			[F,D] = matfaust.fact.fgft_givens(Lap, 'nGivens', J, 'nGivens_per_fac', t); %, 'verbosity', 2);
+			[F,D] = matfaust.fact.fgft_givens(Lap, 'nGivens', J, 'nGivens_per_fac', t, 'enable_large_Faust', true); %, 'verbosity', 2);
 			this.verifyEqual(size(F), size(Lap))
 			%disp('norm F: ')
 			%norm(F, 'fro')
@@ -214,7 +214,7 @@ classdef FaustFactoryTest < matlab.unittest.TestCase
 			% misc/test/src/C++/GivensFGFTParallel.cpp.in
 			this.verifyEqual(err,0.0398154, 'AbsTol', 0.00001)
 			% verify it works the same using the eigtj() alias function
-			[F2,D2] = matfaust.fact.eigtj(Lap, 'nGivens', J, 'nGivens_per_fac', t); %, 'verbosity', 2);
+			[F2,D2] = matfaust.fact.eigtj(Lap, 'nGivens', J, 'nGivens_per_fac', t, 'enable_large_Faust', true); %, 'verbosity', 2);
 			this.verifyEqual(full(F2),full(F))
 			this.verifyEqual(D,D2)
 		end
@@ -225,7 +225,7 @@ classdef FaustFactoryTest < matlab.unittest.TestCase
 			load([this.faust_paths{1} '../../../misc/data/mat/test_GivensDiag_Lap_U_J.mat'])
 			% Lap and J available
 			t = size(Lap,1)/2;
-			[F,D] = matfaust.fact.fgft_givens(sparse(Lap), 'nGivens', J, 'nGivens_per_fac', t); %, 'verbosity', 2);
+			[F,D] = matfaust.fact.fgft_givens(sparse(Lap), 'nGivens', J, 'nGivens_per_fac', t, 'enable_large_Faust', true); %, 'verbosity', 2);
 			this.verifyEqual(size(F), size(Lap))
 			%disp('norm F: ')
 			%norm(F, 'fro')
@@ -235,7 +235,7 @@ classdef FaustFactoryTest < matlab.unittest.TestCase
 			% misc/test/src/C++/GivensFGFTParallel.cpp.in
 			this.verifyEqual(err, 0.0398154, 'AbsTol', 0.00001)
 			% verify it works the same using the eigtj() alias function
-			[F2,D2] = matfaust.fact.eigtj(Lap, 'nGivens', J, 'nGivens_per_fac', t); %, 'verbosity', 2);
+			[F2,D2] = matfaust.fact.eigtj(Lap, 'nGivens', J, 'nGivens_per_fac', t, 'enable_large_Faust', true); %, 'verbosity', 2);
 			this.verifyEqual(full(F2),full(F))
 			this.verifyEqual(D,D2)
 		end
diff --git a/misc/test/src/Python/test_FaustPy.py b/misc/test/src/Python/test_FaustPy.py
index 3ddb1bfbed89ccbc1aeef51e9ba2da8cfbf09cd1..f34c882391fc66cd45c84dfb196c827dae67e798 100644
--- a/misc/test/src/Python/test_FaustPy.py
+++ b/misc/test/src/Python/test_FaustPy.py
@@ -925,7 +925,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'])
-        D, F = pyfaust.fact.fgft_givens(L, J, nGivens_per_fac=1, verbosity=0)
+        D, F = pyfaust.fact.fgft_givens(L, J, nGivens_per_fac=1, verbosity=0, enable_large_Faust=True)
         D = spdiags(D, [0], L.shape[0], L.shape[0])
         print("Lap norm:", norm(L, 'fro'))
         err = norm((F*D.todense())*F.T.todense()-L,"fro")/norm(L,"fro")
@@ -946,7 +946,8 @@ 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)
-        D, F = fgft_givens(L, J, nGivens_per_fac=t, verbosity=0)
+        D, F = fgft_givens(L, J, nGivens_per_fac=t, verbosity=0,
+                           enable_large_Faust=True)
         D = spdiags(D, [0], L.shape[0], L.shape[0])
         print("Lap norm:", norm(L, 'fro'))
         err = norm((F*D.todense())*F.T.todense()-L,"fro")/norm(L,"fro")
@@ -954,7 +955,8 @@ class TestFaustFactory(unittest.TestCase):
         # the error reference is from the C++ test,
         # misc/test/src/C++/GivensFGFTParallel.cpp.in (_double version)
         self.assertAlmostEqual(err, 0.0398154, places=7)
-        D2, F2 = eigtj(L, J, nGivens_per_fac=t, verbosity=0)
+        D2, F2 = eigtj(L, J, nGivens_per_fac=t, verbosity=0,
+                       enable_large_Faust=True)
         D2 = spdiags(D, [0], L.shape[0], L.shape[0])
         print("Lap norm:", norm(L, 'fro'))
         err2 = norm((F2*D.todense())*F2.T.todense()-L,"fro")/norm(L,"fro")
@@ -977,7 +979,8 @@ 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'])
-        D, F = pyfaust.fact.fgft_givens(csr_matrix(L), J, nGivens_per_fac=0, verbosity=0)
+        D, F = pyfaust.fact.fgft_givens(csr_matrix(L), J, nGivens_per_fac=0,
+                                        verbosity=0,enable_large_Faust=True)
         D = spdiags(D, [0], L.shape[0], L.shape[0])
         print("Lap norm:", norm(L, 'fro'))
         err = norm((F*D.todense())*F.T.todense()-L,"fro")/norm(L,"fro")
@@ -999,7 +1002,7 @@ 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)
-        D, F = fgft_givens(csr_matrix(L), J, nGivens_per_fac=t, verbosity=0)
+        D, F = fgft_givens(csr_matrix(L), J, nGivens_per_fac=t, verbosity=0, enable_large_Faust=True)
         D = spdiags(D, [0], L.shape[0], L.shape[0])
         print("Lap norm:", norm(L, 'fro'))
         err = norm((F*D.todense())*F.T.todense()-L,"fro")/norm(L,"fro")
@@ -1007,7 +1010,8 @@ class TestFaustFactory(unittest.TestCase):
         # the error reference is from the C++ test,
         # misc/test/src/C++/GivensFGFTParallel.cpp.in (_double version)
         self.assertAlmostEqual(err, 0.0398154, places=7)
-        D2, F2 = eigtj(csr_matrix(L), J, nGivens_per_fac=t, verbosity=0)
+        D2, F2 = eigtj(csr_matrix(L), J, nGivens_per_fac=t, verbosity=0,
+                       enable_large_Faust=True)
         D2 = spdiags(D, [0], L.shape[0], L.shape[0])
         print("Lap norm:", norm(L, 'fro'))
         err2 = norm((F2*D.todense())*F2.T.todense()-L,"fro")/norm(L,"fro")
diff --git a/src/algorithm/factorization/faust_GivensFGFT.h b/src/algorithm/factorization/faust_GivensFGFT.h
index 255100f3243241e557a1f8387ccf67d15af09826..08a0c22daa1507adfb2fdd9e838e8adc6e2eac59 100644
--- a/src/algorithm/factorization/faust_GivensFGFT.h
+++ b/src/algorithm/factorization/faust_GivensFGFT.h
@@ -56,8 +56,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, 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);
+				GivensFGFT(Faust::MatSparse<FPP,DEVICE>& Lap, int J, unsigned int verbosity = 0, const double stoppingError = 0.0, const bool errIsRel = true, const bool enable_large_Faust = false);
+				GivensFGFT(Faust::MatDense<FPP,DEVICE>& Lap, int J, unsigned int verbosity = 0, const double stoppingError = 0.0, const bool errIsRel = true, const bool enable_large_Faust = false);
 				/** 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 ec920ad6609413c3cb2553d3204f2e5f8f0568a0..19e7a90fe999920b1ae00c8fdb4b7761270975c1 100644
--- a/src/algorithm/factorization/faust_GivensFGFT.hpp
+++ b/src/algorithm/factorization/faust_GivensFGFT.hpp
@@ -379,7 +379,7 @@ void GivensFGFT<FPP,DEVICE,FPP2>::update_err()
 }
 
 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 */, const double stoppingError /* default to 0.0 */, const bool errIsRel) :  GivensFGFTGen<FPP,DEVICE,FPP2>(Lap, J, verbosity, stoppingError, errIsRel), C(Lap.getNbRow(), Lap.getNbCol()), always_theta2(false)
+GivensFGFT<FPP,DEVICE,FPP2>::GivensFGFT(Faust::MatSparse<FPP,DEVICE>& Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError /* default to 0.0 */, const bool errIsRel, const bool enable_large_Faust/* deft to false */) :  GivensFGFTGen<FPP,DEVICE,FPP2>(Lap, J, verbosity, stoppingError, errIsRel, enable_large_Faust), C(Lap.getNbRow(), Lap.getNbCol()), always_theta2(false)
 {
 	//see parent ctor
 	C.setOnes();
@@ -387,7 +387,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 */, const double stoppingError, const bool errIsRel) : GivensFGFTGen<FPP,DEVICE,FPP2>(Lap, J, verbosity, stoppingError, errIsRel), C(Lap.getNbRow(), Lap.getNbCol()), always_theta2(false)
+GivensFGFT<FPP,DEVICE,FPP2>::GivensFGFT(Faust::MatDense<FPP,DEVICE>& Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError, const bool errIsRel, const bool enable_large_Faust/* deft to false */) : GivensFGFTGen<FPP,DEVICE,FPP2>(Lap, J, verbosity, stoppingError, errIsRel, enable_large_Faust), C(Lap.getNbRow(), Lap.getNbCol()), always_theta2(false)
 {
 	// see parent ctor
 	C.setOnes();
diff --git a/src/algorithm/factorization/faust_GivensFGFTComplex.h b/src/algorithm/factorization/faust_GivensFGFTComplex.h
index 3c8078a74397f1341048662e18ffbec2d4e6f17c..15f91b1ec72ce6752a6ed58ba8a89f5982180d37 100644
--- a/src/algorithm/factorization/faust_GivensFGFTComplex.h
+++ b/src/algorithm/factorization/faust_GivensFGFTComplex.h
@@ -53,8 +53,8 @@ namespace Faust {
 			public:
 
 				const static unsigned int ERROR_CALC_PERIOD = 100;
-				GivensFGFTComplex(Faust::MatSparse<FPP,DEVICE>& Lap, int J, unsigned int verbosity = 0, const double stoppingError = 0.0, const bool errIsRel = true);
-				GivensFGFTComplex(Faust::MatDense<FPP,DEVICE>& Lap, int J, unsigned int verbosity = 0, const double stoppingError = 0.0, const bool errIsRel = true);
+				GivensFGFTComplex(Faust::MatSparse<FPP,DEVICE>& Lap, int J, unsigned int verbosity = 0, const double stoppingError = 0.0, const bool errIsRel = true, const bool enable_large_Faust = false);
+				GivensFGFTComplex(Faust::MatDense<FPP,DEVICE>& Lap, int J, unsigned int verbosity = 0, const double stoppingError = 0.0, const bool errIsRel = true, const bool enable_large_Faust = false);
 
 				virtual void next_step();
 
diff --git a/src/algorithm/factorization/faust_GivensFGFTComplex.hpp b/src/algorithm/factorization/faust_GivensFGFTComplex.hpp
index b44cb3719df12a36b4f384a3fb101e3bb22656f1..fd513a9c28c5d2a20c8231b85413084247e9f9c7 100644
--- a/src/algorithm/factorization/faust_GivensFGFTComplex.hpp
+++ b/src/algorithm/factorization/faust_GivensFGFTComplex.hpp
@@ -445,14 +445,14 @@ void GivensFGFTComplex<FPP,DEVICE,FPP2>::update_err()
 }
 
 template<typename FPP, Device DEVICE, typename FPP2>
-GivensFGFTComplex<FPP,DEVICE,FPP2>::GivensFGFTComplex(Faust::MatSparse<FPP,DEVICE>& Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError /* default to 0.0 */, const bool errIsRel) : Faust::GivensFGFTGen<typename FPP::value_type, DEVICE, FPP2, FPP>(Lap, J, verbosity, stoppingError, errIsRel), C(Lap.getNbRow(), Lap.getNbCol())
+GivensFGFTComplex<FPP,DEVICE,FPP2>::GivensFGFTComplex(Faust::MatSparse<FPP,DEVICE>& Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError /* default to 0.0 */, const bool errIsRel, const bool enable_large_Faust) : Faust::GivensFGFTGen<typename FPP::value_type, DEVICE, FPP2, FPP>(Lap, J, verbosity, stoppingError, errIsRel, enable_large_Faust), C(Lap.getNbRow(), Lap.getNbCol())
 {
 	//see parent ctor
 	this->C.setZeros();
 }
 
 template<typename FPP, Device DEVICE, typename FPP2>
-GivensFGFTComplex<FPP,DEVICE,FPP2>::GivensFGFTComplex(Faust::MatDense<FPP,DEVICE>& Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError, const bool errIsRel) : Faust::GivensFGFTGen<typename FPP::value_type, DEVICE, FPP2, FPP>(Lap, J, verbosity, stoppingError, errIsRel), C(Lap.getNbRow(), Lap.getNbCol())
+GivensFGFTComplex<FPP,DEVICE,FPP2>::GivensFGFTComplex(Faust::MatDense<FPP,DEVICE>& Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError, const bool errIsRel, const bool enable_large_Faust) : Faust::GivensFGFTGen<typename FPP::value_type, DEVICE, FPP2, FPP>(Lap, J, verbosity, stoppingError, errIsRel, enable_large_Faust), C(Lap.getNbRow(), Lap.getNbCol())
 {
 	// see parent ctor
 	this->C.setZeros();
diff --git a/src/algorithm/factorization/faust_GivensFGFTGen.h b/src/algorithm/factorization/faust_GivensFGFTGen.h
index 7aa756da7fcf92cb0fc7b6b09cd88e3d81052a8c..903d3f8da0b409c1abd068f852c4fb45e24bcce0 100644
--- a/src/algorithm/factorization/faust_GivensFGFTGen.h
+++ b/src/algorithm/factorization/faust_GivensFGFTGen.h
@@ -107,7 +107,8 @@ namespace Faust {
 				double stoppingError;
 				/** \brief true if the stopping error is taken as relative error (absolute otherwise). */
 				bool errIsRel;
-
+				/** \brief (false to default) true to force the computation of the transform even if it doesn't worth it in term of complexity */
+				bool enable_large_Faust;
 
 			public:
 
@@ -117,10 +118,10 @@ namespace Faust {
 				 * \param J The number of iterations, Givens rotations factors.
 				 * TODO: complete argument list
 				 * */
-				GivensFGFTGen(Faust::MatGeneric<FPP4,DEVICE>* Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError, const bool errIsRel);
+				GivensFGFTGen(Faust::MatGeneric<FPP4,DEVICE>* Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError, const bool errIsRel, const bool enable_large_Faust = false);
 
-				GivensFGFTGen(Faust::MatSparse<FPP4, DEVICE> & Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError, const bool errIsRel);
-				GivensFGFTGen(Faust::MatDense<FPP4, DEVICE> & Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError, const bool errIsRel);
+				GivensFGFTGen(Faust::MatSparse<FPP4, DEVICE> & Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError, const bool errIsRel, const bool enable_large_Faust = false);
+				GivensFGFTGen(Faust::MatDense<FPP4, DEVICE> & Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError, const bool errIsRel, const bool enable_large_Faust = false);
 
 				/** Destructor */
 				virtual ~GivensFGFTGen() {delete[] q_candidates; delete L;};
diff --git a/src/algorithm/factorization/faust_GivensFGFTGen.hpp b/src/algorithm/factorization/faust_GivensFGFTGen.hpp
index f4550e2cf6b0ddbbb5f624cc6fa11340d989a699..f8f75935e187cedae31391cff1a49552242c6589 100644
--- a/src/algorithm/factorization/faust_GivensFGFTGen.hpp
+++ b/src/algorithm/factorization/faust_GivensFGFTGen.hpp
@@ -446,12 +446,18 @@ void GivensFGFTGen<FPP,DEVICE,FPP2,FPP4>::compute_facts()
 	is_D_ordered = false; // facts (re)computed then D must be reordered
 	ite = 0;
 	bool stopping = false;
-	while(J == 0 || ite < facts.size()) // when J == 0 the stopping criterion is the error against Lap
+	if(stopping = !enable_large_Faust && ! stoppingCritIsError && dim_size*dim_size <= J*4)
+	{
+		cerr << "WARNING: the eigtj algorithm stopped because the transform to be computed doesn't worth it according to its complexity (in space and time) relatively to the size of the matrix to decompose. Still, if you want to force the computation, please use the enable_large_Faust flag." << endl;
+		cerr << "nGivens: " << J << " matrix size: " << dim_size << endl;
+		facts.resize(0);
+	}
+	while(! stopping && (J == 0 || ite < facts.size())) // when J == 0 the stopping criterion is the error against Lap
 	{
 		next_step();
 		ite++;
 		if(stopping = (ite > 1 && stoppingCritIsError && errs.size() > 2 && errs[ite-1]-errs[ite-2] > FLT_EPSILON))
-			if(verbosity>0) cerr << "warning: the algorithm stopped because the last error is greater than the previous one." << endl;
+			/*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)
 		{
 			facts.resize(ite);
@@ -529,8 +535,8 @@ void GivensFGFTGen<FPP,DEVICE,FPP2,FPP4>::compute_facts()
 //
 
 template<typename FPP, Device DEVICE, typename FPP2, typename FPP4>
-GivensFGFTGen<FPP,DEVICE,FPP2,FPP4>::GivensFGFTGen(MatGeneric<FPP4,DEVICE>* Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError, const bool errIsRel) :
-Lap(*Lap),  D(Lap->getNbRow()), errs(0), coord_choices(0), q_candidates(new int[Lap->getNbRow()]), is_D_ordered(false), verbosity(verbosity), stoppingCritIsError(stoppingError != 0.0), stoppingError(stoppingError), errIsRel(errIsRel), Lap_squared_fro_norm(0), facts(J>0?J:1), last_fact_permuted(false), J(J), dim_size(Lap->getNbRow())
+GivensFGFTGen<FPP,DEVICE,FPP2,FPP4>::GivensFGFTGen(MatGeneric<FPP4,DEVICE>* Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError, const bool errIsRel, const bool enable_large_Faust /* deft to false */) :
+Lap(*Lap),  D(Lap->getNbRow()), errs(0), coord_choices(0), q_candidates(new int[Lap->getNbRow()]), is_D_ordered(false), verbosity(verbosity), stoppingCritIsError(stoppingError != 0.0), stoppingError(stoppingError), errIsRel(errIsRel), Lap_squared_fro_norm(0), facts(J>0?(J*4<Lap->getNbRow()*Lap->getNbRow()||enable_large_Faust?J:0):1 /* don't allocate if the complexity doesn't worth it and enable_large_Faust is false*/), last_fact_permuted(false), J(J), dim_size(Lap->getNbRow()), enable_large_Faust(enable_large_Faust)
 {
 	if(Lap->getNbCol() != Lap->getNbRow())
 		handleError("Faust::GivensFGFTComplex", "Laplacian must be a square matrix.");
@@ -549,13 +555,13 @@ Lap(*Lap),  D(Lap->getNbRow()), errs(0), coord_choices(0), q_candidates(new int[
 }
 
 template<typename FPP, Device DEVICE, typename FPP2, typename FPP4>
-GivensFGFTGen<FPP,DEVICE,FPP2,FPP4>::GivensFGFTGen(Faust::MatSparse<FPP4, DEVICE> & Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError, const bool errIsRel) :  GivensFGFTGen<FPP,DEVICE,FPP2,FPP4>(&Lap, J, verbosity, stoppingError, errIsRel)
+GivensFGFTGen<FPP,DEVICE,FPP2,FPP4>::GivensFGFTGen(Faust::MatSparse<FPP4, DEVICE> & Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError, const bool errIsRel, const bool enable_large_Faust/* deft to false */) :  GivensFGFTGen<FPP,DEVICE,FPP2,FPP4>(&Lap, J, verbosity, stoppingError, errIsRel, enable_large_Faust)
 {
 	L = new MatSparse<FPP4,DEVICE>(Lap);
 }
 
 template<typename FPP, Device DEVICE, typename FPP2, typename FPP4>
-GivensFGFTGen<FPP,DEVICE,FPP2,FPP4>::GivensFGFTGen(Faust::MatDense<FPP4, DEVICE> & Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError, const bool errIsRel) : GivensFGFTGen<FPP,DEVICE,FPP2,FPP4>(&Lap, J, verbosity, stoppingError, errIsRel)
+GivensFGFTGen<FPP,DEVICE,FPP2,FPP4>::GivensFGFTGen(Faust::MatDense<FPP4, DEVICE> & Lap, int J, unsigned int verbosity /* deft val == 0 */, const double stoppingError, const bool errIsRel, const bool enable_large_Faust/* deft to false */) : GivensFGFTGen<FPP,DEVICE,FPP2,FPP4>(&Lap, J, verbosity, stoppingError, errIsRel,  enable_large_Faust)
 {
 	L = new MatDense<FPP4,DEVICE>(Lap);
 }
@@ -695,6 +701,8 @@ template<typename FPP, Device DEVICE, typename FPP2, typename FPP4>
 Faust::Transform<FPP4,DEVICE> GivensFGFTGen<FPP,DEVICE,FPP2,FPP4>::get_transform(int ord)
 {
 	//TODO: an optimization is possible by changing type of facts to vector<MatGeneric*> it would avoid copying facts into Transform and rather use them directly. It will need a destructor that deletes them eventually if they weren't transfered to a Transform object before.
+	if(facts.size() == 0)
+		throw out_of_range("The transform is empty. The algorithm stopped before computing any factor.");
 	MatSparse<FPP4,DEVICE> & last_fact = *(facts.end()-1);
 	MatSparse<FPP4,DEVICE> P(last_fact.getNbCol(), last_fact.getNbCol()); //last_fact permutation matrix
 	// (to reorder eigenvector transform according to ordered D)
diff --git a/src/algorithm/factorization/faust_GivensFGFTParallel.h b/src/algorithm/factorization/faust_GivensFGFTParallel.h
index 8e51c7f426973ebef5dbbb52cca5ed309b88a1ac..cc65eceb1aade7d211b8f8d04c4189fa7214a52c 100644
--- a/src/algorithm/factorization/faust_GivensFGFTParallel.h
+++ b/src/algorithm/factorization/faust_GivensFGFTParallel.h
@@ -84,8 +84,8 @@ namespace Faust {
 		 * \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, 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);
+		GivensFGFTParallel(Faust::MatDense<FPP,DEVICE>& Lap, int J, int t, unsigned int verbosity = 0, const double stoppingCritIsError = 0.0,  const bool errIsRel = true, const bool enable_large_Faust = false);
+		GivensFGFTParallel(Faust::MatSparse<FPP,DEVICE>& Lap, int J, int t, unsigned int verbosity = 0, const double stoppingCritIsError = 0.0, const bool errIsRel = true, const bool enable_large_Faust = false);
 
 	};
 
diff --git a/src/algorithm/factorization/faust_GivensFGFTParallel.hpp b/src/algorithm/factorization/faust_GivensFGFTParallel.hpp
index 55b7454615c54fd6998535d496999ead7c4b43f8..dcf12a2d91a0452ca379baed21f7f955a20ab0c2 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, const double stoppingError, const bool errIsRel) : GivensFGFT<FPP,DEVICE,FPP2>(Lap, J, verbosity, stoppingError, errIsRel), 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, const bool enable_large_Faust) : GivensFGFT<FPP,DEVICE,FPP2>(Lap, J, verbosity, stoppingError, errIsRel, enable_large_Faust), t(t), fact_nrots(0)
 {
 	if(J > 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, const double stoppingError, const bool errIsRel) : GivensFGFT<FPP,DEVICE,FPP2>(Lap, J, verbosity, stoppingError, errIsRel), 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, const bool enable_large_Faust) : GivensFGFT<FPP,DEVICE,FPP2>(Lap, J, verbosity, stoppingError, errIsRel, enable_large_Faust), t(t), fact_nrots(0)
 {
 	if(J > 0) this->facts.resize(round(J/(float)t));
 	this->always_theta2 = true;
diff --git a/src/algorithm/factorization/faust_GivensFGFTParallelComplex.h b/src/algorithm/factorization/faust_GivensFGFTParallelComplex.h
index 67c121d0234be1cab2e9707e9ac4d178a903e9b3..4bb116b56e2b17374673e9c48d691682fcde25ef 100644
--- a/src/algorithm/factorization/faust_GivensFGFTParallelComplex.h
+++ b/src/algorithm/factorization/faust_GivensFGFTParallelComplex.h
@@ -86,8 +86,8 @@ namespace Faust {
 		 * \param stoppingError defines a stopping criterion based on error (absolute relative error).
 		 *
 		 */
-		GivensFGFTParallelComplex(Faust::MatDense<FPP,DEVICE>& Lap, int J, int t, unsigned int verbosity = 0, const double stoppingCritIsError = 0.0,  const bool errIsRel = true);
-		GivensFGFTParallelComplex(Faust::MatSparse<FPP,DEVICE>& Lap, int J, int t, unsigned int verbosity = 0, const double stoppingCritIsError = 0.0, const bool errIsRel = true);
+		GivensFGFTParallelComplex(Faust::MatDense<FPP,DEVICE>& Lap, int J, int t, unsigned int verbosity = 0, const double stoppingCritIsError = 0.0,  const bool errIsRel = true, const bool enable_large_Faust = false);
+		GivensFGFTParallelComplex(Faust::MatSparse<FPP,DEVICE>& Lap, int J, int t, unsigned int verbosity = 0, const double stoppingCritIsError = 0.0, const bool errIsRel = true, const bool enable_large_Faust = false);
 
 	};
 
diff --git a/src/algorithm/factorization/faust_GivensFGFTParallelComplex.hpp b/src/algorithm/factorization/faust_GivensFGFTParallelComplex.hpp
index 7c71d4cf39c2a5cd5080795668d07ce6c97f4b63..a460341ce0daadea98f4e2c00d26296877aa2849 100644
--- a/src/algorithm/factorization/faust_GivensFGFTParallelComplex.hpp
+++ b/src/algorithm/factorization/faust_GivensFGFTParallelComplex.hpp
@@ -6,16 +6,16 @@ using namespace Faust;
 #endif
 
 template<typename FPP, Device DEVICE, typename FPP2>
-GivensFGFTParallelComplex<FPP,DEVICE,FPP2>::GivensFGFTParallelComplex(Faust::MatDense<FPP,DEVICE>& Lap, int J, int t, unsigned int verbosity, const double stoppingError, const bool errIsRel) : GivensFGFTComplex<FPP,DEVICE,FPP2>(Lap, J, verbosity, stoppingError, errIsRel), t(t), fact_nrots(0)
+GivensFGFTParallelComplex<FPP,DEVICE,FPP2>::GivensFGFTParallelComplex(Faust::MatDense<FPP,DEVICE>& Lap, int J, int t, unsigned int verbosity, const double stoppingError, const bool errIsRel,  const bool enable_large_Faust) : GivensFGFTComplex<FPP,DEVICE,FPP2>(Lap, J, verbosity, stoppingError, errIsRel, enable_large_Faust), t(t), fact_nrots(0)
 {
-	this->facts.resize(round(J/(float)t));
+	if(J > 0) this->facts.resize(round(J/(float)t));
 	this->coord_choices.resize(0);
 }
 
 template<typename FPP, Device DEVICE, typename FPP2>
-GivensFGFTParallelComplex<FPP,DEVICE,FPP2>::GivensFGFTParallelComplex(Faust::MatSparse<FPP,DEVICE>& Lap, int J, int t, unsigned int verbosity, const double stoppingError, const bool errIsRel) : GivensFGFTComplex<FPP,DEVICE,FPP2>(Lap, J, verbosity, stoppingError, errIsRel), t(t), fact_nrots(0)
+GivensFGFTParallelComplex<FPP,DEVICE,FPP2>::GivensFGFTParallelComplex(Faust::MatSparse<FPP,DEVICE>& Lap, int J, int t, unsigned int verbosity, const double stoppingError, const bool errIsRel, const bool enable_large_Faust) : GivensFGFTComplex<FPP,DEVICE,FPP2>(Lap, J, verbosity, stoppingError, errIsRel, enable_large_Faust), t(t), fact_nrots(0)
 {
-	this->facts.resize(round(J/(float)t));
+	if(J > 0) this->facts.resize(round(J/(float)t));
 	this->coord_choices.resize(0);
 }
 
diff --git a/wrapper/matlab/+matfaust/+fact/eigtj.m b/wrapper/matlab/+matfaust/+fact/eigtj.m
index 655708662aebdc1411817831859f83dc01c381da..0b7a957035e102a1ce4d26aa59f8e66457a5fce2 100644
--- a/wrapper/matlab/+matfaust/+fact/eigtj.m
+++ b/wrapper/matlab/+matfaust/+fact/eigtj.m
@@ -75,7 +75,7 @@
 %==========================================================================================
 function [V,D] = eigtj(M, varargin)
 	import matfaust.Faust
-	nGivens_per_fac = 1; % default value
+	nGivens_per_fac = floor(size(M,1)/2); % default value
 	verbosity = 0; % default value
 %	if(~ ismatrix(M) || ~ isreal(M))
 %		error('M must be a real matrix.')
@@ -88,11 +88,18 @@ function [V,D] = eigtj(M, varargin)
 	tol = 0;
 	relerr = true;
 	verbosity = 0;
+	enable_large_Faust = false;
 	argc = length(varargin);
 	order = 1; % ascending order
 	if(argc > 0)
 		for i=1:argc
 			switch(varargin{i})
+				case 'enable_large_Faust'
+					if(argc == i || ~ islogical(varargin{i+1}))
+						error('enable_large_Faust keyword argument is not followed by a logical')
+					else
+						enable_large_Faust = varargin{i+1};
+					end
 				case 'nGivens'
 					if(argc == i || ~ isnumeric(varargin{i+1}) || varargin{i+1}-floor(varargin{i+1}) > 0 || varargin{i+1} <= 0 )
 						error('nGivens keyword arg. is not followed by an integer')
@@ -149,7 +156,7 @@ function [V,D] = eigtj(M, varargin)
 	if(nGivens > 0)
 		nGivens_per_fac = min(nGivens_per_fac, nGivens);
 	end
-	[core_obj, D] = mexfgftgivensReal(M, nGivens, nGivens_per_fac, verbosity, tol, relerr, order);
+	[core_obj, D] = mexfgftgivensReal(M, nGivens, nGivens_per_fac, verbosity, tol, relerr, order, enable_large_Faust);
 	D = sparse(diag(real(D)));
 	V = Faust(core_obj, isreal(M));
 end
diff --git a/wrapper/matlab/src/mexfgftgivens.cpp.in b/wrapper/matlab/src/mexfgftgivens.cpp.in
index d2ac822e27a21e5898bb87a504d4348d443a8099..898dab5107ae4c5fc7eb8d26f8db81ecfec78629 100644
--- a/wrapper/matlab/src/mexfgftgivens.cpp.in
+++ b/wrapper/matlab/src/mexfgftgivens.cpp.in
@@ -55,9 +55,9 @@ typedef @FACT_FPP@ FPP2;
 
 using namespace Faust;
 
-void fgft_givens(const mxArray* matlab_matrix, int J, int t, double tol, unsigned int verbosity, bool relErr, int order, mxArray **plhs);
+void fgft_givens(const mxArray* matlab_matrix, int J, int t, double tol, unsigned int verbosity, bool rel_err,  int order, const bool enable_large_Faust, mxArray **plhs);
 
-void fgft_givens_cplx(const mxArray* matlab_matrix, int J, int t, double tol, unsigned int verbosity, bool relErr, int order, mxArray **plhs);
+void fgft_givens_cplx(const mxArray* matlab_matrix, int J, int t, double tol, unsigned int verbosity, bool rel_err, int order, const bool enable_large_Faust, mxArray **plhs);
 
 void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
 {
@@ -66,11 +66,12 @@ void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
 	int t = 1; // default value for non-parallel Givens FGFT
 	unsigned int verbosity = 0; //default verbosity (no info. displayed)
 	double tol = 0;
-	bool relErr = true;
+	bool rel_err = true;
+	bool enable_large_Faust = false;
 	int order;
 
-	if(nrhs < 2 || nrhs > 7)
-		mexErrMsgTxt("Bad Number of inputs arguments");
+	if(nrhs < 2 || nrhs > 8)
+		mexErrMsgTxt("Bad number of input arguments");
 
 	J = (int) mxGetScalar(prhs[1]);
 	if(nrhs >= 3)
@@ -80,22 +81,24 @@ void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
 	if(nrhs >= 5)
 		tol = (double) mxGetScalar(prhs[4]);
 	if(nrhs >= 6)
-		relErr = (bool) mxGetScalar(prhs[5]);
+		rel_err = (bool) mxGetScalar(prhs[5]);
 	if(nrhs >= 7)
 		order = (int) mxGetScalar(prhs[6]); //eigenvalues order
+	if(nrhs >= 8)
+		enable_large_Faust = (bool) mxGetScalar(prhs[7]);
 
 	tol *= tol; // C++ backend works with squared norm error
 
 	const mxArray* matlab_matrix = prhs[0]; // Laplacian
 
 	if(mxIsComplex(matlab_matrix))
-		fgft_givens_cplx(matlab_matrix, J, t, tol, verbosity, relErr, order, plhs);
+		fgft_givens_cplx(matlab_matrix, J, t, tol, verbosity, rel_err, order, enable_large_Faust, plhs);
 	else
-		fgft_givens(matlab_matrix, J, t, tol, verbosity, relErr, order, plhs);
+		fgft_givens(matlab_matrix, J, t, tol, verbosity, rel_err,  order, enable_large_Faust, plhs);
 
 }
 
-void fgft_givens(const mxArray* matlab_matrix, int J, int t, double tol, unsigned int verbosity, bool relErr, int order, mxArray **plhs)
+void fgft_givens(const mxArray* matlab_matrix, int J, int t, double tol, unsigned int verbosity, bool rel_err, int order, const bool enable_large_Faust, mxArray **plhs)
 {
 	// initialization of the matrix that will be factorized
 	Faust::MatGeneric<SCALAR,Cpu>* Lap;
@@ -113,18 +116,18 @@ void fgft_givens(const mxArray* matlab_matrix, int J, int t, double tol, unsigne
 			Lap = &sLap;
 			//TODO: blas handles to pass later ? (if you want to do the calculation with blas instead of eigen)
 			if(t <= 1)
-				algo = new GivensFGFT<SCALAR,Cpu,FPP2>(sLap, J, verbosity, tol, relErr);
+				algo = new GivensFGFT<SCALAR,Cpu,FPP2>(sLap, J, verbosity, tol, rel_err, enable_large_Faust);
 			else
-				algo = new GivensFGFTParallel<SCALAR,Cpu,FPP2>(sLap, J, t, verbosity, tol, relErr);
+				algo = new GivensFGFTParallel<SCALAR,Cpu,FPP2>(sLap, J, t, verbosity, tol, rel_err, enable_large_Faust);
 
 		}else
 		{
 			mxArray2FaustMat(matlab_matrix,dLap);
 			Lap = &dLap;
 			if(t <= 1)
-				algo = new GivensFGFT<SCALAR,Cpu,FPP2>(dLap, J, verbosity, tol, relErr);
+				algo = new GivensFGFT<SCALAR,Cpu,FPP2>(dLap, J, verbosity, tol, rel_err, enable_large_Faust);
 			else
-				algo = new GivensFGFTParallel<SCALAR,Cpu,FPP2>(dLap, J, t, verbosity, tol, relErr);
+				algo = new GivensFGFTParallel<SCALAR,Cpu,FPP2>(dLap, J, t, verbosity, tol, rel_err, enable_large_Faust);
 		}
 
 
@@ -142,16 +145,18 @@ void fgft_givens(const mxArray* matlab_matrix, int J, int t, double tol, unsigne
 
 		plhs[0] = convertPtr2Mat<Faust::TransformHelper<SCALAR, Cpu>>(th);
 
-		delete algo;
 	}
 	catch (const std::exception& e)
 	{
 		mexErrMsgTxt(e.what());
 	}
 
+	if(algo != nullptr)
+		delete algo;
+
 }
 
-void fgft_givens_cplx(const mxArray* matlab_matrix, int J, int t, double tol, unsigned int verbosity, bool relErr, int order, mxArray **plhs)
+void fgft_givens_cplx(const mxArray* matlab_matrix, int J, int t, double tol, unsigned int verbosity, bool rel_err, int order, const bool enable_large_Faust, mxArray **plhs)
 {
 	// initialization of the matrix that will be factorized
 	Faust::MatGeneric<complex<SCALAR>,Cpu>* Lap;
@@ -160,7 +165,7 @@ void fgft_givens_cplx(const mxArray* matlab_matrix, int J, int t, double tol, un
 	Faust::BlasHandle<Cpu> blas_handle;
 	Faust::SpBlasHandle<Cpu> spblas_handle;
 
-	Faust::GivensFGFTComplex<complex<SCALAR>, Cpu, FPP2>* algo;
+	Faust::GivensFGFTComplex<complex<SCALAR>, Cpu, FPP2>* algo = nullptr;
 
 	try{
 		if (mxIsSparse(matlab_matrix))
@@ -169,22 +174,20 @@ void fgft_givens_cplx(const mxArray* matlab_matrix, int J, int t, double tol, un
 			Lap = &sLap;
 			//TODO: blas handles to pass later ? (if you want to do the calculation with blas instead of eigen)
 			if(t <= 1)
-				algo = new GivensFGFTComplex<complex<SCALAR>,Cpu,FPP2>(sLap, J, verbosity, tol, relErr);
+				algo = new GivensFGFTComplex<complex<SCALAR>,Cpu,FPP2>(sLap, J, verbosity, tol, rel_err, enable_large_Faust);
 			else
-				algo = new GivensFGFTParallelComplex<complex<SCALAR>,Cpu,FPP2>(sLap, J, t, verbosity, tol, relErr);
+				algo = new GivensFGFTParallelComplex<complex<SCALAR>,Cpu,FPP2>(sLap, J, t, verbosity, tol, rel_err, enable_large_Faust);
 
 		}else
 		{
 			mxArray2FaustMat(matlab_matrix,dLap);
 			Lap = &dLap;
 			if(t <= 1)
-				algo = new GivensFGFTComplex<complex<SCALAR>,Cpu,FPP2>(dLap, J, verbosity, tol, relErr);
+				algo = new GivensFGFTComplex<complex<SCALAR>,Cpu,FPP2>(dLap, J, verbosity, tol, rel_err, enable_large_Faust);
 			else
-				algo = new GivensFGFTParallelComplex<complex<SCALAR>,Cpu,FPP2>(dLap, J, t, verbosity, tol, relErr);
+				algo = new GivensFGFTParallelComplex<complex<SCALAR>,Cpu,FPP2>(dLap, J, t, verbosity, tol, rel_err, enable_large_Faust);
 		}
 
-
-
 		algo->compute_facts();
 
 		Faust::Vect<SCALAR,Cpu> D(Lap->getNbRow());
@@ -198,13 +201,15 @@ void fgft_givens_cplx(const mxArray* matlab_matrix, int J, int t, double tol, un
 
 		plhs[0] = convertPtr2Mat<Faust::TransformHelper<complex<SCALAR>, Cpu>>(th);
 
-		delete algo;
 	}
 	catch (const std::exception& e)
 	{
+		mexErrMsgTxt(e.what());
 
 	}
 
+	if(algo != nullptr)
+		delete algo;
 }
 
 
diff --git a/wrapper/python/pyfaust/fact.py b/wrapper/python/pyfaust/fact.py
index 2d2f23e5c4f6ca1fb710acd10fd4af270fad08eb..72602d333efaeb34e12268883066a9443089b39e 100644
--- a/wrapper/python/pyfaust/fact.py
+++ b/wrapper/python/pyfaust/fact.py
@@ -170,7 +170,7 @@ def svdtj(M, nGivens=None, tol=0, order='ascend', relerr=True,
     return U, S, V
 
 def eigtj(M, nGivens=None, tol=0, order='ascend', relerr=True,
-          nGivens_per_fac=None, verbosity=0):
+          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 left eigenvectors (as the columns of the Faust object V).
 
@@ -255,12 +255,13 @@ def eigtj(M, nGivens=None, tol=0, order='ascend', relerr=True,
         svdtj
     """
     D, core_obj = _FaustCorePy.FaustFact.eigtj(M, nGivens, tol, relerr,
-            nGivens_per_fac, verbosity, order)
+                                               nGivens_per_fac, verbosity, order,
+                                               enable_large_Faust)
     return D, Faust(core_obj=core_obj)
 
 # experimental block start
 def fgft_givens(Lap, nGivens=None, tol=0, order='ascend', relerr=True,
-          nGivens_per_fac=None, verbosity=0):
+          nGivens_per_fac=None, verbosity=0, enable_large_Faust=False):
     """
     Computes the FGFT of the Laplacian matrix Lap (using fact.eigtj).
 
@@ -292,7 +293,8 @@ def fgft_givens(Lap, nGivens=None, tol=0, order='ascend', relerr=True,
     See also:
         eigtj, fgft_palm
     """
-    return eigtj(Lap, nGivens, tol, order, relerr, nGivens_per_fac, verbosity)
+    return eigtj(Lap, nGivens, tol, order, relerr, nGivens_per_fac, verbosity,
+                enable_large_Faust)
 # experimental block end
 
 def _check_fact_mat(funcname, M):
diff --git a/wrapper/python/src/FaustCoreCy.pxd b/wrapper/python/src/FaustCoreCy.pxd
index 92dcc73680ff1b9d5c331a5c724dd56a4f50a2de..8e8ca85bf5c174ae3192772031b009bbc01e7852 100644
--- a/wrapper/python/src/FaustCoreCy.pxd
+++ b/wrapper/python/src/FaustCoreCy.pxd
@@ -188,20 +188,22 @@ cdef extern from "FaustFact.h":
 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 num_cols, unsigned int J,
+                                                       unsigned int t, FPP* D, unsigned int verbosity,
                                                        const double stoppingError,
                                                        const bool errIsRel,
-                                                       const int order)
+                                                       const int order,
+                                                       const bool enable_large_Faust)
 
     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,
+                                                              int* id_col, int nnz, unsigned int num_rows,
+                                                              unsigned int num_cols, unsigned int J,
                                                               unsigned int t, FPP* D,
                                                               unsigned int verbosity,
                                                               const double stoppingError,
                                                               const bool errIsRel,
-                                                              const int order)
+                                                              const int order,
+                                                              const bool enable_large_Faust)
 
     cdef FaustCoreCpp[FPP]* fact_givens_fgft_cplx[FPP,FPP2](const
                                                             FPP* Lap, unsigned int num_rows,
@@ -210,17 +212,19 @@ cdef extern from "FaustFactGivensFGFT.h":
                                                             FPP2* D, unsigned int verbosity,
                                                             const double stoppingError,
                                                             const bool errIsRel,
-                                                            const int order)
+                                                            const int order,
+                                                            const bool enable_large_Faust)
 
     cdef FaustCoreCpp[FPP]* fact_givens_fgft_sparse_cplx[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,
-                                 FPP2* D,
-                                 unsigned int verbosity,
-                                 const double stoppingError,
-                                 const bool errIsRel,
-                                 const int order)
+                                                                    int* id_col, int nnz, unsigned int num_rows,
+                                                                    unsigned int num_cols, unsigned int J,
+                                                                    unsigned int t,
+                                                                    FPP2* D,
+                                                                    unsigned int verbosity,
+                                                                    const double stoppingError,
+                                                                    const bool errIsRel,
+                                                                    const int order,
+                                                                    const bool enable_large_Faust)
 
     cdef void svdtj[FPP, FPP2](FaustCoreCpp[FPP]** U, FaustCoreCpp[FPP] **V, FPP* S,
                          const FPP* M_data, unsigned int num_rows, unsigned int
diff --git a/wrapper/python/src/FaustFactGivensFGFT.h b/wrapper/python/src/FaustFactGivensFGFT.h
index 701f4a87b35bdf4c21025b55c3b10eede713efa9..d8174ae51faa0c0c54d9352db2f951eb716b2cbc 100644
--- a/wrapper/python/src/FaustFactGivensFGFT.h
+++ b/wrapper/python/src/FaustFactGivensFGFT.h
@@ -7,22 +7,22 @@ FaustCoreCpp<FPP>* fact_givens_fgft(const FPP* Lap, unsigned int num_rows, unsig
 
 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 double stoppingError = 0.0, const bool errIsRel = true, const int order = 1 /* ascendant order for eigenvalues */);
+        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, const int order = 1 /* ascendant order for eigenvalues */, const bool enable_large_Faust=false);
 
 
 template<typename FPP, typename FPP2 = float>
-FaustCoreCpp<FPP>* fact_givens_fgft_generic(GivensFGFT<FPP, Cpu, FPP2>* algo, FPP* D, const int order);
+FaustCoreCpp<FPP>* fact_givens_fgft_generic(GivensFGFT<FPP, Cpu, FPP2>* algo, FPP* D, const int order, const bool enable_large_Faust=false);
 
 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 double stoppingError = 0.0, const bool errIsRel = true, const int order = 1 /* ascendant order for eigenvalues */);
+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 double stoppingError = 0.0, const bool errIsRel = true, const int order = 1 /* ascendant order for eigenvalues */, const bool enable_large_Faust=false);
 
 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 double stoppingError = 0.0, const bool errIsRel = true, const int order = 1 /* ascendant order for eigenvalues */);
+        unsigned int J, unsigned int t, FPP2* D, unsigned int verbosity = 0, const double stoppingError = 0.0, const bool errIsRel = true, const int order = 1 /* ascendant order for eigenvalues */, const bool enable_large_Faust=false);
 
 
 template<typename FPP, typename FPP2 = float>
-FaustCoreCpp<FPP>* fact_givens_fgft_generic_cplx(GivensFGFTComplex<FPP, Cpu, FPP2>* algo, FPP2* D, const int order);
+FaustCoreCpp<FPP>* fact_givens_fgft_generic_cplx(GivensFGFTComplex<FPP, Cpu, FPP2>* algo, FPP2* D, const int order, const bool enable_large_Faust=false);
 
 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 J, unsigned int t, unsigned int verbosity = 0, const double stoppingError = 0.0, const bool errIsRel = true /* end of input parameters*/);
diff --git a/wrapper/python/src/FaustFactGivensFGFT.hpp b/wrapper/python/src/FaustFactGivensFGFT.hpp
index 10dc2eb8f5634c97a6cb916f84efe1578ccab7f4..8a358c680e654075d3e33f363c3ee3410fd3c680 100644
--- a/wrapper/python/src/FaustFactGivensFGFT.hpp
+++ b/wrapper/python/src/FaustFactGivensFGFT.hpp
@@ -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 double stoppingError, const bool errIsRel, const int order)
+        unsigned int J, unsigned int t /* end of input parameters*/, FPP* D, unsigned int verbosity, const double stoppingError, const bool errIsRel, const int order, const bool enable_large_Faust)
 {
       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);
+          algo = new GivensFGFT<FPP, Cpu, FPP2>(mat_Lap, (int)J, verbosity, stoppingError, errIsRel, enable_large_Faust);
       }
       else
       {
-          algo = new GivensFGFTParallel<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity, stoppingError, errIsRel);
+          algo = new GivensFGFTParallel<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity, stoppingError, errIsRel, enable_large_Faust);
       }
       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 double stoppingError, const bool errIsRel, const int order)
+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, const int order, const bool enable_large_Faust)
 {
     //TODO: optimization possible here by avoiding Lap copy in MatDense (by
     //just using the data in Lap as underlying pointer of MatDense)
@@ -33,53 +33,60 @@ 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);
+        algo = new GivensFGFT<FPP, Cpu, FPP2>(mat_Lap, (int)J, verbosity, stoppingError, errIsRel, enable_large_Faust);
     }
     else
     {
-        algo = new GivensFGFTParallel<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity, stoppingError, errIsRel);
+        algo = new GivensFGFTParallel<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity, stoppingError, errIsRel, enable_large_Faust);
     }
     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)
+FaustCoreCpp<FPP>* fact_givens_fgft_generic(GivensFGFT<FPP, Cpu, FPP2>* algo, FPP* D, const int order, const bool enable_large_Faust)
 {
 
+    FaustCoreCpp<FPP>* fc = nullptr;
 
     algo->compute_facts();
 
-    Faust::Transform<FPP,Cpu> trans = std::move(algo->get_transform(order));
-    TransformHelper<FPP,Cpu> *th = new TransformHelper<FPP,Cpu>(trans, true); // true is for moving and not copying the Transform object into TransformHelper (optimization possible cause we know the original object won't be used later)
+    try {
+        Faust::Transform<FPP,Cpu> trans = std::move(algo->get_transform(order));
+        TransformHelper<FPP,Cpu> *th = new TransformHelper<FPP,Cpu>(trans, true); // true is for moving and not copying the Transform object into TransformHelper (optimization possible cause we know the original object won't be used later)
 
 
-    //copy ordered diagonal in buffer D (allocated from the outside)
-    algo->get_D(D, order);
+        //copy ordered diagonal in buffer D (allocated from the outside)
+        algo->get_D(D, order);
 
-    delete algo;
+        return new FaustCoreCpp<FPP>(th);
+    }
+    catch(out_of_range e)
+    {
 
-    return new FaustCoreCpp<FPP>(th);
+    }
+    delete algo;
+    return fc;
 }
 
 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 double stoppingError, const bool errIsRel, const int order)
+        unsigned int J, unsigned int t /* end of input parameters*/, FPP2* D, unsigned int verbosity, const double stoppingError, const bool errIsRel, const int order, const bool enable_large_Faust)
 {
       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);
+          algo = new GivensFGFTComplex<FPP, Cpu, FPP2>(mat_Lap, (int)J, verbosity, stoppingError, errIsRel, enable_large_Faust);
       }
       else
       {
-          algo = new GivensFGFTParallelComplex<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity, stoppingError, errIsRel);
+          algo = new GivensFGFTParallelComplex<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity, stoppingError, errIsRel, enable_large_Faust);
       }
       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 double stoppingError, const bool errIsRel, const int order)
+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 double stoppingError, const bool errIsRel, const int order, const bool enable_large_Faust)
 {
     //TODO: optimization possible here by avoiding Lap copy in MatDense (by
     //just using the data in Lap as underlying pointer of MatDense)
@@ -88,32 +95,38 @@ 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);
+        algo = new GivensFGFTComplex<FPP, Cpu, FPP2>(mat_Lap, (int)J, verbosity, stoppingError, errIsRel, enable_large_Faust);
     }
     else
     {
-        algo = new GivensFGFTParallelComplex<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity, stoppingError, errIsRel);
+        algo = new GivensFGFTParallelComplex<FPP, Cpu, FPP2>(mat_Lap, (int)J, (int) t, verbosity, stoppingError, errIsRel, enable_large_Faust);
     }
     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)
+FaustCoreCpp<FPP>* fact_givens_fgft_generic_cplx(GivensFGFTComplex<FPP, Cpu, FPP2>* algo, FPP2* D, const int order, const bool enable_large_Faust)
 {
 
-
+    FaustCoreCpp<FPP>* fc = nullptr;
     algo->compute_facts();
 
-    Faust::Transform<FPP,Cpu> trans = std::move(algo->get_transform(order));
-    TransformHelper<FPP,Cpu> *th = new TransformHelper<FPP,Cpu>(trans, true); // true is for moving and not copying the Transform object into TransformHelper (optimization possible cause we know the original object won't be used later)
-
+    try
+    {
+        Faust::Transform<FPP,Cpu> trans = std::move(algo->get_transform(order));
+        TransformHelper<FPP,Cpu> *th = new TransformHelper<FPP,Cpu>(trans, true); // true is for moving and not copying the Transform object into TransformHelper (optimization possible cause we know the original object won't be used later)
 
-    //copy ordered diagonal in buffer D (allocated from the outside)
-    algo->get_D(D, order);
 
+        //copy ordered diagonal in buffer D (allocated from the outside)
+        algo->get_D(D, order);
+        fc =new FaustCoreCpp<FPP>(th);
+    }
+    catch(out_of_range e)
+    {
+    }
     delete algo;
 
-    return new FaustCoreCpp<FPP>(th);
+    return fc;
 }
 
 template<typename FPP, typename FPP2>
diff --git a/wrapper/python/src/_FaustCorePy.pyx b/wrapper/python/src/_FaustCorePy.pyx
index 4e4e3667d78999cf1905234925b0ff9e5567f040..41112d894f3633be17073306d39497968c310a66 100644
--- a/wrapper/python/src/_FaustCorePy.pyx
+++ b/wrapper/python/src/_FaustCorePy.pyx
@@ -1435,7 +1435,7 @@ cdef class FaustFact:
 
     @staticmethod
     def fact_givens_fgft_sparse(Lap, J, t, verbosity=0, stoppingError = 0.0,
-                                errIsRel=True, order='ascend'):
+                                errIsRel=True, order='ascend', enable_large_Faust=False):
         from scipy.sparse import spdiags
         cdef double [:] data1d #only for csr mat factor
         cdef int [:] indices # only for csr mat
@@ -1467,16 +1467,23 @@ cdef class FaustFact:
                                                                    verbosity,
                                                                    stoppingError,
                                                                    errIsRel,
-                                                                   int(order))
+                                                                   int(order),
+                                                                   enable_large_Faust)
 
         core._isReal = True
         #D_spdiag = spdiags(D, [0], Lap.shape[0], Lap.shape[0])
         #return core, D_spdiag
+        if(core._isReal and core.core_faust_dbl == NULL or \
+           not core._isReal and core.core_faust_cplx == NULL):
+            raise Exception("Empty transform (nGivens is too big ? Set"
+                            " enable_large_Faust to True to force the computation).")
+
+
         return core, D
 
     @staticmethod
     def fact_givens_fgft(Lap, J, t, verbosity=0, stoppingError = 0.0,
-                         errIsRel=True, order=1):
+                         errIsRel=True, order=1, enable_large_Faust=False):
         isReal = Lap.dtype in [ 'float', 'float128',
                              'float16', 'float32',
                              'float64', 'double']
@@ -1506,17 +1513,24 @@ cdef class FaustFact:
                                                                            &D_view[0], verbosity,
                                                                            stoppingError,
                                                                            errIsRel,
-                                                                           int(order))
+                                                                           int(order),
+                                                                           enable_large_Faust)
 
         core._isReal = True
         #from scipy.sparse import spdiags
         #D_spdiag = spdiags(D, [0], Lap.shape[0], Lap.shape[0])
         #return core, D_spdiag
+        if(core._isReal and core.core_faust_dbl == NULL or \
+           not core._isReal and core.core_faust_cplx == NULL):
+            raise Exception("Empty transform (nGivens is too big ? Set"
+                            " enable_large_Faust to True to force the computation).")
+
+
         return core, D
 
     @staticmethod
     def fact_givens_fgft_cplx(Lap, J, t, verbosity=0, stoppingError = 0.0,
-                         errIsRel=True, order=1):
+                         errIsRel=True, order=1, enable_large_Faust=False):
         isReal = Lap.dtype in [ 'float', 'float128',
                              'float16', 'float32',
                              'float64', 'double']
@@ -1541,22 +1555,29 @@ cdef class FaustFact:
 
         core = FaustCore(core=True)
         core.core_faust_cplx = FaustCoreCy.fact_givens_fgft_cplx[complex,double](&Lap_view[0,0],
-                                                                    Lap_num_rows,
-                                                                    Lap_num_cols, J, t,
-                                                                    &D_view[0], verbosity,
-                                                                    stoppingError,
-                                                                    errIsRel,
-                                                                    int(order))
+                                                                                 Lap_num_rows,
+                                                                                 Lap_num_cols, J, t,
+                                                                                 &D_view[0], verbosity,
+                                                                                 stoppingError,
+                                                                                 errIsRel,
+                                                                                 int(order),
+                                                                                 enable_large_Faust)
 
         core._isReal = False
         #from scipy.sparse import spdiags
         #D_spdiag = spdiags(D, [0], Lap.shape[0], Lap.shape[0])
         #return core, D_spdiag
+        if(core._isReal and core.core_faust_dbl == NULL or \
+           not core._isReal and core.core_faust_cplx == NULL):
+            raise Exception("Empty transform (nGivens is too big ? Set"
+                            " enable_large_Faust to True to force the computation).")
+
         return core, D
 
     @staticmethod
     def fact_givens_fgft_sparse_cplx(Lap, J, t, verbosity=0, stoppingError = 0.0,
-                                errIsRel=True, order='ascend'):
+                                     errIsRel=True, order='ascend',
+                                     enable_large_Faust=False):
         from scipy.sparse import spdiags
         cdef complex[:] data1d #only for csr mat factor
         cdef int [:] indices # only for csr mat
@@ -1577,39 +1598,46 @@ cdef class FaustFact:
         D_view = D
 
         core = FaustCore(core=True)
-        core.core_faust_cplx = FaustCoreCy.fact_givens_fgft_sparse_cplx[complex, double](&data1d[0],
-                                                                   &indices[0],
-                                                                   &indptr[0],
-                                                                   Lap.nnz,
-                                                                   Lap_num_rows,
-                                                                   Lap_num_cols, J, t,
-                                                                   &D_view[0],
-                                                                   verbosity,
-                                                                   stoppingError,
-                                                                   errIsRel,
-                                                                   int(order))
+        core.core_faust_cplx = \
+        FaustCoreCy.fact_givens_fgft_sparse_cplx[complex, double](&data1d[0],
+                                                                  &indices[0],
+                                                                  &indptr[0],
+                                                                  Lap.nnz,
+                                                                  Lap_num_rows,
+                                                                  Lap_num_cols, J, t,
+                                                                  &D_view[0],
+                                                                  verbosity,
+                                                                  stoppingError,
+                                                                  errIsRel,
+                                                                  int(order),
+                                                                  enable_large_Faust)
 
         core._isReal = False
         #D_spdiag = spdiags(D, [0], Lap.shape[0], Lap.shape[0])
         #return core, D_spdiag
+        if(core._isReal and core.core_faust_dbl == NULL or \
+           not core._isReal and core.core_faust_cplx == NULL):
+            raise Exception("Empty transform (nGivens is too big ? Set"
+                            " enable_large_Faust to True to force the computation).")
+
         return core, D
 
     @staticmethod
-    def eigtj(M, maxiter=None, tol=0, relerr=True,  nGivens_per_fac=None, verbosity=0,
-          order='ascend'):
-        if(maxiter == None): 
+    def eigtj(M, nGivens=None, tol=0, relerr=True,  nGivens_per_fac=None, verbosity=0,
+          order='ascend', enable_large_Faust=False):
+        if(nGivens == None):
             if(tol == 0):
-                raise Exception("You must specify maxiter or tol argument"
+                raise Exception("You must specify nGivens or tol argument"
                         " (to define a stopping  criterion)")
-            maxiter = 0 
+            nGivens = 0
         if(nGivens_per_fac == None): nGivens_per_fac = int(M.shape[0]/2)
         if(isinstance(M, np.ndarray) and \
             not np.allclose(np.matrix(M, copy=False).H, M) or M.shape[0] != M.shape[1]):
             raise ValueError(" the matrix/array must be symmetric or hermitian.")
-        if(not isinstance(maxiter, int)): raise TypeError("maxiter must be a int")
+        if(not isinstance(nGivens, int)): raise TypeError("nGivens must be a int")
         if(not isinstance(nGivens_per_fac, int)): raise TypeError("nGivens_per_fac must be a int")
         nGivens_per_fac = max(nGivens_per_fac, 1)
-        if(maxiter > 0): nGivens_per_fac = min(nGivens_per_fac, maxiter)
+        if(nGivens > 0): nGivens_per_fac = min(nGivens_per_fac, nGivens)
         tol *= tol # the C++ impl. works on squared norms to measure errors
 
         M_is_real = M.dtype in [ 'float', 'float128',
@@ -1618,26 +1646,26 @@ cdef class FaustFact:
         if(isinstance(M, np.ndarray)):
             M = np.asfortranarray(M)
             if(M_is_real):
-                core_obj,D = FaustFact.fact_givens_fgft(M, maxiter, nGivens_per_fac,
+                core_obj,D = FaustFact.fact_givens_fgft(M, nGivens, nGivens_per_fac,
                         verbosity, tol,
-                        relerr, order)
+                        relerr, order, enable_large_Faust)
             else: #complex
-                core_obj,D = FaustFact.fact_givens_fgft_cplx(M, maxiter, nGivens_per_fac,
+                core_obj,D = FaustFact.fact_givens_fgft_cplx(M, nGivens, nGivens_per_fac,
                         verbosity, tol,
-                        relerr, order)
+                        relerr, order, enable_large_Faust)
 
         elif(isinstance(M, csr_matrix)):
             if(M_is_real):
-                core_obj,D = FaustFact.fact_givens_fgft_sparse(M, maxiter,
+                core_obj,D = FaustFact.fact_givens_fgft_sparse(M, nGivens,
                         nGivens_per_fac,
                         verbosity,
                         tol,
                         relerr,
-                        order)
+                        order, enable_large_Faust)
             else: #complex
-                core_obj,D = FaustFact.fact_givens_fgft_sparse_cplx(M, maxiter, nGivens_per_fac,
+                core_obj,D = FaustFact.fact_givens_fgft_sparse_cplx(M, nGivens, nGivens_per_fac,
                         verbosity, tol,
-                        relerr, order)
+                        relerr, order, enable_large_Faust)
         else:
             raise TypeError("The matrix to diagonalize must be a"
                             " scipy.sparse.csr_matrix or a numpy array.")