Overload of slice in TransformHelperPoly.

It's needed because the factors must be created (in lazy instantiation) before the slice can happen.

Other minor changes.

src/faust_linear_operator/CPU/faust_MatSparse.hpp

@@ -1249,7 +1249,7 @@ void Faust::MatSparse<FPP,Cpu>::vstack(MatSparse<FPP, Cpu>& top, MatSparse<FPP,
 		resize(nnz, nrows, ncols);
 	typedef Eigen::Triplet<FPP> T;
 	std::vector<T> tripletList;
-	tripletList.reserve(top.getNonZeros()+bottom.getNonZeros());
+	tripletList.reserve(tnnz+bnnz);
 	for(faust_unsigned_int i=0; i < top.getNbRow() ; i++)
 		for(typename Eigen::SparseMatrix<FPP,Eigen::RowMajor>::InnerIterator it(top.mat,i); it; ++it)
 		{

src/faust_linear_operator/CPU/faust_TransformHelperPoly.h

@@ -104,7 +104,8 @@ namespace Faust
 			TransformHelper<FPP,Cpu>* optimize_time_Fv(const bool transp=false, const bool inplace=false, const int nsamples=1);
 			TransformHelper<FPP,Cpu>* swap_cols(const faust_unsigned_int id1, const faust_unsigned_int id2, const bool permutation=false, const bool inplace=false, const bool check_transpose=true);
 			TransformHelper<FPP,Cpu>* swap_rows(const faust_unsigned_int id1, const faust_unsigned_int id2, const bool permutation=false, const bool inplace=false, const bool check_transpose=true);
-
+            TransformHelper<FPP, Cpu>* slice(faust_unsigned_int start_row_id, faust_unsigned_int end_row_id,
+					                    faust_unsigned_int start_col_id, faust_unsigned_int end_col_id);
 			Vect<FPP,Cpu> multiply(const Vect<FPP,Cpu> &x, const bool transpose=false, const bool conjugate=false);
 			Vect<FPP,Cpu> multiply(const FPP* x, const bool transpose=false, const bool conjugate=false);
 			void multiply(const FPP* x, FPP* y, const bool transpose=false, const bool conjugate=false);

src/faust_linear_operator/CPU/faust_TransformHelperPoly.hpp

@@ -933,5 +933,16 @@ namespace Faust
 			return ret;
 		}
 
+	template<typename FPP>
+		TransformHelper<FPP, Cpu>* TransformHelperPoly<FPP>::slice(faust_unsigned_int start_row_id, faust_unsigned_int end_row_id,
+				faust_unsigned_int start_col_id, faust_unsigned_int end_col_id)
+		{
+			auto this_ = const_cast<TransformHelperPoly<FPP>*>(this);
+			this_->basisChebyshev_fact_all();
+			auto ret = TransformHelperGen<FPP, Cpu>::slice(start_row_id, end_row_id, start_col_id, end_col_id);
+			if(this_->laziness == INSTANTIATE_COMPUTE_AND_FREE)
+				this_->basisChebyshev_free_fact_all();
+			return ret;
+		}
 
 }

src/faust_linear_operator/faust_MatGeneric.hpp

@@ -164,7 +164,7 @@ template<typename FPP,FDevice DEVICE>
 std::string Faust::MatGeneric<FPP,DEVICE>::to_string(int32_t nrows, int32_t ncols, bool transpose, Real<FPP> density, int32_t nnz, bool is_identity, MatType type)
 {
 	std::ostringstream str;
-	str << " (" << MatGeneric<FPP,DEVICE>::get_scalar_type_str() << ")";
+	str << " (" << MatGeneric<FPP,DEVICE>::get_scalar_type_str() << ") ";
 	if(type == Dense)
 		str << "DENSE,";
 	else if(type == Sparse)

src/faust_linear_operator/faust_TransformHelperGen.h

@@ -91,7 +91,7 @@ namespace Faust
 			int get_Fv_mul_mode() const;
 			void eval_sliced_Transform();
 			void eval_fancy_idx_Transform();
-			TransformHelper<FPP, DEV>* slice(faust_unsigned_int start_row_id, faust_unsigned_int end_row_id,
+			virtual TransformHelper<FPP, DEV>* slice(faust_unsigned_int start_row_id, faust_unsigned_int end_row_id,
 					faust_unsigned_int start_col_id, faust_unsigned_int end_col_id);
 			TransformHelper<FPP, DEV>* fancy_index(faust_unsigned_int* row_ids, faust_unsigned_int num_rows, faust_unsigned_int* col_ids, faust_unsigned_int num_cols);
 			virtual TransformHelper<FPP,DEV>* optimize_storage(const bool time=true);

wrapper/python/pyfaust/poly.py

@@ -526,7 +526,7 @@ def expm_multiply(A, B, t, K=10, dev='cpu'):
     else:
         return Y
 
-def invm_multiply(A, B, rel_err=1e-6, max_K=2048, dev='cpu'):
+def invm_multiply(A, B, rel_err=1e-6, max_K=np.inf, dev='cpu'):
     """
     Computes an approximate of the action of the matrix inverse of A on B using Chebyshev polynomials.
 
@@ -572,7 +572,6 @@ def invm_multiply(A, B, rel_err=1e-6, max_K=2048, dev='cpu'):
     g = abs(1 / c + sqrt(1/c**2 - 1))
     K = min(max_K, int(((log(1/eps) + log(2/(m*sqrt(1-c**2))) - log(g-1)) /
                         log(g))))
-    print("K=", K)
     Abar = 2*A/(b-a) - (b+a)*Id/(b-a)
     T = basis(Abar, K, 'chebyshev')
     coeffs = array([ 2 / (m*sqrt(1-c**2)) * (-1)**k * g**(-k) for k in